Planetary Nebulae and How to Find Them: Color Identification in Big Broadband Surveys

TL;DR

This paper develops color-color space criteria using synthetic models to identify planetary nebulae in large broadband surveys, aiding in their census and study across the galaxy and beyond.

Contribution

It introduces a grid of synthetic magnitudes for PNe at various evolutionary stages and evaluates their detectability in upcoming surveys like LSST.

Findings

Color-color diagrams effectively distinguish PNe from other sources.

Synthetic magnitudes match well with real PNe colors from SDSS.

LSST can detect most PNe in the Magellanic Clouds and some in Andromeda.

Abstract

Planetary nebulae (PNe) provide tests of stellar evolution, can serve as tracers of chemical evolution in the Milky Way and other galaxies, and are also used as a calibrator of the cosmological distance ladder. Current and upcoming large scale photometric surveys have the potential to complete the census of PNe in our galaxy and beyond, but it is a challenge to disambiguate partially or fully unresolved PNe from the myriad other sources observed in these surveys. Here we carry out synthetic observations of nebular models to determine $ugrizy$ color-color spaces that can successfully identify PNe among billions of other sources. As a primary result we present a grid of synthetic absolute magnitudes for PNe at various stages of their evolution, and we make comparisons with real PNe colors from the Sloan Digital Sky Survey. We find that the $r-i$ versus $g-r$, and the $r-i$ versus $u-g$, color-color diagrams show the greatest promise for cleanly separating PNe from stars, background galaxies, and quasars. Finally, we consider the potential harvest of PNe from upcoming large surveys. For example, for typical progenitor host star masses of $\sim$3 M$_\odot$, we find that the Large Synoptic Survey Telescope (LSST) should be sensitive to virtually all PNe in the Magellanic Clouds with extinction up to $A_{\rm V}$ of $\sim$5 mag; out to the distance of Andromeda, LSST would be sensitive to the youngest PNe (age less than $\sim$6800 yr) and with $A_{\rm V}$ up to 1 mag.

Tables

Table 1. Table 1 Parameters used for the 13 Cloudy models that change from model to model. These parameters represent the evolution of the central star ( T ⋆ subscript 𝑇 ⋆ T_{\star} , L ⋆ subscript 𝐿 ⋆ L_{\star} ) and nebula ( R i ​ n subscript 𝑅 𝑖 𝑛 R_{in} , n e subscript 𝑛 𝑒 n_{e} ) while M s ​ t ​ o ​ p subscript 𝑀 𝑠 𝑡 𝑜 𝑝 M_{stop} is used to stop the simulation once the nebula has accumulated enough mass.

$Age$	$\log T_{\star }$	$\log \frac{L_{\star }}{L_{\odot }}$	$R_{in}$	$\log n_e$	$M_{stop}$
(yrs)	(K)	($L_{\odot }$)	(pc)	($c{m}^{-3}$)	($M_{\odot }$)
3502	4.786	3.782	0.035	4.386	0.071
4154	4.887	3.766	0.044	4.159	0.088
4860	4.988	3.736	0.056	3.911	0.113
5663	5.089	3.671	0.073	3.652	0.146
6328	5.155	3.571	0.088	3.460	0.178
6745	5.183	3.468	0.098	3.349	0.199
7031	5.194	3.363	0.106	3.277	0.201
7224	5.195	3.257	0.111	3.230	0.201
7387	5.181	3.056	0.115	3.192	0.202
7453	5.165	2.912	0.117	3.176	0.202
7552	5.140	2.73	0.119	3.153	0.203
7751	5.114	2.566	0.125	3.108	0.202
8351	5.095	2.452	0.142	2.980	0.202

Table 2. Table 2 Abundances used for PN models within Cloudy (Aller & Czyzak, 1983 ; Khromov, 1989 ) . Values are relative to solar values from Grevesse & Sauval ( 1998 ) , Holweger ( 2001 ) , Allende Prieto et al. ( 2001 , 2002 ) . Elements not mentioned are assumed to be depleted enough to be of no consequence. For more information on available abundances within Cloudy refer to “Hazy” the Cloudy manual (Ferland et al., 2013 ) .

Atom	[X/H]
	(dex)
H	0.0000
He	-1.0000
C	-3.1079
N	-3.7447
O	-3.3565
F	-6.5229
Ne	-3.9586
Na	-5.7212
Mg	-5.7959
Al	-6.5686
Si	-5.0000
P	-6.6990
S	-5.0000
Cl	-6.7696
Ar	-5.5686
K	-6.9208
Ca	-7.9208
Fe	-6.3010
Ni	-7.7447

Table 3. Table 3 Properties of the u ​ g ​ r ​ i ​ z ​ y 𝑢 𝑔 𝑟 𝑖 𝑧 𝑦 ugrizy filter set used in our analysis. All magnitudes are normalized to the AB magnitude system (Oke & Gunn, 1983 ) .

Filter	Range	${\lambda }_{\mathrm{eff}}$	$W_{\mathrm{eff}}$	$m_{sat}$	$m_{5\sigma }$
	(nm)	(nm)	(nm)	(mags)	(mags)
$u$	320–400	373.2	54.6	14.7	23.9
$g$	400–552	473.0	133.2	15.7	25.0
$r$	552–691	613.8	133.7	15.8	24.7
$i$	691–818	748.7	832.5	15.8	24.0
$z$	818–922	866.8	937.5	15.3	23.3
$y$aafootnotemark:	950–1080	967.6	81.0	13.9	22.1

Table 4. Table 4 Predicted absolute magnitudes of our synthetic nebular spectra for the u ​ g ​ r ​ i ​ z ​ y 𝑢 𝑔 𝑟 𝑖 𝑧 𝑦 ugrizy filter set. No extinction is applied. Models are fully unresolved with Ω aperture / Ω neb = subscript Ω aperture subscript Ω neb absent \Omega_{\rm aperture}/\Omega_{\rm neb}= 1.3 and f n ​ e ​ b = subscript 𝑓 𝑛 𝑒 𝑏 absent f_{neb}= 1. For varying Ω aperture / Ω neb subscript Ω aperture subscript Ω neb \Omega_{\rm aperture}/\Omega_{\rm neb} and f n ​ e ​ b subscript 𝑓 𝑛 𝑒 𝑏 f_{neb} values see Table LABEL:table:_neb_fraction .

Variant I:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=1.3$			$f_{neb}$ = $1$
Age	$u$	$g$	$r$	$i$	$z$	$y$
(yrs)	(mags)	(mags)	(mags)	(mags)	(mags)	(mags)
3502	-2.59	-3.63	-3.38	-2.20	-2.18	-2.09
4154	-2.70	-3.91	-3.47	-2.21	-2.20	-2.15
4860	-2.46	-3.94	-2.98	-1.87	-1.98	-1.94
5663	-2.11	-3.66	-2.53	-1.44	-1.62	-1.66
6328	-1.93	-3.41	-2.36	-1.21	-1.40	-1.49
6745	-1.95	-3.23	-2.43	-1.14	-1.29	-1.39
7031	-1.98	-3.00	-2.53	-1.06	-1.13	-1.28
7224	-2.00	-2.72	-2.64	-0.98	-0.98	-1.21
7387	-1.71	-2.10	-2.49	-0.63	-0.54	-0.95
7453	-1.44	-1.60	-2.27	-0.33	-0.18	-0.72
7552	-1.09	-0.90	-1.98	0.08	0.28	-0.45
7751	-0.75	-0.22	-1.68	0.45	0.71	-0.24
8351	-0.51	0.24	-1.46	0.73	1.00	-0.11

Table 5. Table 5 Absolute magnitudes for all 13 PN models for corresponding Ω aperture / Ω neb subscript Ω aperture subscript Ω neb \Omega_{\rm aperture}/\Omega_{\rm neb} and f n ​ e ​ b subscript 𝑓 𝑛 𝑒 𝑏 f_{neb}

Variant I:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=1.3$			$f_{neb}$ = $1$
Age	u	g	r	i	z	y
(yrs)	(mags)	(mags)	(mags)	(mags)	(mags)	(mags)
3502	-2.59	-3.63	-3.38	-2.20	-2.18	-2.09
4154	-2.70	-3.91	-3.47	-2.21	-2.20	-2.15
4860	-2.46	-3.94	-2.98	-1.87	-1.98	-1.94
5663	-2.11	-3.66	-2.53	-1.44	-1.62	-1.66
6328	-1.93	-3.41	-2.36	-1.21	-1.40	-1.49
6745	-1.95	-3.23	-2.43	-1.14	-1.29	-1.39
7031	-1.98	-3.00	-2.53	-1.06	-1.13	-1.28
7224	-2.00	-2.72	-2.64	-0.98	-0.98	-1.21
7387	-1.71	-2.10	-2.49	-0.63	-0.54	-0.95
7453	-1.44	-1.60	-2.27	-0.33	-0.18	-0.72
7552	-1.09	-0.90	-1.98	0.08	0.28	-0.45
7751	-0.75	-0.22	-1.68	0.45	0.71	-0.24
8351	-0.51	0.24	-1.46	0.73	1.00	-0.11
Variant II:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.9$			$f_{neb}$ = $4\times {10}^{-1}$
3502	-1.88	-2.78	-2.52	-1.39	-1.35	-1.24
4154	-1.90	-3.03	-2.58	-1.35	-1.33	-1.27
4860	-1.62	-3.04	-2.09	-0.99	-1.10	-1.05
5663	-1.25	-2.76	-1.63	-0.56	-0.73	-0.77
6328	-1.06	-2.51	-1.46	-0.32	-0.50	-0.59
6745	-1.06	-2.33	-1.54	-0.24	-0.39	-0.49
7031	-1.09	-2.10	-1.63	-0.16	-0.24	-0.38
7224	-1.10	-1.82	-1.74	-0.08	-0.08	-0.31
7387	-0.82	-1.20	-1.59	0.27	0.36	-0.05
7453	-0.55	-0.70	-1.37	0.57	0.72	0.19
7552	-0.20	0.00	-1.07	0.97	1.18	0.45
7751	0.14	0.67	-0.78	1.35	1.60	0.66
8351	0.38	1.13	-0.56	1.62	1.89	0.79
Variant III:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.7$			$f_{neb}$ = $2\times {10}^{-1}$
3502	-1.45	-2.17	-1.90	-0.83	-0.77	-0.65
4154	-1.35	-2.38	-1.93	-0.74	-0.70	-0.63
4860	-1.03	-2.38	-1.43	-0.36	-0.45	-0.39
5663	-0.63	-2.09	-0.97	0.10	-0.07	-0.10
6328	-0.41	-1.84	-0.79	0.34	0.17	0.08
6745	-0.41	-1.66	-0.86	0.43	0.29	0.18
7031	-0.42	-1.42	-0.95	0.51	0.44	0.30
7224	-0.44	-1.14	-1.06	0.59	0.60	0.37
7387	-0.15	-0.53	-0.91	0.94	1.04	0.63
7453	0.12	-0.02	-0.69	1.24	1.39	0.86
7552	0.47	0.67	-0.39	1.64	1.85	1.13
7751	0.81	1.34	-0.10	2.02	2.27	1.34
8351	1.04	1.79	0.12	2.29	2.56	1.47
Variant IV:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.5$			$f_{neb}$ = $1\times {10}^{-1}$
3502	-1.08	-1.53	-1.24	-0.27	-0.17	-0.03
4154	-0.82	-1.65	-1.20	-0.08	-0.02	0.07
4860	-0.42	-1.61	-0.68	0.35	0.28	0.35
5663	0.05	-1.31	-0.20	0.83	0.69	0.67
6328	0.31	-1.05	-0.01	1.10	0.94	0.87
6745	0.34	-0.86	-0.07	1.20	1.07	0.97
7031	0.34	-0.62	-0.16	1.29	1.22	1.09
7224	0.34	-0.34	-0.26	1.38	1.39	1.17
7387	0.63	0.27	-0.11	1.72	1.82	1.42
7453	0.90	0.76	0.11	2.03	2.18	1.66
7552	1.24	1.46	0.41	2.42	2.64	1.93
7751	1.58	2.11	0.70	2.80	3.05	2.14
8351	1.81	2.55	0.92	3.07	3.34	2.27
Variant V:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.3$			$f_{neb}$ = $4\times {10}^{-2}$
3502	-0.77	-0.83	-0.47	0.27	0.45	0.63
4154	-0.31	-0.73	-0.27	0.67	0.80	0.93
4860	0.23	-0.58	0.29	1.20	1.21	1.31
5663	0.82	-0.23	0.82	1.77	1.69	1.71
6328	1.19	0.05	1.06	2.11	1.99	1.94
6745	1.30	0.25	1.03	2.24	2.15	2.07
7031	1.35	0.49	0.96	2.36	2.31	2.20
7224	1.38	0.77	0.86	2.45	2.48	2.28
7387	1.68	1.37	1.02	2.80	2.92	2.54
7453	1.94	1.86	1.24	3.10	3.26	2.78
7552	2.29	2.52	1.54	3.49	3.71	3.05
7751	2.61	3.14	1.83	3.85	4.11	3.25
8351	2.83	3.54	2.05	4.11	4.38	3.39
Variant VI:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.15$			$f_{neb}$ = $9\times {10}^{-3}$
3502	-0.62	-0.32	0.10	0.61	0.85	1.08
4154	-0.01	0.11	0.59	1.23	1.45	1.65
4860	0.66	0.51	1.24	1.90	2.06	2.23
5663	1.41	1.01	1.91	2.63	2.72	2.82
6328	1.96	1.38	2.29	3.13	3.16	3.20
6745	2.23	1.62	2.36	3.38	3.40	3.39
7031	2.40	1.88	2.35	3.56	3.60	3.56
7224	2.50	2.15	2.30	3.70	3.79	3.67
7387	2.84	2.73	2.48	4.06	4.22	3.95
7453	3.10	3.17	2.70	4.35	4.54	4.18
7552	3.42	3.75	3.00	4.71	4.95	4.45
7751	3.71	4.25	3.28	5.03	5.30	4.66
8351	3.90	4.55	3.50	5.25	5.54	4.80
Variant VII:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.1$			$f_{neb}$ = $4\times {10}^{-3}$
3502	-0.59	-0.21	0.24	0.68	0.94	1.18
4154	0.06	0.35	0.83	1.36	1.61	1.84
4860	0.76	0.89	1.54	2.09	2.30	2.52
5663	1.55	1.50	2.28	2.88	3.05	3.22
6328	2.17	1.96	2.77	3.47	3.60	3.70
6745	2.52	2.24	2.94	3.80	3.90	3.96
7031	2.75	2.51	3.00	4.03	4.14	4.16
7224	2.91	2.79	2.99	4.20	4.35	4.30
7387	3.27	3.33	3.20	4.57	4.77	4.60
7453	3.52	3.73	3.42	4.85	5.07	4.84
7552	3.83	4.22	3.71	5.18	5.44	5.11
7751	4.09	4.64	3.99	5.47	5.76	5.32
8351	4.27	4.88	4.20	5.67	5.96	5.45
Variant VIII:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.075$			$f_{neb}$ = $2\times {10}^{-3}$
3502	-0.58	-0.16	0.30	0.71	0.98	1.23
4154	0.08	0.47	0.96	1.42	1.69	1.93
4860	0.80	1.10	1.70	2.18	2.43	2.67
5663	1.62	1.81	2.50	3.01	3.24	3.45
6328	2.28	2.35	3.08	3.66	3.86	4.02
6745	2.67	2.68	3.34	4.04	4.21	4.35
7031	2.95	2.97	3.47	4.31	4.49	4.59
7224	3.15	3.24	3.52	4.52	4.72	4.76
7387	3.52	3.74	3.76	4.90	5.13	5.09
7453	3.77	4.10	3.99	5.17	5.42	5.32
7552	4.07	4.52	4.28	5.48	5.75	5.59
7751	4.32	4.86	4.55	5.74	6.03	5.80
8351	4.48	5.07	4.74	5.92	6.21	5.94
Variant IX:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.05$			$f_{neb}$ = $1\times {10}^{-3}$
3502	-0.57	-0.13	0.34	0.72	1.00	1.25
4154	0.10	0.53	1.02	1.45	1.73	1.98
4860	0.82	1.22	1.78	2.22	2.49	2.74
5663	1.65	1.98	2.61	3.08	3.34	3.57
6328	2.33	2.58	3.24	3.76	4.00	4.20
6745	2.75	2.95	3.57	4.17	4.39	4.57
7031	3.05	3.25	3.77	4.47	4.69	4.85
7224	3.28	3.53	3.88	4.70	4.93	5.05
7387	3.66	4.00	4.16	5.09	5.34	5.40
7453	3.91	4.32	4.39	5.34	5.61	5.64
7552	4.19	4.69	4.68	5.64	5.92	5.90
7751	4.43	4.98	4.93	5.89	6.18	6.12
8351	4.59	5.16	5.11	6.05	6.35	6.26
Variant X:		${\Omega }_{\mathrm{aperture}}/{\Omega }_{\mathrm{neb}}=0.01$			$f_{neb}$ = $3\times {10}^{-5}$
3502	-0.57	-0.10	0.37	0.74	1.02	1.28
4154	0.11	0.60	1.10	1.48	1.77	2.03
4860	0.85	1.36	1.88	2.27	2.57	2.83
5663	1.69	2.22	2.75	3.16	3.45	3.72
6328	2.39	2.93	3.47	3.88	4.18	4.45
6745	2.84	3.38	3.92	4.34	4.64	4.91
7031	3.17	3.72	4.26	4.68	4.98	5.25
7224	3.44	3.98	4.52	4.95	5.25	5.52
7387	3.84	4.39	4.91	5.34	5.65	5.91
7453	4.08	4.63	5.15	5.59	5.89	6.15
7552	4.36	4.91	5.43	5.86	6.16	6.42
7751	4.58	5.13	5.65	6.08	6.38	6.64
8351	4.72	5.27	5.79	6.22	6.52	6.78

Full text

1. Introduction

Planetary Nebulae (PNe) are the shells of gas ionized by hot central stars. The PN forms from previously ejected material lost during the asymptotic giant branch (AGB) phase of a low-to-intermediate mass star ($\rm 1-8~{}M_{\odot}$). The central star of a planetary nebula (CSPN) plays a key role in the PN characteristics since its fast stellar wind plows into the AGB wind to form the nebula while the CSPN’s high surface temperature photoionizes the gas in the newly-formed nebula. A PN will expand and fade over time, while the CSPN rises in temperature, until it eventually cools towards the white dwarf (WD) cooling track (Kwok et al., 1978; Vassiliadis & Wood, 1994; Bloecker, 1995; Miller Bertolami, 2016). Compared to the lifetime of the star, the PN phase is short-lived remaining visible for only $\sim 10^{4}$ years (Iben, 1995; Jacob et al., 2013).

For over 60 years the formation process of PNe has been questioned favoring two contending processes. The first consists of a single star in which heavy mass loss occurs during the AGB phase while stellar rotation and magnetic fields shape the expanding nebula (Gurzadian, 1962; Garcia-Segura et al., 1999; Garcia‐Segura et al., 2005; Matt et al., 2004; Blackman, 2004). The second process favors a wide range of interactions between the evolving star and a binary companion for shaping the nebula (Fabian & Hansen, 1979; Soker, 1997; De Marco, 2009). Theoretical and observational considerations suggest that a single AGB star is unlikely to produce a strong enough magnetic field to dramatically shape the nebula, favouring the possibility of binary interactions being mainly responsible for the formation and shaping of non-spherical PN (Soker, 2006; Nordhaus et al., 2007). It is likely that there are still many PNe left undiscovered as our best estimates place the total number of galactic PNe anywhere between $\sim 6600$ and $4.6\pm 1.3\times 10^{4}$ depending on the formation process (De Marco, 2005; Moe & De Marco, 2006).

Naturally, because they represent a specific and short-lived phase of stellar evolution, PNe can be difficult to study but they are important for improving our understanding of late-stage low- to intermediate-mass stars (Iben, 1995; Frew, 2008). PNe that have been found in the Milky Way and in neighboring galaxies have been extremely valuable for a variety of studies. Because of their bright emission lines, PNe are identifiable across the Galaxy and in nearby stellar systems. PNe can be used as tracers of galaxy kinematics (Hurley-Keller & Morrison, 2004; Peng et al., 2004; Merrett et al., 2006; Coccato & Coccato, 2016); as a rung on the distance ladder via the PN luminosity function (Ciardullo, 2003; Gesicki et al., 2018); and as potential tracers of the chemical evolution of the Milky Way and other galaxies (Walsh et al., 1997; Richer et al., 1998; Kniazev et al., 2008; Saviane et al., 2009; Magrini & Gonçalves, 2009; Kwitter et al., 2012; Sanders et al., 2012; Cavichia et al., 2017).

There is ongoing interest in techniques for readily identifying more PNe efficiently and reliably. A large number of PNe have been discovered through techniques that take advantage of their bright nebular emission lines. $\rm H\alpha$ surveys such as the SuperCOSMOS $\rm H\alpha$ Survey (SHS) (Parker et al., 2005; Frew et al., 2014), the INT Photometric $\rm H\alpha$ Survey (IPHAS) (Drew et al., 2005), and the VST Photometric $\rm H\alpha$ Survey (VPHAS+) (Drew et al., 2014) have been very successful in identifying galactic PNe and are cataloged in the Hong Kong/AAO/Strasbourg/$\rm H\alpha$ (HASH) (Parker et al., 2016) database which contains $\sim 3500$ objects. Integral field spectroscopy of [O III] $\lambda 5007$ has been used to find PNe in crowded areas (Pastorello et al., 2013). Dust can make it difficult for optical surveys to detect PNe but they have also been found in the UKIRT Wide-field Imaging Survey for $\rm H_{2}$ survey (UWISH2) through their $\rm H_{2}$ emission (Gledhill et al., 2018). More galactic PNe have recently been discovered through their multi-wavelength characteristics ranging from optical to radio emission (Fragkou et al., 2018).

Large, all-sky surveys, like the Sloan Digital Sky Survey (SDSS) and, in the future, the Large Synoptic Survey Telescope (LSST), can be used to potentially identify thousands of new PN based on broadband photometry where PNe are poorly characterized. Although the observed broadband colors of PNe have been studied (e.g., in 2MASS and WISE colors Schmeja & Kimeswenger, 2001; Frew & Parker, 2010; Iłkiewicz & Mikołajewska, 2017), and although the theoretically expected broadband colors of PN central stars have previously been calculated (e.g., Weston et al., 2009; Morris & Montez, 2015), the theoretically expected broadband colors of PN nebular emission have yet to be characterized and compared with observations. Existing PNe catalogs are far from complete, and upcoming all-sky photometric surveys have great potential to uncover troves of additional PNe, assuming that it will be possible to efficiently distinguish true PNe from the large numbers of false positives.

We focus on the optical region of the PN spectrum because of the prominent nebular emission lines at these wavelengths. A vast amount of $ugriz(y)$ survey data is or will become available through SDSS and LSST, therefore we attempt to characterize PNe in this photometric system. We consider the broadband $ugrizy$ characteristics of a synthetic PN as a function of evolutionary age, for both resolved and unresolved PNe. Our methodology for creating our synthetic PN models and calculating the synthesized $ugrizy$ observations are presented in §2. A grid of broadband absolute magnitudes and the efficacy of different color-color diagrams to reliably identify PNe are provided in §3. In §4 we compare our results to existing observational SDSS studies and consider future applications with LSST. Finally, §5 presents a summary of our conclusions.

2. Methodology

In this section, we describe our methodology for creating synthetic broadband $ugrizy$ observations of PNe as functions of evolutionary age for both resolved and unresolved PNe. We begin by describing the physical ingredients and assumptions adopted. Because we seek to undertake a more comprehensive study of the broadband behavior of PNe than previously attempted, we intentionally construct our methodology to fully and self-consistently include the effects of both the central star and the nebular evolution. Next we describe the production of the nebular spectra followed by the synthetic $ugrizy$ magnitudes and colors from these spectra, and we discuss our treatment of resolved and unresolved PNe. Finally, we include the effects of interstellar reddening for resolved and unresolved spectra.

2.1. Planetary Nebula Ingredients

To build our synthetic PN we start with the central star and nebular properties. In particular, throughout this study we work with a 3 $M_{\odot}$ progenitor star, because it is for such a star that nebular evolution models exist to which we can self-consistently couple the evolution of the central star and thereby produce synthetic spectra and colors over the course of the evolution of the star+PN. Figure 1 depicts the HR diagram for a $3\ M_{\odot}$ central star evolutionary track from the models of Vassiliadis & Wood (1994); Bloecker (1995), along with $1\ M_{\odot}$ and $5\ M_{\odot}$ progenitor masses for context. The final mass of the white dwarf in this model is ${0.605~{}M_{\odot}}$, similar to the $\approx$0.6 M⊙ mass that is often assumed as a fiducial mass for PN central stars (Perinotto et al., 2004). Thirteen specific age positions on the evolutionary track are indicated, representing the discrete ages that we have chosen to include in our study, spanning the evolution of the system from emergence of the PN to exposing the white dwarf. Figure 2 (top panel) indicates the temporal behavior of the CSPN luminosity ($L$) and temperature ($T_{\rm eff}$) corresponding to these 13 ages.

For every age position on the central star’s evolutionary track (Figure 1) we require a nebular emission model corresponding to that evolutionary age and with the appropriate central star mass. Fortunately, 1D radiation-hydrodynamics simulations of a PN have been performed by Perinotto et al. (2004); Schönberner et al. (2005). Because those simulations were done for the expansion of a PN illuminated by a central remnant of 0.595 $M_{\odot}$, it is essential to use a central-star evolutionary model that produces a remnant as close to this same final mass as possible. That is why our choice above is of a central star whose final mass is most similar to this value111The slight discrepancy between the central star final masses of 0.595 $M_{\odot}$ versus 0.605 $M_{\odot}$ in the nebular versus central star evolutionary tracks, respectively, is negligible for compact nebula with radii $<1{\rm~{}pc}$ (Jacob et al., 2013)., corresponding to a progenitor mass of 3 $M_{\odot}$.

We calculate the nebular size at each time step using a simple approximation based on the Perinotto et al. (2004); Schönberner et al. (2005) nebular model. Based on the radiation-hydrodynamic simulations, the densest part of the nebula is coincident with the inner radius (or rim), $R_{\rm in}$, which follows the relationship $R=R_{0}+c(t-t_{0})^{b}$, where $R_{0}=0.031$ pc, $t_{0}=3016$ yrs, $c=6.7\times 10^{-7}{\rm~{}pc~{}yr}^{-1}$, and $b=1.4$. For each time step in the evolutionary track, we determine the inner radius ($R_{\rm in}$) based on this prescription. Since the prescription is only applicable for $3016<t<10000$ yrs, we only considered the evolutionary track between these ages. As a result, the radii of our synthetic nebula model expands from $0.035$ pc to $0.142$ pc.

The outer radius $R_{\rm out}$ of the nebula is determined by Cloudy calculations (see next section), specifically by a mass stopping criteria for each model. Phillips (2007) describe an empirical relationship to the total ionized mass, $M_{\rm ionized}$. $M_{\rm ionized}$ increases linearly for $R\leq 0.1{\rm~{}pc}$ according to $M_{\rm ionized}=M_{0}(R/0.1{\rm~{}pc})$ until reaching a constant total mass, $M_{0}$, for $R>0.1{\rm~{}pc}$. In our radiative transfer calculations, we restrict the total ionized mass of the nebula according to this empirical relationship using $M_{0}\approx 0.2\ M_{\odot}$. As a result the final mass of our synthetic nebula are designed to closely match the final masses of the ionized nebula described in Phillips (2007). Also as a consequence of this stopping mass criteria the thickness of the nebular shell remains thin, specifically, $(R_{\rm out}-R_{\rm in})/R_{\rm in}$ is always $<$ 0.25.

In addition to the inner radius, for each considered position on the evolutionary track we must calculate the electron density, $n_{e}$, of the nebula. Based on the behavior of $\sim 240$ Galactic and Magellanic PNe Frew (2008) derived a nebular density relation with nebular radius of $\log n_{e}=-2.31(\pm 0.04)\log R+1.02(\pm 0.04)$, where $n_{e}$ is the electron density of the nebula and $R$ is the nebular radius, in our case $R_{\rm in}$. Figure 2 (bottom panel) shows the evolution of $n_{e}$ and $R_{\rm in}$ over all ages of the nebula.

2.2. Cloudy, With a Chance of Photons

The central star and nebular properties described in the previous section are used as input to Cloudy radiative transfer calculations. Cloudy is a plasma simulation software that simulates non-equilibrium gas conditions to predict an observable spectrum(version 13.03; Ferland et al. (2013)).

In Table 1, we list the stellar and nebular ingredients needed to run our Cloudy radiative transfer models. For the nebular composition we use chemical abundances typical of a planetary nebula (see Table 2 Ferland et al., 2013; Aller & Czyzak, 1983; Khromov, 1989). We assume spherical geometry for the nebula, a uniform filling factor of unity, and no stopping temperature criteria. Distance between the observer and the nebula is assumed to be 1 kpc. The resulting coarse (R = 200) and high-resolution (R = 2000) spectra at the youngest nebular age considered ($3502$ yrs, $61094$ K, $6053~{}L_{\odot}$) are shown in Figure 3 for our region of interest ($300{\rm~{}to~{}}1100{\rm~{}nm}$). The high-resolution spectrum is used for display purposes. Since the high-resolution spectrum does not preserve flux, we use the coarse spectrum for our photometric calculations.

2.3. Synthesizing Broadband Photometric Observations

We consider observations of our synthetic nebulae performed with a standard $ugrizy$ photometric system222We used the LSST filters as defined on 2016-12-07 obtained from the Spanish Virtual Observatory (SVO) Filter Profile Service. The $y$ filter is shown for completeness but is not used in this study.. The set of photometric filters cover wavelengths between 320 nm and 1080 nm and their properties are provided in Table 3 with their transmission curves shown in Figure 3. The filter transmission curve is interpolated onto the wavelength grid of the coarse resolution Cloudy spectra described in §2.2. For each age of our synthetic nebula we convolved the nebular and central star flux spectra with the filter transmission curves. Then we normalized the resulting flux to the AB magnitude system to determine the $ugrizy$ magnitudes. We scaled these magnitudes to 10 pc to determine the absolute magnitudes, which are presented in Table 4.

2.4. Accounting for Spatially Resolved Extended Emission

The extended nature of a PN requires additional consideration. For aperature photometry, when the nebula and central star are resolved only a fraction of the nebular flux will be measured. To understand the impact of resolved nebula and central stars, we constructed a toy model of the synthetic nebulae. For each stage, we assumed a three-dimensional spherical shell with the radius, thickness, and density used in our synthetic nebulae. We projected this shell onto the plane of the sky. Assuming aperture photometry, we considered a range of aperture diameters on the sky, $\Omega_{\rm aperture}$, and nebular diameter on the sky, $\Omega_{\rm neb}$. For a range of aperture to nebular diameter ratios, $\Omega_{\rm aperture}/\Omega_{\rm neb}$, we calculated the fraction of nebular emission, $f_{\rm neb}$, that is measured by the aperture. We recalculated the magnitudes after scaling the nebular flux by this fraction and adding it to the central star flux. We consider models with the same nebular fraction as a fractional variant and label them with roman numerals in increasing order as $f_{\rm neb}$ decreases.

This method of estimating resolved nebular magnitudes provides the most flexibility with regards to survey resolution characteristics. Note that because in reality a given PN, at a given distance, will evolve in its apparent angular size, it will not in general evolve along a single variant “track”. Rather, the collection of variant tracks represent the overall parameter space through which individual PNe may evolve, for a range of distances, photometric apertures, and states of evolution.

2.5. Reddening

Intervening dust poses a serious problem for observing and identifying objects with photometric colors. The dust will dim (extinction) and redden an object which can limit the distance at which an object can be detected and change its position in a color-color diagram. To determine the reddening vector in the ugrizy system, we applied the reddening curve of Cardelli et al. (1989) to all ages of the synthetic nebular spectra then calculated the reddened magnitudes. Figure 4 shows how varying levels of extinction ($\rm A_{\rm V}$) can affect the spectrum of a PN. Figures 5 and 6 show the reddening vector that all 13 models would follow with $A_{\rm V}=1{\rm~{}mag}$.

3. Results

In this section, we present the main results of our study. First, we provide a grid of synthetic $ugrizy$ photometry representing our model PN at different evolutionary ages and different observed angular sizes. Second, we describe the evolution of prominent emission lines in each filter and their effect on the calculated magnitudes. We then identify the $ugrizy$ color-color spaces that are most effective at differentiating PNe from other celestial objects.

3.1. Synthesized Nebula Absolute Magnitudes

A key result of this study is the collection of the absolute magnitudes for our synthetic nebulae in the $ugrizy$ photometric system (see Table 4, Table LABEL:table:_neb_fraction and Figure 5 for all variants). These absolute magnitudes span $\sim 4800$ years of PN evolution for a $\rm 3~{}M_{\odot}$ progenitor star evolving into a $\rm 0.605~{}M_{\odot}$ white dwarf. The contributions from continuum and line emission dictate the behavior of the broad band magnitudes as a function of age. In each band, the continuum drops by about an order of magnitude as the models approach older stages resulting in the overall decline in the broadband magnitudes as the nebula evolves. We do not display variants beyond VII, as variants VIII–X become too faint to be practical.

The line emission as a function of age in each band is often dictated by one or two species that dominate a given bandpass. In the $u$ band, the prominent emission lines are: [Ne III] ($\lambda 3869$, $\lambda 3968$), [O II] ($\lambda 3726$, $\lambda 3729$), H II ($\lambda 3835$, $\lambda 3889$, $\lambda 3970$), and He III ($\lambda 3203$). In the $g$ band, the prominent emission lines are: [O III] ($\lambda 5007$, $\lambda 4959$), H I ($\rm H_{\rm\beta}$ $\lambda 4861$, $\rm H_{\rm\gamma}$ $\lambda 4340$, $\rm H_{\rm\delta}$ $\lambda 4102$), He II ($\lambda 4686$), He I ($\lambda 4471$), [S II] ($\lambda 4074$, $\lambda 4070$, $\lambda 4078$), and N I ($\lambda 5200$, $\lambda 5198$). In the $r$ band, the prominent emission lines are: H I ($\rm H_{\rm\alpha}$ $\lambda 6563$), [N II] ($\lambda 6584$, $\lambda 6548$), He I ($\lambda 5876$), [S II] ($\lambda 6720$, $\lambda 6731$, $\lambda 6716$), S III ($\lambda 6312$), [Cl III] ($\lambda 5538$) and [O I] ($\lambda 6300$). In the $i$ band, the prominent emission lines are: [Ar III] ($\lambda 7135$, $\lambda 7751$), [O II] ($\lambda 7323$, $\lambda 7332$), He I ($\lambda 7065$), and Cl IV ($\lambda 8047$). Fluxes from the aforementioned lines exhibit three common behaviors: (1) decreasing throughout the evolution, (2) double-peaked with an initial rise, a decrease, then another rise as the CSPN $T_{\rm eff}$ approaches its hottest temperature before decreasing again as $T_{\rm eff}$ cools, and, (3) single-peak with a steady increase in flux as the as CSPN $T_{\rm eff}$ approaches its hottest temperature before decreasing as $T_{\rm eff}$ cools. Lines with decreasing behavior are: [Ne III], H I, [O III], He I, S III, [Cl III], Cl IV, [Ar III]. Lines with a double peak behavior are: [O II], [S II], N I, [N II], [O I]. Only He II lines show a single peak behavior. These behaviors of the line emission are due to the ionization and physical parameters (radius and density) of the modeled nebula.

When convolved with the bandpass transmission curves the broadband magnitudes can mimic the age-related behavior of the strongest emission lines. The $u$ band magnitude behavior is double-peaked with a steady decrease largely influenced by the [O II] emission and the overall decreasing flux from the other lines in this bandpass. The [Ne III] lines are very close to the edge of the $u$ filter bandpass, where the transmission is reduced, therefore, despite being the most prominent line in the $u$ band, [Ne III] does not contribute as much flux as the [O II] lines. The $g$ band magnitude exhibits decreasing behavior and is largely influenced by [O III]. The $r$ band magnitude is double-peaked, which is the result of the H I emission dominating at earlier ages ($t_{\rm age}<7000~{}{\rm yr}$) until [N II] begins to dominate. The $i$ band magnitude decreases according to the evolution of [Ar III] with some influence from [O II].

In the CMDs (Figure 5) there exists a hook feature before the eighth model ($t_{\rm age}\leq 7224$ yrs) (also somewhat visible in the CCDs) due to the evolution of these lines. Once the models have reached the empirically observed maximum mass we imposed (at the eighth model $t_{\rm age}=7224$ yrs) they become matter bounded and transparent to ionizing radiation resulting in a somewhat steady decline in brightness and color evolution.

The absolute magnitudes given in Table 4 are only valid for the total integrated flux from a PN. When a PN is resolved, only a fraction of the nebula will be observed, thus reducing the nebular contribution. This effect changes the measured flux in each filter and the measured color for any two filters. Figures 5, 6, and 8 show how magnitudes and colors change as the aperture encloses less of the nebula.

3.2. Identifying PNe in $ugrizy$ Color-Color Diagrams

We compare our PN models and variants to SDSS objects in Figure 6. Contours of stars, galaxies, and quasars are used to illustrate their locations in the various diagrams. Data for these object types come from SDSS DR7 (Abazajian et al., 2009) and were queried using the SDSS SkyServer. We only included objects with SDSS spectroscopic observations to allow for further filtering on $z$ for the star, galaxy, and quasar object types.

The top panel of Figure 6 shows the $r-i$ vs $g-r$ color-color diagram. Here, models occupy a separated space from the SDSS data with $r-i$ between 0 and $-2.5$ and $g-r$ between $\sim$1.75 and $-1.5$. The unresolved synthetic nebular models (variant I), as well as the oldest resolved models (variants VI through IX) remain separated from stars, quasars, and galaxies through out all or most stages of evolution in this diagram for $r-i<-0.75$. (We do not show variants II–V as these are effectively equivalent to variant I.) As the central star and nebula are resolved, the colors converge towards the locus of blue evolved stars, as expected. Similarly, the bottom panel of Figure 6 shows that our models are also separated from most of the SDSS data in the $r-i$ vs $u-g$ color-color diagram for $r-i<-0.75$. These two diagrams will be essential in identifying new PN candidates of all ages.

In the middle panel of Figure 6, the PN models occupy the left region of the diagram. They fall between values of $\sim$1.75 to $-1.5$ in $g-r$ and $\sim$1.6 to $-1.0$ in $u-g$. Our synthetic nebulae cross through the loci of early-type MS stars, quasars, and blue evolved stars in this diagram. The youngest ($t_{\rm age}<7000$ yrs) and oldest ($t_{\rm age}>7500$ yrs) stages in our fully unresolved models (variant I) are somewhat separated from these loci making these areas good places to look for young and old unresolved PN candidates. The change in color values for all CCDs is a result of the nebular emission as seen in Figure 5 and is described in the previous section.

4. Discussion

4.1. Comparison to Prior Work

A number of studies have shown that emission line objects like PNe can be identified or characterized by broadband filters. In the UV, optical, and near-infrared (NIR) Veyette et al. (2014) used Hubble Space Telescope (HST) Advanced Camera for Surveys and the Panchromatic Hubble Andromeda Treasury (PHAT) survey data to identify broadband detections of known PNe and were able to roughly estimate the excitation classifications of groups of PNe with the addition of archival m5007 narrow-band magnitudes. In the optical, Kniazev et al. (2014) used the Sloan Digital Sky Survey (SDSS) colors to separate PNe from other point sources in the outskirts of M31. Parker et al. (2016) developed a method using mid-infrared colors to identify galactic PN candidates using available data in the Galactic Legacy Infrared Mid-Plane Survey Extraordinaire I (GLIMPSE I) point source archive. Data from the Infrared Astronomical Satellite (IRAS) was used to inform optical spectroscopic follow-up by Suárez et al. (2006) to confirm PNe candidates along with post-AGB stars and sources transitioning from AGB to PN stages of evolution. Corradi et al. (2008) and following papers in the series (Corradi et al., 2010; Rodríguez-Flores et al., 2014) describe a similar method for selecting symbiotic star candidates using the INT Photometric $\rm H\alpha$ Survey (IPHAS) and the Two Micron All Sky Survey (2MASS) colors showing that symbiotic stars can be distinguished from PNe through colors but ultimately spectroscopic follow-up is necessary for verification.

We compare our predicted PNe colors and magnitudes (for variant I) to confirmed PNe in M31 (Kniazev et al., 2014). Kniazev et al. (2014) used 29 known PNe in M31 (Jacoby & Ford, 1986; Nolthenius & Ford, 1987) to construct magnitude and color criteria to identify additional PNe with SDSS photometry. 70 were confirmed to be PNe through follow-up spectroscopy (Kniazev et al., 2014). In Figure 7 the color-magnitude and color-color diagrams presented in Kniazev et al. (2014) are reproduced along with some of their selection criteria. We queried SDSS (DR14; Abolfathi et al., 2018) for clean sources within 1 degree of M31 and with point source criteria outlined in Kniazev et al. (2014). The spectroscopically verified PNe from Nolthenius & Ford (1987) along with our 3 $M_{\odot}$ PN evolutionary models scaled to the distance of M31 (760 kpc; van den Bergh, 1999) are overlaid onto the CMD and CCD. At such a distance, all of our nebular models are unresolved. In the CMD, the sample identified by Kniazev et al. (2014) is consistent with the youngest stages of our PN models ($t_{\rm age}<7300$ years). We find similar consistency in the CCD, except for youngest stage ($t_{\rm age}\sim 4200$ years), which is outside of the $\rm(g-r\geq-0.4)$ color cut used by Kniazev et al. (2014) to reduce contamination from high redshift objects. As suggested by our PN models, this color cut will limit the identification of the youngest and oldest PNe in M31, however, these PNe might be present and identifiable using the CCD.

We also compared our results for resolved and unresolved colors and magnitudes to six faint PNe (PN G094.0+27.4, PN G211.4+18.4, PN G049.3+88.1, PN G025.3+40.8, PN G047.0+42.4, and PN G158.8+37.1) recovered by Yuan & Liu (2013) through identification of excess $\rm[OIII]~{}\lambda 5007$ or $\rm[OIII]~{}\lambda 4959$ emission in the SDSS spectra of objects who’s photometry was affected by the nebula. For each PN, we used the object’s coordinates to obtain an SDSS g band image of the region (DR14; Abolfathi et al., 2018) and visually identified the PN in the image (using ds9). We then visually located the closest object to the CS position of the PN and obtained the SDSS ugriz magnitudes of that source. In Figure 8 we compare the colors of these six sources in relation to our models. The colors for PN G094.0+27.4, PN G211.4+18.4, PN G047.0+42.4, and PN G158.8+37.1 are consistent with that expected from a central star with very little nebular emission as they appear close to our model variants VIII-X (marked as stars in Figure 8). The nebular radii of these three objects are > $90^{\prime\prime}$, which is much larger than the typical $\sim 1^{\prime\prime}$ seeing resolution limit of SDSS, hence only a small fraction of the nebular emission is detected in these objects. PN G025.3+40.8 (marked as a triangle in Figure 8), with a nebular radius of $13^{\prime\prime}$, is consistent with a higher fraction of nebular flux detected as it appears near our variant VI models indicating a fractional nebular flux of $\sim 9\times 10^{-3}$ was observed with the CS. PN G049.3+88.1 (marked as a circle in Figure 8), with a nebular radius of $2\farcs 7$, is consistent with nearly all of the nebular flux detected as it appears near our variant I models. While this is not a rigorous comparison of these sources to our models, generally, we see that CSPNe with a larger sized nebula contain less nebular flux than CSPNe with a smaller sized nebula as predicted by our models.

4.2. Applications to Upcoming Surveys

Upcoming surveys such as LSST can potentially greatly expand the number of known PNe at different stages of evolution, both within the Milky Way and in other nearby galaxies. As an example, other authors have estimated LSST’s distance limit for RR Lyrae variables to be 400 kpc (LSST Science Collaboration et al., 2009); LSST’s potential reach for PNe may also be similarly impressive. LSST will commence operations in 2023 and is scheduled to be a 10-year long survey of the southern sky aiming at addressing questions from scales of the solar system to dark energy. The telescope harbors an 8.4 meter mirror and a 3,200 megapixel camera that will image 37 billion stars and galaxies. The camera is equipped with the $ugrizy$ filter set with sensitivity information detailed in Table 3.

The first consideration for LSST’s reach of PNe is spatial resolution of PNe flux. The optimal seeing limit for LSST is projected to be $\sim 0\farcs 7$ (LSST Science Collaboration et al., 2009), if we take this as an estimate of the size of aperture ($\Omega_{\rm LSST}$), then for any $\Omega_{\rm neb}>\Omega_{\rm LSST}$ the nebula will potentially be resolved. In Figure 9, we calculate the ratio $\Omega_{\rm LSST}/\Omega_{\rm neb}$ for a range of distances and include the locations of our synthetic nebular models shown in Figures 5 and 6. At a distance to the Magellanic Clouds (49.97 and 62.1 kpc; Pietrzyński et al., 2013; Graczyk et al., 2014) young models ($\rm t_{\rm age}<6300~{}yrs$) are unresolved while older models are marginally resolved. For our oldest models ($\rm t_{\rm age}>6300~{}yrs$) the nebular diameters are $\sim 1\farcs 0$. With worse seeing or by utilizing a slightly adjustable aperture in the imaging analysis, it is possible to observe all the nebular flux for marginally-resolved cases.

The reach for a given survey is also dependent on the survey saturation and limiting magnitudes (see Table 3; LSST Science Collaboration et al., 2009), as well as extinction to the source. Using Table 4, it is possible to study unresolved nebular fluxes for specific combinations of distance and extinction to determine the reach for LSST. Because the identification of potential PNe requires detection in at least three filters ($g$, $r$, and $i$) we considered the saturation and 5$\sigma$ limiting magnitudes, $m_{\rm sat}$ and $m_{5\sigma}$, respectively, for these three filters. In Figure 10, we determined the LSST detection limit for our unresolved nebular models as a function of distance and extinction. For modest extinction values ($A_{\rm V}<\rm 5~{}mag$) at the LMC/SMC, all of our nebular models predict detectable emission in $g$, $r$, and $i$, however, we potentially reach the saturation limit in $g$ and $r$, for our youngest models when extinction drops below 1 mag (see Figure 10). For another example, if LSST observed a system like Andromeda (M31) at a distance of 760 kpc (van den Bergh, 1999), in the areas of low extinction ($A_{\rm V}\lesssim 2{\rm~{}mag}$), similar to those found in the outskirts of that galaxy (Kniazev et al., 2014; Dalcanton et al., 2015), the youngest nebular models can be detected in these three bands, but fewer of the older nebular models would be detected in $i$ ($\rm t_{\rm age}>7300~{}yrs$). Hence, surveys like LSST can be used to discover many new potential PNe beyond our galaxy.

We can also consider LSST’s reach for PNe specifically within the Milky Way. Here, estimating the reach is further complicated by the numerous combinations of nebular age/size, distance, and extinction. The contribution of the nebular flux to the aperture is a strong function of $\Omega_{\rm LSST}/\Omega_{\rm neb}$ (see Figure 5). To consider the reach within the Milky Way, we provide variants of the absolute magnitudes in Table LABEL:table:_neb_fraction for a range of $\Omega_{\rm LSST}/\Omega_{\rm neb}$. For example, we considered the distance-$A_{\rm V}$ based reach for LSST at a distance consistent with the galactic center (8.3 kpc; Majaess et al., 2018). Young PNe ($t_{\rm age}\approx 3500{\rm~{}yr}$, $r_{\rm neb}\sim 0.035{\rm~{}pc}$) will be resolved at the galactic center with only $\sim$1.6% of the nebular emission measured by a $0\farcs 7$ aperture ($\Omega_{\rm LSST}/\Omega_{\rm neb}=0.075$). As a result, such a PN would be detectable in $g$, $r$, and $i$ for lines of sight where $A_{\rm V}<10{\rm~{}mags}$. However, for lines of sight where $A_{\rm V}$ drops below $2{\rm~{}mags}$, the model-predicted flux would reach the saturation limits. For older models ($\rm t_{\rm age}>6300~{}yrs$) at the distance of the galactic center, since the nebula grows in size, the nebular flux contribution reduces further. For these older models, saturation is unlikely with $A_{V}$ serving as the only limit on detection with maximum line of sight $A_{V}$ values in the range of $5-10{\rm~{}mags}$.

4.3. Limitations of the Current Work

The methodology and results presented in this work demonstrate the potential of using synthetic optical PNe spectra and a set of filters to produce accurate magnitudes and colors. However, we have only considered the evolutionary track for a $\rm 3~{}M_{\odot}$ progenitor with chemical abundances similar to those of galactic PNe. Constructing these models for a range of progenitor masses will give a better understanding of the distribution of PNe in color-color space and their identification in future surveys.

In developing our methodology, we discovered that a crucial ingredient to construct accurate models is the hydrodynamical expansion over time of the nebular radius as the central star evolves. In particular, self-consistent radiation-hydrodynamic models of a PN (e.g., Jacob et al., 2013; Schönberner et al., 2005; Toalá & Arthur, 2014) are important ingredients necessary to track to evolutionary behavior of the nebula and its radiation.

As mentioned in section 2.4, this method can only identify PNe candidates. All potential PNe candidates will require supporting spectroscopic or narrowband observations to verify their nature. This is also necessary to distinguish true PNe from other emission line objects such as Wolf-Rayet (WR) stars, HII regions, symbiotic stars, and CVs which may contaminate the same color-color space.

Finally, our study highlights an immediate outcome and discovery space made possible by a survey like LSST. We have not carefully considered the impact of repeated observations, which can be used to search for variability of the central star and/or potential companions and push the 5-$\sigma$ limiting magnitude fainter when observations are combined.

5. Conclusions

In this work we presented synthetic absolute magnitudes in the $ugrizy$ filters of 13 PN models representing the evolution of a PN for a $\rm 3~{}M_{\odot}$ progenitor star over $\sim$5,000 years (Table 4). We have calculated the colors for various photometric aperture sizes to explore spatially resolved and unresolved cases. We showed that our model magnitudes and colors are consistent with real observations for both resolved and unresolved cases. We also showed that LSST will allow for the identification of many PNe in the Milky Way and neighboring galaxies.

Color-color diagrams will be useful in differentiating PNe from other point sources in the upcoming LSST era. We showed that the PNe model colors are entirely separate from most other SDSS objects cataloged in SDSS; indeed, we can simply provide a cutoff above $r-i=-0.75$ to find viable PNe candidates. Moreover, we showed that colors can also differentiate younger PNe ($\rm t_{\rm age}<7000~{}yrs$) from older PNe ($\rm t_{\rm age}>7500~{}yrs$). This will likely still be possible even in the presence of large amounts of extinction. There is still work to be done on how colors of PNe with various progenitor masses will change.

Given the LSST seeing limited resolution and magnitude limits, PNe detected at various distances will be spatially resolved to varying degrees. For distances to the galactic center, the youngest PNe ($\rm t_{\rm age}\approx 3500~{}yrs$) will saturate with $\rm A_{\rm V}\leq 2~{}mags$, be visible up to $\rm A_{\rm V}\approx 10~{}mags$, and be partially resolved. The oldest PNe ($\rm t_{\rm age}\approx 8300~{}yrs$) will be visible up to $\rm A_{\rm V}\approx 5~{}mags$ and be mostly resolved. Young PNe ($\rm t_{\rm age}<6300~{}yrs$) at the LMC/SMC will be fully unresolved and visible for moderate $\rm A_{\rm V}$ between $\rm\sim 1~{}and~{}5~{}mags$ while older models will be resolved containing at least 20% of the nebular emission. For distances similar to M31, all PNe will be unresolved but only younger PNe ($\rm t_{\rm age}<7300~{}yrs$) will be visible with $\rm A_{\rm V}$ between $\rm\sim 0.6~{}and~{}2.4~{}mags$. This and future works will help prepare the astronomical community for the massive amounts of data that will be provided by LSST and enable the discovery of unprecedented numbers of new PNe, enabling detailed tests of our current theories of this phase of stellar evolution.

Acknowledgements

The specific use of Cloudy in this work is through the Python wrapper pyCloudy which allows the user to interact with the Cloudy input and output files through Python (Morisset, 2013). This research has made use of the Spanish Virtual Observatory (http://svo.cab.inta-csic.es) supported from the Spanish MINECO/FEDER through grant AyA2014-55216. This research has made use of the VizieR catalogue access tool, CDS, Strasbourg, France. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This work made use of the IPython package (Pérez & Granger, 2007). This research made use of SciPy (Jones et al., 2001). This research made use of NumPy (Van Der Walt et al., 2011) This research made use of matplotlib, a Python library for publication quality graphics (Hunter et al., 2007). This research made use of ds9, a tool for data visualization supported by the Chandra X-ray Science Center (CXC) and the High Energy Astrophysics Science Archive Center (HEASARC) with support from the JWST Mission office at the Space Telescope Science Institute for 3D visualization. This research made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration et al., 2013). Funding for the Sloan Digital Sky Survey IV has been provided by the Alfred P. Sloan Foundation, the U.S. Department of Energy Office of Science, and the Participating Institutions. SDSS-IV acknowledges support and resources from the Center for High-Performance Computing at the University of Utah. The SDSS web site is www.sdss.org. SDSS-IV is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS Collaboration including the Brazilian Participation Group, the Carnegie Institution for Science, Carnegie Mellon University, the Chilean Participation Group, the French Participation Group, Harvard-Smithsonian Center for Astrophysics, Instituto de Astrofísica de Canarias, The Johns Hopkins University, Kavli Institute for the Physics and Mathematics of the Universe (IPMU) / University of Tokyo, Lawrence Berkeley National Laboratory, Leibniz Institut für Astrophysik Potsdam (AIP), Max-Planck-Institut für Astronomie (MPIA Heidelberg), Max-Planck-Institut für Astrophysik (MPA Garching), Max-Planck-Institut für Extraterrestrische Physik (MPE), National Astronomical Observatories of China, New Mexico State University, New York University, University of Notre Dame, Observatário Nacional / MCTI, The Ohio State University, Pennsylvania State University, Shanghai Astronomical Observatory, United Kingdom Participation Group, Universidad Nacional Autónoma de México, University of Arizona, University of Colorado Boulder, University of Oxford, University of Portsmouth, University of Utah, University of Virginia, University of Washington, University of Wisconsin, Vanderbilt University, and Yale University. Funding for the SDSS and SDSS-II has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, the U.S. Department of Energy, the National Aeronautics and Space Administration, the Japanese Monbukagakusho, the Max Planck Society, and the Higher Education Funding Council for England. The SDSS is managed by the Astrophysical Research Consortium for the Participating Institutions. The Participating Institutions are the American Museum of Natural History, Astrophysical Institute Potsdam, University of Basel, University of Cambridge, Case Western Reserve University, University of Chicago, Drexel University, Fermilab, the Institute for Advanced Study, the Japan Participation Group, Johns Hopkins University, the Joint Institute for Nuclear Astrophysics, the Kavli Institute for Particle Astrophysics and Cosmology, the Korean Scientist Group, the Chinese Academy of Sciences (LAMOST), Los Alamos National Laboratory, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, Ohio State University, University of Pittsburgh, University of Portsmouth, Princeton University, the United States Naval Observatory, and the University of Washington.

George Vejar thanks the LSSTC Data Science Fellowship Program, his time as a Fellow has benefited this work.