A SPITZER IRS SURVEY OF NGC 1333: INSIGHTS INTO DISK EVOLUTION FROM A VERY YOUNG CLUSTER

As part of a Spitzer GTO program (PID 40525) we observed 78 objects in NGC 1333 on 2007 October 5–6 and 9–10, 2008 February 21 to March 2, and 2008 March 30, in Infrared Spectrograph (IRS; Houck et al. 2004) campaigns 44, 48, and 49, respectively. Our sample of YSOs was chosen based on their excess in photometric observations from both the Spitzer Infrared Array Camera (IRAC; Fazio et al. 2004) 8.0 μm band and Multi-band Imaging Photometer (MIPS; Rieke et al. 2004) 24.0 μm band, published in Gutermuth et al. (2008). We observed all objects with detections (σ < 0.2 mag) in the IRAC 8.0 μm band and MIPS 24 μm band that are brighter than 11.54 mag at 8.0 μm and 8.5 mag at 24 μm. For this reason we adopt their identifiers for our sample (these identifiers were also used in Winston et al. 2009, 2010). The basic data for this sample, namely, [GMM08] SED index, R.A., decl., AOR ID, and alternate names, are listed in Table 1.

This selection excluded a total of 55 objects from the Gutermuth et al. (2008) YSO sample: 39 had no MIPS 24 μm photometry, 2 had no IRAC 8.0 μm photometry, and 14 were too faint at either 8.0 μm or 24 μm. Of the Gutermuth et al. (2008) sample, four objects were already scheduled to be observed with IRS, and therefore these objects were not observed in program 40525. One of these conflicts was the Class 0 protostar IRAS 2A, which was observed in Spitzer program 2. The other three, ASR 118, SSTc2d J032924.1+311958, and SSTc2d J032929.3+311835, were observed in program 30843 on 2006 September 18 and 2007 March 9 and 21. These objects were selected from their photometry as possible "cold disks" that have inner holes and are presented in Merín et al. (2010), indexed there by #2, #4, and #5, respectively. Gutermuth et al. (2008) identify these objects as #67, #122, and #137, and we include these three observations with IRS in our sample using the Gutermuth et al. (2008) identifiers. This addition brings the total number of observations to 81.

Selecting for objects with excess in the 8.0 μm IRAC and 24 μm MIPS bands biases our sample toward earlier star formation stages (Classes 0/I and II), where massive dusty envelopes and disks are present. Gutermuth et al. (2008) find the 137 cluster members identified with IRAC and MIPS to represent 83% ± 11% of the total cluster membership. Our sample consists of over half of these sources (N = 81), which is ∼49% of the extrapolated total cluster membership.

We observed each object from 5 to 38 μm using the Staring Mode of the Spitzer IRS in both orders of the short–low (SL, 5.2–14.5 μm) and long–low (LL, 14–38 μm) resolution modules (λ/Δλ ∼ 60–120). From the Spitzer Science Center's (SSC) basic-calibrated data products using pipeline S18.7, we extracted a spectrum using the IRS instrument team's Spectral Modeling, Analysis and Reduction Tool (SMART; Higdon et al. 2004). We corrected for bad or "rogue" pixels by interpolating the flux from two adjacent good pixels in the dispersion direction.

For each observation, we took two images with the object at different positions, or nods, in the slit, each one-third of the width of the slit away from the edge. For most images in our sample, we followed the standard procedure for extracting spectra as is done in Furlan et al. (2006, 2008), Kim et al. (2009), McClure et al. (2010). We removed background emission by subtracting the image with the source at one nod from the image with the source at the opposite nod. We then extracted the spectrum from the sky-subtracted two-dimensional (2D) image using a column that tapers as the width of the point-spread function does. Spectral extraction and removal of background emission for 60 of these objects proved troublesome because of a combination of bright cloud emission, multiple sources in the slit, and/or the extinguished, faint nature of the source. For these observations with particularly difficult background emission, we individually inspected and carefully extracted a spectrum using the SMART optimal extraction tool, AdOpt (Lebouteiller et al. 2010). AdOpt allows the user to narrow the extraction window and fit background emission in that window with an nth order polynomial. Using AdOpt, we fit the background emission in most cases with a polynomial of zeroth or first order, though a few cases required as high as third order.

With AdOpt, tapered column extraction, or a combination of the two methods, we extracted a spectrum for each object from each module, order, and nod. To calibrate the flux of the spectra, we created relative spectral response functions (RSRFs) by dividing a template spectrum, α Lac (A1V; Cohen et al. 2003) for SL (5.2–14.5 μm) and ξDra (K2III) for LL (14.0–36 μm), by its tapered column or AdOpt-extracted IRS spectrum. We multiplied the target spectra by the similarly prepared RSRF for each module, order, and nod to calibrate the flux from the IRS. For each nod, we combined the calibrated spectra from each module and order. To produce a final spectrum, we averaged these spectra from the two nods and calculated uncertainty as half of the nod difference. Occasionally, slight mispointing resulted in a decrease in flux in the modules that were mispointed, causing a mismatch between the SL and LL spectra. We corrected for this by scaling the fainter module, which presumably lacks flux due to mispointing, up by a scalar factor. A typical scaling factor was 1.1, with a maximum value of 3.5 (#34). As an additional check on our reduction, we extracted a spectrum for each object with the optimal extraction software, UR Optimal (W. Forrest and C. Tayrien 2011, private communication), which extracts a spectrum using data from 5, 7, or 9 pixels around the nod position. The spectrum extracted with tapered column extraction and/or AdOpt had smaller nod-difference uncertainties than that of UR Optimal; therefore, we use the tapered column and/or AdOpt-extracted spectrum for our analysis. We estimate the absolute photometric uncertainty in flux in the final spectra to be ∼5%.

We supplement the IRS spectra with broadband photometry from the literature in the J, H, and K_s bands of the Two Micron All Sky Survey (2MASS; Skrutskie et al. 2006). NGC 1333 was also observed as part of the 2MASS extended mission, 2MASS 6X. Data from this mission were taken from observations with six times the integration time of standard 2MASS observations. This resulted in ∼1 mag improvement in sensitivity limit and reduced the uncertainty for faint objects near the detection limit of regular 2MASS data (Carpenter et al. 2001). Where available, we use the 2MASS 6X data over the standard 2MASS from the literature. We also supplement our data, where available with R and I broadband photometry from Aspin et al. (1994). We used 2MASS K_s broadband photometry in our SEDs, and where 2MASS K_s-band photometry was not available we used K-band photometry from Aspin et al. (1994) or Lada et al. (1996). Every object was observed photometrically with Spitzer, and we supplement our data with broadband photometry from Gutermuth et al. (2008) in the IRAC 3.6, 4.5, 5.8, and 8.0 μm bands and the MIPS 24 μm band.

2.2.1. Notes on Individual Objects

Spatially variable, bright sky emission resulted in careful reduction of 60 of the 81 images using off-nod and off-order sky subtraction, as well as polynomial fit sky subtraction with AdOpt. Below are the details of particularly troublesome objects. We omit objects from our sample that have a spectrum that is not representative of a YSO, namely, #8 and #38.

#5. Bright emission lines in the SL spectrum could not be removed through various sky subtraction methods; this results in a seemingly discontinuous spectrum.

#8. Identified as IRAS 4B in Gutermuth et al. (2008), #8 is not identical with the Class 0 object but corresponds instead to one of IRAS 4B's outflow lobes. The observation for #8 was pointed at the bright emission originating from one of IRAS 4B's outflow cavities to the south of the embedded object's submillimeter coordinates. The SED is shown in Figure 1, and it is left out of our sample. The SL spectrum is dominated by the emission from the outflow as evidenced by the collisionally excited molecular hydrogen emission lines. The LL spectrum contains a forest of attenuated, blended water emission lines from IRAS 4B itself (Watson et al. 2007) though it lies about 5'' away, at the edge of the slit.

Zoom In Zoom Out Reset image size

#23. Careful extraction was unable to remove H₂ and polycyclic aromatic hydrocarbon (PAH) emission features, which do not originate from #23.

#34. The SL spectrum was very difficult to disentangle from the sky emission and had to be scaled up by 3.5 to match the LL spectrum. Also, the SL spectrum possibly suffers from an offset in the dispersion direction.

#37. The LL spectrum is dominated by a bright source near the target coordinates. The flux contribution of this nearby bright source dominated the spectrum and could not be removed. For this reason we only present the SL spectrum of #37 and only use the SL spectrum in our analysis.

#38. Based on the presence of redshifted PAH emission, #38 is a galaxy. Using the redshifted 6.2, 7.7, 8.6, and 11.2 μm PAH emission features, we calculate a redshift from each feature and obtain an average redshift of 0.29. The SED of the galaxy labeled here as #38 is shown in Figure 2, and it is left out of our sample.

Zoom In Zoom Out Reset image size

#52. The second module of SL nod 1 had data missing shortward of 5.43 μm as a result of the SSC pipeline trying to correct for rows in the detector that were altered by peak-up camera saturation. This resulted in an unreliable nod 1 spectrum, and for this reason we use only the nod 2 spectrum for SL2.

#58. Two close sources in the LL slit provided some difficultly. Narrow window extraction at nod position with AdOpt yielded the best result. LL was scaled by 0.7 to match SL. We include #58 in our analysis with the caveat that the LL spectrum is likely not solely representative of #58.

#63. Both nods of the first module of SL displayed a striping effect in the dispersion direction, deemed "jail-bars," likely due to correction for a bright source in the peak-up array. Careful extraction using AdOpt was able to recover the source emission, and we use this spectrum in our analysis. This effect, as well as the faint nature of the object, is the likely cause of the noise in the spectrum of #63.

#67. Both nods of the first module of LL were severely affected by stray emission in the detector from a very bright source just outside the LL slit; for this reason we exclude the LL1 data from the SED. This effect is exemplified in its spectrum in Merín et al. (2010) as a wide false 30 μm feature, also seen in the IRS spectrum of Hn 10 E in Chamaeleon I (Manoj et al. 2011). Merín et al. (2010) classify this source (called therein ASR 118) as a "cold disk," a term equivalent to our TD. However, due to lack of reliable data in LL1, we exclude #67 from our analysis and cannot definitively address its radial structure.

#68. The emission lines in the SL module at 5.5, 6.1, 6.9, 8.0, and 9.6 μm are highly spatially variable molecular hydrogen lines from the sky, which could not be removed in the reduction process.

#75, #94, #110, and #133. The spectrum of each of these objects contains one or more spikes due to difficult-to-remove background emission and/or possible remnant bad pixels.

#84. Sky subtraction in LL proved difficult; the strong dip around 30 μm is not a real feature but could not be removed.

#92 and #102. Variable, bright sky emission and the faint nature of the source cause these carefully extracted spectra to be noisy and uncertain.

Classification of a YSO is typically accomplished using the slope of its SED spanning from the near- to mid-infrared (Lada 1987). We find the logarithmic slope, n, or spectral index, in a given wavelength range, henceforth expressed as $n_{\lambda _{1}-\lambda _{2}}$, which is given as

$$\begin{equation} n_{\lambda _1 - \lambda _2} = \frac{\log (\lambda _2 F_{\lambda _2})- \log (\lambda _1 F_{\lambda _1})}{\log (\lambda _2) - \log (\lambda _1)}. \end{equation} \tag{ 2 }$$

An IR slope, typically from 2 to 25 μm, greater than 0.3 corresponds to that of a Class I object, n_{2 − 25} between −0.3 and 0.3 is that of a Flat Spectrum object, n_{2 − 25} between −1.6 and −0.3 defines a Class II object, and n_{2 − 25} less than −1.6 classifies a Class III object. The empirical classification of the SED of YSOs has become inextricably tied with an evolutionary sequence (Adams et al. 1987) as discussed in Section 1; therefore, proper classification is important in order to understand the evolutionary state of these systems.

We empirically classify the observed SEDs of the objects in our sample with 2MASS K_s-band photometry using n_{2 − 25} and report these values and the associated classification in Table 2. The flux at 2 μm is the 2MASS K_s photometric flux (λ = 2.17 μm), and the flux at 25 μm is the integrated flux of the observed IRS spectrum from 23.5 to 26.5 μm. SED classification from 2 to 25 μm could not be done for eight objects because they lack detection in the K_s band of 2MASS and for one object (#67) because of unreliable data near 25 μm. We find 16 Class I objects, 14 Flat Spectrum (henceforth FS) objects, and 40 Class II objects based on n_{2 − 25} SED classification. No objects were classified as Class III based on their n_{2 − 25}. This is not surprising as Class III objects would not be bright in the mid-infrared and we selected our sample based on their flux at 8 and 24 μm.

When interpreted as an indication of evolutionary state, the n_{2 − 25} index is at times known to be misleading (Robitaille et al. 2006; Crapsi et al. 2008; McClure et al. 2010; Manoj et al. 2011). Edge-on disks are known to have very red spectral indices (Pontoppidan et al. 2005) and are frequently identified as Class I objects based on n_{2 − 25}. The Spitzer IRS spectra reveal silicate emission and ice absorption features, which are not accounted for in an empirical n_{2 − 25} classification and allow us to gain a better understanding of the evolutionary state of these objects. Also important to note is that K_s-band measurements from the extended 2MASS mission were taken in 2000 and the IRS spectra in 2007 and 2008. As these T Tauri stars are known to be highly variable, this difference in epochs of observation can affect the 2–25 μm slope and subsequent SED classification. The empirical SED classification from 2 to 25 μm gives us one look into the evolutionary state of these objects.

To better understand the evolutionary state of our sample, including those objects without known 2MASS K_s-band measurements, we examine the "extinction-free" spectral indices of the observed spectra for the entire sample (McClure et al. 2010). These indices, n_{5 − 12} and n_{12 − 20}, characterize the SED slope of the spectra before extinction correction between wavelengths with roughly the same value of A_λ/A_V according to McClure (2009), namely, 5.3, 12.9, and 19.8 μm. The flux at each of these wavelengths is calculated by averaging the flux over three channels, i.e., 0.3 μm, centered at the anchor point. This index is truly extinction free for objects with A_V > 8 and varies slightly with wavelength for objects with A_V < 8. Based on the n_{5 − 12} spectral index, we classify the evolutionary state of the sample (N = 79) according to the following criteria: n_{5 − 12} ⩾ −0.2—the spectrum is dominated by envelope emission; −0.2 > n_{5 − 12} ⩾ −2.25—dominated by disk emission; and n_{5 − 12} < −2.25—dominated by photosphere emission (McClure et al. 2010). Using n_{5 − 12}, 23 objects are classified as envelope dominated, 55 as disk dominated, and 1 as photosphere dominated. Large values of n_{12 − 20}, greater than 2.0, indicate typical values for TDs in Taurus and Ophiuchus (McClure et al. 2010).

Figure 6 plots these two extinction-free indices against each other, allowing us to gain another perspective of the evolutionary state of the objects in our sample. Three objects (#7, #35, and #137) fall within the area of this plot where TDs are found as they are characterized by a steeply rising SED beyond 12 μm (McClure et al. 2010; Furlan et al. 2011). Upon examination of the SEDs of #7 and #35, the presence of deep silicate and CO₂ ice absorption features in their spectrum and absence of detection in the three bands of 2MASS indicate that these objects are deeply embedded in an envelope. Their decreasing SED from 5 to 12 μm and steeply rising SED longer than ∼10 μm suggest that these are deeply embedded protostars (Furlan et al. 2008) and are not TDs. The extinction-free indices are misleading here as the steep excess at longer wavelengths is a result of an envelope rather than, as would be the case for TDs, the excess originating at the outer disk wall (Espaillat 2009). Furlan et al. (2011) found that embedded protostars seen edge-on are also found in this area of the plot. However, examination of the SED of #137 shows no evidence of a significant envelope affecting the spectrum; therefore, the steep rise picked out in a high n_{12 − 20} is likely from emission originating in the outer layers of a heated outer disk wall. We will return to #137 in our discussion of radial structure in Section 5.2.

Zoom In Zoom Out Reset image size

We examine the SEDs of our sample with respect to their n_{5 − 12} classification and identify objects that require more careful consideration and possible reclassification. Figure 7 displays the SED classes (based on n_{2 − 25}) with respect to the extinction-free index n_{5 − 12}. Of the 23 envelope-dominated objects, 12 are Class I, 3 are FS, and 4 are Class II. The other four lack data either at 25 μm, as in the case of object #67, or in the K_s band and therefore could not be classified. Ice and silicate absorption features, as a result of local extinction, are evidence for a significant envelope affecting the SED (Boogert et al. 2008; Furlan et al. 2008; Zasowski et al. 2009). We use the presence or absence of these features to get an idea of how strongly an envelope dominates the emission of these objects, as well as whether their SED might truly be considered as dominated by disk emission. The four Class II objects that are envelope dominated according to n_{5 − 12} (#68, #69, #82, and #114) plus object #67, which has no n_{2 − 25} but is envelope dominated based on n_{5 − 12}, do not show ice or silicate absorption features in their observed or dereddened spectra; for these reasons we reclassify these objects as disk dominated though we do not include them in our disk analysis in Section 5. If an envelope were present, we would derive a large A_V for these objects in Section 3. The relatively low extinction values calculated for #68, #69, #82, and #114 (A_V = 5.4, 3.7, 5.7, and 3.0, respectively) signal that an envelope does not strongly affect these systems. These objects are also close to the disk/envelope boundary in n_{5 − 12}, especially #67, #68, #69, and #82, which have values between −0.20 and −0.09. This "boundary" is tentative and an approximation based on the properties and locations of objects in Taurus and Ophiuchus in this plot. The uncertainties of #68, #69, and #82 span the disk/envelope boundary. Five Class I envelope-dominated objects (#16, #19, #20, #113, and #115) display significant 10 μm silicate emission features in their dereddened spectrum. Silicate emission features can be seen from the dusty disk in a Class II object, from the warm envelope in a face-on Class I object (Watson et al. 2004), and from both a disk and a warm residual envelope in an FS object. These objects appear more evolved than the canonical Class I stage, and though they were not classified observationally by n_{2 − 25} as FS, these objects resemble the historical definition of an FS object as one that is in a transition phase between Class I and Class II, whose envelope is present yet tenuous and disk is beginning to be seen in the SED (Casali & Matthews 1992; Greene & Lada 1996). This explanation is consistent with the envelope-dominated classification using n_{5 − 12}.

Zoom In Zoom Out Reset image size

Worthy of note in the envelope-dominated sample is object #122, or SSTc2d J032924.1+311958, classified in Merín et al. (2010) as a likely edge-on disk due to ice and silicate absorption features in the IRS spectrum. Our analysis classifies #122 as a Class I envelope. As stated in Section 4.1, the n_{2 − 25} classification of edge-on disks is known to be misleading. However, the extinction-free spectral index n_{5 − 12} classifies six of seven edge-on disks in the Taurus sample of Furlan et al. (2011) as disks, with the one outlier being HARO 6-5B with an n_{5 − 12} = 0.26. Object #122 has an n_{5 − 12} of 0.69, as well as ice and silicate absorption features that suggest that a significant envelope is still present in this system. Further, the spectrum also resembles that of the famous edge-on disk, HH 30 in Taurus (Furlan et al. 2011).

Of the 55 disks, 5 are Class I, 11 are FS, and 36 are Class II according to n_{2 − 25}. Three disk-dominated objects (#3, #7, and #35) have no SED class as they were not detected in the K_s band of 2MASS; as mentioned before, this is a strong sign that these objects are deeply embedded protostars. The SEDs of #3, #7, and #35 display strong 10 μm silicate and 15 μm CO₂ ice absorption features, as well as a steeply rising SED beyond ∼10 μm. For these reasons we reclassify #3, #7, and #35 as envelope dominated. Four more objects classified as dominated by disk emission require closer examination (#13, #33, #36, and #63). Class I based on n_{2 − 25}, #63 has strong ice absorption features, a steep drop in flux from 8 to 9 μm, and a flat SED >12 μm that is reminiscent of the most deeply embedded (A_V > 30) Class II objects in Ophiuchus (McClure et al. 2010). The other objects, #13, #33, and #36, are Class I and show strong silicate absorption and CO₂ ice absorption features. These objects also seem to suffer from high extinction and demonstrate how the extinction-free indices can give you a more accurate idea of evolutionary state, especially for objects with large A_V. Based on a closer examination of the n_{5 − 12} classification of objects in our sample with respect to their SEDs, our YSO sample in NGC 1333 consists of 21 objects dominated by envelope emission, 57 by disk emission, and 1 by photospheric emission; the SEDs of these objects are shown in Figures 8–10, respectively.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

According to Aspin et al. (1994), the photosphere-dominated object #135 (ASR 54) has K = 15.74, for which we would assign it to Class I. The SED of ASR 54 (Figure 10) shows 15 μm CO₂ ice absorption, as well as ice absorption around 10 μm, likely a result of foreground absorption and not the presence of an envelope. ASR 54 was noted by Gutermuth et al. (2008) as a candidate for a star that appears to be exciting PAH emission near the IRAS 4 filament in the southern NGC 1333. From IRAC photometry and the K-band detection by Aspin et al. (1994), Gutermuth et al. (2008) suggest that this is a heavily extinguished (A_K ∼ 10) late B-type star with some excess at 24 μm, implying that it might be a transitional or debris disk. The extinction-free spectral index n_{5 − 12} agrees with Gutermuth et al. (2008) categorizing ASR 54 as a photosphere. There is certainly an excess at longer wavelengths, suggestive of a circumstellar accretion disk, though we lack reliable extinction information that would allow us to determine whether this object is a TD or debris disk. We also cannot conclusively say whether this is the object exciting the PAH emission seen in the southern end of NGC 1333.

From our sample of full IRS spectra and 2MASS K_s-band photometry, the empirical n_{2 − 25} Class I fraction is 23% ± 12% (16/70) and the fraction of sources dominated by envelope emission, as defined in the n_{5 − 12} analysis in Section 4.2, is 27% ± 11% (21/79). Gutermuth et al. (2008) found that 39 of the 137 members of NGC 1333 were protostars resulting in a protostar class fraction of 28%. Our findings are consistent with Gutermuth et al. (2008) and reaffirm that the protostar fraction is quite high in NGC 1333, which implies the presence of a substantial young population.

In comparison to a cluster of similar age analyzed in the same way, such as Ophiuchus (∼0.8 Myr, envelope fraction 9% (10/106), Class I fraction 25%; McClure et al. 2010), NGC 1333 has a significantly larger envelope fraction and a similar Class I fraction. Using Mathis's extinction curve with R_V = 5, i.e., A_V/A_J = 3.06, the median A_V for Ophiuchus is 10.7 mag (McClure et al. 2010) and the median A_V for NGC 1333 is 5.8 mag (see Section 3). While Ophiuchus is heavily extinguished, which skews age estimates from class fractions to the younger side, as extinction can make a disk have the SED class of an envelope (Pontoppidan et al. 2005), NGC 1333 is not as affected by extinction and therefore its large envelope fraction demonstrates that NGC 1333 is genuinely young. Given this, NGC 1333's larger envelope fraction would imply that it is even younger than Ophiuchus, whose age is estimated to be ∼0.8 Myr (McClure et al. 2010).

The fraction of our sample with disk-dominated emission, as defined in the n_{5 − 12} analysis in Section 4.2, is 72% ± 11% (57/79). Our sample of 70 objects with an SED class yields an FS fraction of 20% ± 12% (14/70) and a Class II fraction of 57% ± 12% (40/70). Gutermuth et al. (2008) identify 98 Class II objects in their 137 YSO sample, yielding a Class II fraction of 72%. Our results are consistent with the findings of Gutermuth et al. (2008) as their Class II sample is comparable to our disk-dominated sample rather than our Class II sample as they do not classify objects as FS. Comparison with Ophiuchus, which has a larger fraction of objects with disk-dominated emission (85%, 90/106), further demonstrates that NGC 1333 is likely to be younger than Ophiuchus.

Our selection for objects with excess in the 8 μm IRAC and 24 μm MIPS bands does little to bias these class fractions as disk- and envelope-dominated objects typically have substantial IR excesses at 8 and 24 μm. In the study of NGC 1333 X-ray sources by Winston et al. (2010), they found 41 Class III objects in their 95-object sample (43%). Winston et al. (2010) selected for likelihood of X-ray activity, and therefore SED class fractions would be biased toward later stages of star formation that have higher X-ray activity. Gutermuth et al. (2008) calculated a disk fraction of NGC 1333 using sources within a circle of radius of 55 or 0.4 pc centered at the median right ascension and declination of the YSOs they identified. Those authors calculate a disk fraction of 83% ± 11% as the fraction of the sample that has K_s < 14 to all objects detected at K_s after accounting for field star contamination. This method helps to mitigate biases against diskless members, which would make up the remainder of the cluster, specifically 17% ± 11%. In a study of low-mass stars and substellar objects, Wilking et al. (2004) estimate a disk fraction for the northern cluster of NGC 1333 to be 75% ± 20%, implying a diskless population of 25% ± 20%. They also use K-band excess for this calculation and correct for their incomplete census of members with disks.

When a protostar forms, it is surrounded by an infalling cloud of gas and dust that eventually lands as a flared accretion disk around the star. As the disk evolves, small dust grains coagulate, become larger, and eventually settle out of the upper layers of the disk down to the midplane. This decreases the amount of dust in the upper layers of the disk and reduces the heating and associated flaring, which causes the disk to become flatter (D'Alessio et al. 2006; Dullemond et al. 2007).

The continuum emission from circumstellar accretion disks can be characterized through the use of spectral indices that sample silicate-feature-free emission, such as n_{6 − 13} and n_{13 − 31} (D'Alessio et al. 2006; Furlan et al. 2006, 2009; Watson et al. 2009). The majority of the emission characterized by n_{6 − 13} and n_{13 − 31} is shown by D'Alessio et al. (2006) to originate in two distinct parts of the disk. Emission from the inner parts of a radially continuous disk is characterized by n_{6 − 13}, which generally corresponds to 0.1–3 AU for the least settled disks and 0.1–1 AU for the most settled disks. The middle to outer regions of radially continuous disks can be characterized by n_{13 − 31}, around 3–25 AU in the least settled disks and 1–10 AU for the most settled disks (D'Alessio et al. 2006).

We calculate n_{6 − 13} and n_{13 − 31} for our sample of disks in NGC 1333, as described above, using Equation (2) with fluxes at 6, 13, and 31 μm derived by integrating the extinction-corrected IRS spectra from 5.4–6.0 μm, 12.8–14.0 μm, and 30.3–32.0 μm, respectively. We report these values and their uncertainties in Table 3. The uncertainties in these spectral indices are propagated from the uncertainty in the original spectra. Some objects in the sample have significant uncertainty in their spectrum redward of 30 μm, resulting in a large uncertainty in n_{13 − 31}. The increase in noise seen in some spectra at wavelengths longer than 30 μm can be attributed to the data being taken late in the Spitzer mission (campaigns 44, 48, and 49) after the detector bias was decreased in the LL array. The degree of settling in the disk can be inferred through comparison of these continuum indices with disk models that take into account settling.

Table 3. Continuum Indices

Index	n2 − 6a	$\sigma _{n_{2-6}}$	n6 − 13	$\sigma _{n_{6-13}}$	n13 − 31	$\sigma _{n_{13-31}}$
18	−0.72	0.04	−0.31	0.05	0.08	0.04
44	−2.23	0.12	−1.39	0.14	−1.32	0.44
46	−1.91	0.12	−0.91	0.18	−1.16	0.17
48	−1.66	0.08	−1.18	0.14	−1.00	0.23
49	−1.16	0.13	−0.80	0.15	−0.78	0.24
50	−1.63	0.06	−0.62	0.07	−0.91	0.04
51	−1.91	0.09	−0.94	0.14	−1.08	0.62
52	−1.76	0.11	−0.47	0.14	−0.89	0.18
53	−1.07	0.13	−0.84	0.15	−1.93	0.64
54	−1.93	0.08	−1.07	0.12	−1.12	0.34
55	−1.88	0.09	−1.05	0.18	0.63	0.46
57	−1.89	0.11	−0.50	0.13	−0.49	0.12
58	−1.69	0.14	−0.52	0.15	−0.07	0.26
59	−1.65	0.06	−0.58	0.07	−0.39	0.09
61	−0.90	0.05	−0.59	0.06	−0.69	0.05
64	−1.35	0.19	−0.17	0.21	−0.34	0.07
65	−2.22	0.10	−0.33	0.12	−0.04	0.06
73	−1.32	0.08	−0.21	0.08	−0.20	0.06
74	−1.27	0.06	−0.97	0.07	−0.63	0.09
75	−1.68	0.19	−0.81	0.50	0.40	0.48
78	−0.65	0.09	−0.61	0.12	−0.45	0.18
84	−0.82	0.07	−0.24	0.08	0.38	0.20
88	−1.55	0.07	−0.89	0.09	−1.29	0.07
94	−1.19	0.16	−1.19	0.67	−0.64	0.74
99	−0.72	0.13	−0.58	0.17	−0.88	0.34
101	−1.67	0.17	−1.04	0.36	1.08	0.35
106	−1.32	0.06	−1.46	0.08	−0.43	0.07
110	−1.64	0.23	−0.94	0.59	2.35	0.55
111	−0.86	0.06	−0.90	0.05	−0.91	0.06
116	−0.95	0.07	−0.39	0.08	−0.55	0.08
118	−1.88	0.10	−0.83	0.15	−0.43	0.26
119	−2.01	0.13	−1.25	0.16	−0.67	0.13
120	−1.56	0.16	−0.62	0.19	0.02	0.14
121	−1.60	0.11	−0.85	0.13	−0.77	0.29
123	−1.07	0.06	−0.84	0.09	−2.19	0.15
125	−1.80	0.07	−0.34	0.09	−0.92	0.17
126	−1.76	0.12	−0.70	0.18	0.17	0.19
127	−1.39	0.13	−1.30	0.15	−0.35	0.15
128	−1.71	0.08	−1.34	0.16	−0.18	0.18
131	−1.74	0.09	−0.76	0.14	0.78	0.14
133	−1.29	1.01	−0.71	1.10	0.26	0.08
134	−1.60	0.11	−1.82	0.12	0.38	0.42
136	−2.68	0.08	−1.73	0.15	−0.34	0.31
137	−2.60	0.07	−0.22	0.12	1.79	0.09

Notes. Continuum indices for disk-dominated objects in NGC 1333. ^aObjects with upper limit for K are indicated with ellipses.

Download table as:  ASCIITypeset image

We compare our sample with predictions from models of irradiated accretion disks following the formalism of D'Alessio et al. (1998, 1999, 2001, 2005, 2006). In the models used in this work, the surface density (and hence mass) of the disk is proportional to the mass accretion rate, $\dot{M}$, over the Shakura–Sunyaev viscosity parameter, α. We chose to fix α at a typically adopted value of 0.01 and vary the mass accretion rate. The reverse could have also been performed. Therefore, varying α is equivalent to varying the mass accretion rate. These models probe the effects of turbulence on the heating of the disk and thus the shape of the SED. Turbulent mixing has a non-negligible effect over the vertical structure of the disk (i.e., Zsom et al. 2011; Ciesla 2010). However, studying the role of the turbulence as a mixing mechanism and its effects on the emission of the disk is beyond the scope of this paper. Following these models, a grid of disks around stars with M_* = 0.2 and 0.5 M_☉ was constructed (Espaillat 2009). The disks have α = 0.01, mass accretion rate $\dot{M}$ ranging from 10⁻¹⁰ to 10⁻⁷ M_☉ yr⁻¹, and a dust depletion parameter of 0.001, 0.01, 0.1, and 1 (no settling). The continuum indices n_{6 − 13} and n_{13 − 31} were computed for disk inclinations of 20°, 40°, 60°, and 80°.

We compare these models with our observations in Figure 11. Also in Figure 11 are the observed data for Taurus (middle panel; Furlan et al. 2006) and our disk sample in NGC 1333 (bottom panel). Furlan et al. (2005, 2006) and D'Alessio et al. (2006) demonstrate that as the settling parameter decreases, n_{6 − 13} and n_{13 − 31} also decrease for the models described above; therefore, lower values of n_{6 − 13} and n_{13 − 31} represent more settled disks overall. As discussed previously, it is important to remember that these indices characterize the settling of dust in the upper layers of different spatial regions of the disk. Overall the objects in NGC 1333 appear to be well explained by models with  = 0.01–0.001. Dust depletion factors of 0.01 or 0.001 mean that 99%–99.9% of the small grains, presumably from the ISM, have grown and sunk to the disk midplane, where they cannot be easily seen. This leads to an optically effectively flatter disk, with a bluer spectrum, as discussed by Furlan et al. (2005, 2006, 2009, 2011) and D'Alessio et al. (2006). We find no statistically significant difference in the distribution of these continuum indices, as discussed in Section 5.4.

Zoom In Zoom Out Reset image size

5.2.1. Continuum Indices and Disk Geometry

Another continuum index, n_{2 − 6}, characterizes the excess emission from the innermost region of the disk (≲1 AU; McClure et al. 2010). A disk that displays a large excess at longer wavelengths but has a very negative n_{2 − 6} is similar to the SED of CoKu Tau/4 in Taurus, which has been modeled to have the geometry of an inwardly truncated disk with a completely cleared-out inner hole with radius ∼10 AU (D'Alessio et al. 2005; Calvet et al. 2005). TDs, as defined in this analysis, have little to no excess from a disk at short wavelengths (<8 μm) and show IR excess emission from an optically thick disk at longer wavelengths. PTDs, like UX Tau A, have an inner disk at ≲1 AU that would result in a short-wavelength excess similar to that of a full disk though they possess one or more gaps and an outer optically thick disk. It is important to note that n_{2 − 6} does not necessarily characterize the amount of excess above the photosphere.

The n_{13 − 31} continuum index, as discussed previously, probes properties of dust in a radially continuous disk about 1–10 AU from the star (D'Alessio et al. 2006). For objects with interrupted radial structure, such as TDs, n_{13 − 31} characterizes the excess originating from the outer disk wall, as the SED will show a steep rise at longer wavelengths that is essentially the superposition of a blackbody at the temperature of the disk wall onto the spectrum of the star (D'Alessio et al. 2006; Espaillat et al. 2007; Espaillat 2009). For TDs whose outer disk wall is heated due to stellar irradiation and accretion shocks on the stellar surface, we expect this warmer optically thick dust to manifest itself in the spectrum as a high n_{13 − 31} (D'Alessio et al. 2006; Furlan et al. 2006).

In this way, we can identify objects that have low n_{2 − 6} and high n_{13 − 31} as TDs (McClure et al. 2010; Manoj et al. 2011). PTDs in this analysis are disks with one or more gaps, resulting in at least an inner disk and an outer disk. Like TDs, PTDs with large gaps will have an increase in the continuum from 13 to 31 μm from the wall of the outer optically thick disk. For this reason, a large n_{13 − 31} can identify PTDs with large (tens of AU) gaps, though PTDs with one or more smaller gaps would appear to have an n_{13 − 31} similar to a radially continuous disk (Espaillat et al. 2007; Espaillat 2009; Furlan et al. 2009; McClure et al. 2010).

We calculate n_{2 − 6} for our sample of disks using Equation (2) with fluxes at 2 and 6 μm derived from extinction-correlated 2MASS K_s band photometry and IRS spectra from 5.4–6.0 μm and report these values along with their uncertainties in Table 3. In Figure 12, we plot n_{13 − 31} versus n_{2 − 6} for the sample of disks as done in McClure et al. (2010) to separate out objects whose radial structure is broken up by holes or gaps. As in Ophiuchus, most of the disks in NGC 1333 fall within −2 < n_{2 − 6} < −0.7 and −1.33 < n_{13 − 31} < 0.3. In Ophiuchus and Taurus, PTDs and TDs are found in the upper octile in n_{13 − 31} with TDs in the lower octile of n_{2 − 6} and PTDs between the lower octile and the median (McClure et al. 2010). Applying this classification to NGC 1333, #137, the outlier in the TD range of n_{12 − 20}, also stands out here in n_{13 − 31} and has a very steep slope from 2 to 6 μm, placing it in the TD range as defined by Ophiuchus and Taurus of the n_{13 − 31} versus n_{2 − 6} plot. Six objects (#55, #75, #101, #110, #131, and #134) fall within or very close to the supposed PTD range of n_{2 − 6}, as they are outliers in n_{13 − 31}. All of these objects appear to be quite faint from 1 to 35 μm with νF_ν ∼ 10⁻¹¹ to 10⁻¹² erg s⁻¹ cm⁻². Two objects, #55 and #75, do not appear to have a rising SED at longer wavelengths but instead have an n_{13 − 31} outlier as a result of noise in the spectrum >30 μm. The three objects in the PTD range with the largest values of n_{13 − 31}, #101, #110, and #131, appear to have SEDs similar to TDs, as they show steep slopes from 2 to 6 μm and flux excess from 13 to 31 μm. However, #131 has an excess <8 μm similar to that of a disk and is likely not a TD. These objects do have a steep slope from 2 to 6 μm, with n_{2 − 6} ∼ −1.7. A small amount of dust in a central clearing would cause the slope to be slightly flatter, steeper than a full disk but flatter than a photosphere. For these reasons we consider #101 and #110 as possible TDs.

Zoom In Zoom Out Reset image size

The object with the lowest value of n_{13 − 31} in the PTD range, #134, does not look like any known PTDs (Espaillat 2009; Espaillat et al. 2010), but it shows a short-wavelength excess similar to a disk, a significant dip in flux between 10 and 20 μm, and starts to rise out to about 30 μm. We do not classify #134 as a TD or PTD though note that the strong dip in the middle of the IRS spectrum possibly indicates interrupted radial structure and makes this object a good candidate for further study. Winston et al. (2010) identified #136 as a TD based on its IRAC and MIPS photometry. Though it is not an outlier in n_{13 − 31}, #136 has the steepest n_{2 − 6} (−2.6), which supports its classification as a TD. We will discuss #136 further in Section 5.2.2. As mentioned previously, n_{2 − 6} is a slope and does not necessarily characterize the amount of excess above the photosphere. For this reason it may not be able to identify some types of PTDs and TDs.

There are two blue outliers (#53 and #123) in n_{13 − 31} with an n_{13 − 31} below −4/3, the limit for an infinite geometrically thin, optically thick disk. Low values of n_{13 − 31} have been found in models of outwardly truncated disks (McClure et al. 2008), and for this reason we suggest that #53 and #123 have outwardly truncated disks. Important to note is that objects may also have an n_{13 − 31} below −4/3 if they are no longer optically thick, a sign of a more evolved state. However, we do not believe this is the case for these objects as they are not identified as Class III from Section 4.1 or photosphere dominated from Section 4.2. The presence of a companion is known to truncate a disk (e.g., Artymowicz & Lubow 1994) as in the case of SR20 in Ophiuchus (McClure et al. 2008). The lack of established multiplicity information for the sample does not allow us to address definitively the multiplicity for these truncated disks in NGC 1333.

5.2.2. Silicate Emission and Disk Geometry

The broad silicate emission features in the spectrum of disks at 10 and 20 μm arise from optically thin dust at the disk surface and can be used to probe the dust properties in the upper layers of the disk (Malfait et al. 1998; Natta et al. 2000; D'Alessio et al. 2006). For objects whose structure deviates from being radially continuous, such as TDs, the emission seen in these features comes from optically thin dust within the gap or hole, such as in GM Aur, or from optically thin dust in the surface layers of the outer disk wall, such as in CoKu Tau/4 (D'Alessio et al. 2005). Therefore, measuring the strength of the 10 and 20 μm silicate emission can also be used to probe the radial structure of the disk.

The equivalent width of the silicate emission features in the IRS spectra is a measure of the optically thin emission per unit area of optically thick emission from dust in the disk. The equivalent width is characterized by the sum of the excess flux from the silicate feature at each wavelength, F_λ, with respect to the flux of the underlying continuum, F_{λ, cont}. Equivalent width is given by

$$\begin{equation} {\rm EW}_\lambda = \int _{\lambda _1 }^{\lambda _2} \frac{F_{\lambda } - F_{\lambda ,{\rm cont}}}{F_{\lambda ,{\rm cont}}} d \lambda. \end{equation} \tag{ 3 }$$

In our analysis we adopt λ₁ = 8 μm and λ₂ = 13 μm and λ₁ = 16 μm and λ₂ = 28 μm in order to measure the 10 and 20 μm silicate emission features, respectively. We adopt positive values of EW_λ for emission features.

While the equivalent width expresses the strength of the silicate feature in units of its underlying continuum, the integrated continuum-subtracted flux of the silicate feature allows us to probe how much optically thin dust is actually producing this emission. For this reason, we also calculate the flux of the 10 and 20 μm silicate features using

$$\begin{equation} F_\lambda = \int _{\lambda _1 }^{\lambda _2} (F_{\lambda } - F_{\lambda ,{\rm cont}}) d \lambda \end{equation} \tag{ 4 }$$

with the same limits as above. The EW and integrated flux are calculated in the same way as in Watson et al. (2009), Furlan et al. (2009), McClure et al. (2010), and Manoj et al. (2011).

We determine the continuum for the calculation of both properties by fitting the spectrum with a third- to fifth-order polynomial with the fit anchored at 5.61–7.94 μm, 13.02–13.50 μm, 14.32–14.80 μm, 30.16–32.19 μm, and 35.07–35.92 μm, which Watson et al. (2009) showed to be free of distinct silicate emission features. For objects where a single polynomial did not produce an adequate fit, we used two polynomials to fit the continuum from 5–14 μm and 14–36 μm. For a few cases a smoothly varying fit that did not exceed the extinction-corrected flux of the spectrum could not be found. For these objects the equivalent width and integrated flux from the portion of the spectrum for which a satisfactory fit could not be found are left out of our analysis. We calculate EW₁₀, EW₂₀, F₁₀, and F₂₀, as well as their uncertainties, for 44 of the 57 disk-dominated objects with no evidence of an envelope as discussed in Section 5.1; these values are listed in Table 4. The main source of uncertainty in these calculations results from the continuum fit; the uncertainties noted in Table 4 are calculated assuming an uncertainty in the continuum fit of 10%, though in some cases the uncertainty is likely to be larger. A continuum fit successfully extracting the 33 μm crystalline forsterite feature could not be found for all but a few of the spectra as a result of significant uncertainty at wavelengths ≳30 μm as discussed in Section 5.1; therefore, we do not discuss this feature in our analysis.

Table 4. Dust Processing Indicators

Index	EW10a	${\sigma _{{\rm EW}_{10}}}$a	EW20a	${\sigma _{{\rm EW}_{20}}}$a	F10b	${\sigma _{F_{10}}}$b	F20b	${\sigma _{F_{20}}}$b	F11.3/F9.8	$\sigma _{F_{11.3}/F_{9.8}}$
18	1.92	0.08	1.52	0.14	31.30	0.89	14.04	0.89	0.53	0.04
44	0.60	0.06	...	...	0.21	0.02	0.27	0.02	0.64	0.14
46	1.14	0.08	1.16	0.27	0.13	0.01	0.03	0.01	0.77	0.09
48	1.66	0.08	1.11	0.17	0.89	0.03	0.17	0.03	0.48	0.04
49	1.07	0.07	1.86	0.18	1.25	0.07	0.59	0.07	0.46	0.06
50	4.66	0.11	1.88	0.15	24.97	0.31	3.77	0.31	0.52	0.03
51	1.99	0.09	2.25	0.35	0.23	0.01	0.08	0.01	0.99	0.07
52	6.05	0.13	3.17	0.17	14.97	0.14	3.07	0.14	0.43	0.02
53	1.83	0.08	...	...	0.55	0.02	...	...	0.73	0.06
54	3.20	0.10	1.35	0.22	2.17	0.04	0.25	0.04	0.55	0.03
55	1.19	0.08	...	...	0.10	0.01	...	...	0.54	0.11
57	6.44	0.14	3.35	0.16	13.23	0.11	3.92	0.11	0.34	0.02
58	0.81	0.07	...	...	0.29	0.02	...	...	1.00	0.14
59	2.94	0.09	1.43	0.15	6.83	0.13	1.53	0.13	0.51	0.03
61	2.14	0.08	1.88	0.15	11.87	0.31	4.46	0.31	0.50	0.04
64	2.83	0.09	2.85	0.16	6.75	0.13	3.15	0.13	0.31	0.02
65	3.13	0.12	1.47	0.15	7.99	0.28	2.53	0.28	0.46	0.03
73	8.35	0.16	4.87	0.18	195.40	1.26	77.74	1.26	0.28	0.01
74	2.21	0.08	1.71	0.15	5.56	0.14	1.62	0.14	0.47	0.03
75	3.37	0.17	2.94	0.29	0.54	0.02	0.30	0.02	0.71	0.06
78	3.78	0.10	3.87	0.19	5.05	0.08	2.15	0.08	0.48	0.03
84	1.14	0.07	...	...	0.98	0.05	...	...	0.36	0.04
88	2.35	0.08	1.08	0.14	47.59	1.07	5.60	1.07	0.33	0.03
94	2.31	0.31	1.81	0.34	0.30	0.04	0.10	0.04	0.42	0.06
99	1.14	0.07	4.42	0.24	0.73	0.04	0.91	0.04	0.44	0.05
101	1.93	0.11	...	...	0.41	0.02	...	...	0.46	0.05
106	0.60	0.06	1.18	0.15	2.19	0.20	1.34	0.20	0.39	0.10
110	0.74	0.20	...	...	0.06	0.02	...	...	...	...
111	0.63	0.06	1.00	0.14	10.30	0.91	4.71	0.91	0.49	0.09
116	4.33	0.22	4.93	0.48	40.16	1.87	24.16	1.87	0.40	0.05
118	2.36	0.09	2.97	0.25	0.52	0.02	0.24	0.02	0.75	0.06
119	1.56	0.08	0.60	0.16	2.57	0.10	0.32	0.10	0.56	0.05
120	3.04	0.10	3.25	0.20	1.42	0.03	0.59	0.03	0.48	0.03
121	0.73	0.07	0.93	0.27	0.16	0.01	0.04	0.01	0.60	0.08
123	4.63	0.11	0.44	0.15	21.68	0.26	1.05	0.26	0.55	0.03
125	1.82	0.08	0.97	0.18	1.02	0.03	0.20	0.03	0.85	0.07
126	1.12	0.08	4.46	0.25	0.60	0.03	1.02	0.03	0.30	0.05
127	2.64	0.09	1.77	0.19	1.15	0.03	0.28	0.03	0.45	0.03
128	1.30	0.09	...	...	0.13	0.01	...	...	0.97	0.15
131	2.23	0.09	0.77	0.25	0.15	0.00	0.04	0.00	0.47	0.04
133	2.26	0.09	0.68	0.15	2.28	0.06	0.58	0.06	0.51	0.04
134	2.94	0.09	5.51	0.34	0.26	0.01	0.24	0.01	0.45	0.03
136	1.04	0.08	7.28	0.35	0.43	0.03	0.82	0.03	0.60	0.09
137	3.08	0.10	2.00	0.15	0.67	0.01	1.16	0.01	0.65	0.04

Notes. Objects with unsatisfactory continuum fit are indicated with ellipses. ^aEW_λ values are in units of μm. ^bF_λ values are in units of 10⁻¹² × erg s⁻¹ cm⁻².

Download table as:  ASCIITypeset image

The EW can be used to identify PTDs that have gaps that result in a deficit of optically thick emission (Furlan et al. 2009). Equivalent width measures optically thin emission with respect to optically thick emission; therefore, for the same amount of optically thin emission, a lack of optically thick emission would manifest itself in a large EW. If the gaps have some amount of optically thin dust, this will further enhance their ability to be identified by having a large EW. In Figure 13 we plot n_{13 − 31} versus EW₁₀ for our sample and the Taurus sample from Furlan et al. (2006), along with a polygon indicating the boundaries in these indices that correspond to boundaries found by the models of radially continuous, irradiated disks described in Section 5.1. These boundaries are a function of the parameters varied in the models (e.g., $\dot{M}$, M_☉, and disk inclination). We also note that the EW₁₀ depends on the settling parameter, , which we also varied in the models. Objects with EW₁₀ and n_{13 − 31} that fall outside of this polygon have an interrupted radial structure, meaning that the disk has one or more gaps or is inwardly or outwardly truncated. This is reflected by the positions of the well-modeled TDs and PTDs in Taurus such as CoKu Tau/4, GM Aur, and LkCa 15. Figure 14 demonstrates the areas outside the polygon of radially continuous disk models that generally correspond to disks with altered radial structure in Taurus, Chamaeleon I, and Ophiuchus that have been modeled (Espaillat et al. 2007; McClure et al. 2008; Kim et al. 2009).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

First, we examine the disks that are extreme outliers in EW₁₀ through the criterion that their EW₁₀ is in the top 12.5%. Using this restrictive criterion, we identify #50, #52, #57, #73, and #116 as PTDs. As discussed previously, objects that are likely to be TDs possess a large n_{13 − 31}, greater than that of the radially continuous disk models. While #110 falls within the PTD range in Figure 12, it is certainly not an outlier in EW₁₀, though it is an outlier in n_{13 − 31}, reinforcing its previous classification as a TD. Curiously, #110 has a position in Figure 13 similar to UX Tau A, a PTD that has been modeled to have no optically thin dust in the gap, resulting in a small EW₁₀ (Espaillat et al. 2007). We cannot confirm or exclude this possibility without detailed modeling. All disks that were in the PTD range of Figure 12, except for #110, fall within the polygon of the models in Figure 13. The other TD in the PTD range from Figure 12, #101, falls just within the polygon of radially continuous, irradiated disk models in Figure 13 though its uncertainty in n_{13 − 31} crosses the top boundary. Of the other disks in the PTD range in Figure 12, #75 and #134 have EW₁₀ greater than the median, though they are not in the upper 12.5% and therefore are not large enough to be considered significant outliers. The classification of #137 as a TD from Figure 6 and Figure 12 is bolstered here as its position in Figure 13 places it outside the polygon of radially continuous disk models and in the TD range. Our analysis agrees with that of Merín et al. (2010), which classifies #137 (known therein as source 5 or SSTc2d J032929.3+311835) as a "cold disk."

The TD suggested by Winston et al. (2010) that did not fall in the TD range of Figure 12 is also not unique here in that it falls within the polygon of models. Though this object has a "normal" n_{13 − 31}, examination of its SED shows that it is a TD. The SED of #136 is photospheric <9 μm, and the continuum slightly rises toward the long-wavelength side of the spectrum. We interpret this as a sign of a large central hole size, relatively cool outer disk wall, and a settled low-mass outer disk. Both the large inner hole and the low-mass settled outer disk prevent #136 from being identified as a TD based on n_{13 − 31}.

The case of #136 suggests that n_{13 − 31} can pick out some TDs though others can have an n_{13 − 31} similar to that of a radially continuous disk. The spectral indices n_{13 − 31} and n_{2 − 6} are only good at identifying TDs with a certain range of hole sizes. If the hole is too small (≲1 AU) or too big (≳50 AU), the indices will not identify them as TDs. It is also possible that n_{13 − 31} can improperly identify radially continuous disks as TDs; therefore, individual SED inspection is key for proper classification. Further discussion of the classification TDs and PTDs in NGC 1333 is found in Section 6.

5.2.3. Analysis and Comparison of Silicate Emission Features

With respect to this classification we further analyze the radial structure of disks in NGC 1333 using F₁₀, EW₁₀, and EW₂₀. Figure 15 shows F₁₀ (top) and EW₂₀ (bottom) plotted against EW₁₀. As expected, F₁₀ and EW₁₀ are well correlated, with Pearson's linear correlation coefficient r = 0.65, and a probability p_rand = 0.0024% that this correlation could have been generated from a random distribution (see Watson et al. 2009). A small value of p_rand, namely, ⩽1%, indicates that the quantities are significantly correlated. The F₁₀ versus EW₁₀ plot shows that most of the PTDs have large values of F₁₀, as well as an enhanced EW₁₀. A handful of disks have similar F₁₀ values as the PTDs but have much smaller EW₁₀. Objects with enhanced EW₁₀ in Ophiuchus (McClure et al. 2010) and Chamaeleon I (Manoj et al. 2011) show similar values of F₁₀ to radially continuous disks, which suggests that these objects have a deficit in optically thick emission as a result of a gap in their disk. Interestingly, the disk with the highest F₁₀ also has the highest EW₁₀ in NGC 1333, which is not seen in Ophiuchus and Chamaeleon I. The TDs in NGC 1333 have low F₁₀, also the case for TDs in Ophiuchus and Chamaeleon I. The flux at 10 μm comes from the inner regions of the disk; therefore, a significantly cleared inner disk, as is the case for TDs, would result in a lack of emission at 10 μm.

Zoom In Zoom Out Reset image size

The bottom panel of Figure 15 compares EW₁₀ and EW₂₀ and shows that the two trace each other well, with r = 0.32 and p_rand = 6.1%. All of the PTDs, which have enhanced EW₁₀, fall above the median for EW₂₀. However, most of the PTDs have a lower EW₂₀ than EW₁₀, which is also seen in Chamaeleon I. A large gap would affect the underlying optically thick emission contributing to EW₂₀, as well as EW₁₀. This trend in lower EW₂₀ with EW₁₀ outliers possibly suggests that most of these EW₁₀ outliers have small gaps in the inner disk that only affect the EW₁₀. As discussed in Section 5.2.1, the small nature of these gaps would cause them to be more difficult to identify as PTDs with the n_{2 − 6} versus n_{13 − 31} analysis (Figure 12). As was the case in Ophiuchus and Chamaeleon I, the trends of the continuum indices, fluxes, and flux ratios with EW₂₀ are weaker versions of the trends with EW₁₀. Most of the scatter in EW₂₀ is due to the uncertainty in the continuum fit for the 20 μm feature.

One probe of the degree of crystallization and grain growth in the inner 1–2 AU of the disk is the flux density ratio F_11.3/F_9.8 (Przygodda et al. 2003; van Boekel et al. 2003, 2005; Kessler-Silacci et al. 2006; Honda et al. 2006; Bouwman et al. 2008; Olofsson et al. 2009). F_11.3/F_9.8 quantifies the shape of the 10 μm silicate feature and as such tracks the degree of processing of the dust that creates the silicate emission features. The higher the ratio, the more the dust has been processed (more grain growth and/or crystallization), with a "pristine" profile defined as one similar to the optically thin interstellar-grain profile with F_11.3/F_9.8 = 0.34. We calculate F_11.3 and F_9.8 according to Equation (4) with λ₁ = 10.9 and 9.4 μm and λ₂ = 11.7 and 10.2 μm, respectively. F_11.3/F_9.8 for the sample is listed in Table 4. An F_11.3/F_9.8 representative of the processing of dust in the disk of object #110 was impossible to obtain due to contamination of F_9.8 by residual molecular hydrogen emission not removed by sky subtraction.

In Figure 16 we plot F_11.3/F_9.8 against EW₁₀. The objects with the highest values of F_11.3/F_9.8 correspond to the lowest EW strength, which is characteristic of radially continuous disks. This inverse relationship has Pearson r = −0.34 and p_rand = 2.4%. We also note the trend that as EW₁₀ decreases the frequency of disks with high F_11.3/F_9.8 increases. This trend is seen in Taurus (Furlan et al. 2006), Ophiuchus (McClure et al. 2010), and Chamaeleon I (Manoj et al. 2011), as well as in other regions (van Boekel et al. 2005; Kessler-Silacci et al. 2006; Bouwman et al. 2008; Apai et al. 2005; Pascucci et al. 2009). Also, the PTDs have relatively low F_11.3/F_9.8, with most of them falling below the median and close to the pristine line. This would imply in general that the dust seen in the PTDs has experienced less grain growth and/or crystallization than the radially continuous disks. Therefore, the optically thin dust at 1–2 AU in PTDs is more pristine than radially continuous disks. The trend in low F_11.3/F_9.8 (more pristine shape) in outliers in EW₁₀ is also seen in models of PTDs (Espaillat 2009) and observation of PTDs in Taurus (Sargent et al. 2009), Ophiuchus (McClure et al. 2010), and Chamaeleon I (Kim et al. 2009). Modeling suggests that the emission responsible for the shape and strength of the 10 μm silicate feature in these disks comes from the optically thin dust within the gaps or in the surface layers of the outer disk wall (Espaillat 2009). As also seen in Chamaeleon I, the disks with n_{13 − 31} below −4/3 in NGC 1333 show a slight preference toward higher values of F_11.3/F_9.8, which implies that these truncated or optically thin disks have a higher degree of dust processing than TDs or PTDs.

Zoom In Zoom Out Reset image size

We use a one-dimensional (1D) two-sided Kolmogorov–Smirnov (K-S) test to determine if the distribution of n_{2 − 6}, n_{6 − 13}, and n_{13 − 31} in NGC 1333 is drawn from a significantly different population from that of Taurus (Furlan et al. 2006), Chamaeleon I (Manoj et al. 2011), or the L1688 cloud in Ophiuchus (McClure et al. 2010). To do this, we took our sample of 36 Class II disk-dominated objects and did not include objects that were identified previously as TDs, truncated disks, or PTDs, leaving us with 27 objects representative of the "normal" radially continuous disk population in NGC 1333. We culled the data from Taurus, Chamaeleon I, and L1688 in the same way. The distribution of these indices is shown in Figure 17 and p-values inset from a 1D two-sided K-S test comparing NGC 1333 to each region.

Zoom In Zoom Out Reset image size

We find no statistically significant difference between NGC 1333 and Chamaeleon I, Taurus, or L1688 in the distribution of n_{6 − 13} and n_{13 − 31}. These statistical tests imply that there is a similar dispersion in these continuum indices in all four regions. The distribution implies a heavy concentration in the = 0.001–0.01 range and a tail that reaches 0.1, but very few objects that demand a well-mixed disk. We do, also, find a difference between NGC 1333 and Taurus and L1688 in the distribution of n_{2 − 6} (p-value = 0.07 and 0.01, respectively). This would imply that there are significant differences in the slope of the excess emission in the innermost regions of disks in these regions. NGC 1333 and Chamaeleon I have distributions of n_{2 − 6} skewed toward bluer (smaller) values. In Chamaeleon I, this is attributed to a spectral type difference between these two regions (Manoj et al. 2011). The same argument can be made for NGC 1333, which has many more objects in the M3–M8 than K5–M2 spectral type range and a median spectral type of M3. Taurus and Ophiuchus both have a median spectral type of M0. The median n_{2 − 6} value for M3–M8 objects in NGC 1333 is −1.64, which, when compared to the median n_{2 − 6} for K5–M2 objects of −1.50, shows that the prevalence of M3–M8 stars in our sample skews the n_{2 − 6} distribution of NGC 1333 to the smaller side.

Trends with model parameters within the polygon from Figure 13 are shown with arrows in Figure 14 (Furlan et al. 2009; D'Alessio et al. 2006; Espaillat 2009). To determine if there are differences in the scatter of points within this polygon, a 2D two-sided K-S test was performed for only the disks in Taurus, Chamaeleon I, Ophiuchus, and NGC 1333 that fall within this polygon and have a spectral type later than M0. We find no statistically significant difference in the 2D distribution of points within this polygon between NGC 1333 and Taurus, Chamaeleon I, or L1688, with p-values of 0.36, 0.20, and 0.21, respectively.

We construct a representative median for our sample in NGC 1333 from the 35 objects classified as both Class II and disk dominated observed in program 40525. The objects in our median span a spectral type range from A3 to M8, with 7 between K4 and M2 and 20 between M3 and M8 and a median spectral type of M4. We also include seven objects with unknown spectral type. The extinction measured in A_V (derived from the R_V = 5 version of the Mathis (1990) extinction curve) ranges from 0 to 22.3, with a median of 3.1, and includes seven objects with A_V > 10. We follow the same procedure as Furlan et al. (2009) to construct our median. We normalized the dereddened IRS spectrum of all 35 Class II disks through multiplication by a scalar factor equal to the median dereddened H-band flux of the 35 objects (0.039 Jy) divided by the H-band flux of each object, computed the median of the sample at each wavelength, and converted to νF_ν. This method effectively normalizes the excess emission for differences in stellar luminosity for stars with spectral types K–M (D'Alessio et al. 1999), which compose the majority of our median sample (27/28 objects with known spectral type).

In Figure 18 we compare our median to the medians of K5 to M2 stars in Taurus, L1688, and Chamaeleon I previously presented in Furlan et al. (2009), each scaled to the median H-band flux of Taurus. The shape of the NGC 1333 median is quite similar to that of the other medians, across the IRS range. After scaling to the median H-band flux of Taurus, the NGC 1333 median is still fainter than the others. This can be attributed to the spectral type differences between these samples. The Taurus, L1688, and Chamaeleon I medians are for objects with K5 to M2 spectral types, while our sample has six objects in the K5 to M2 range and 20 have spectral types later than M2. Late-type stars heat their disks less as they are less luminous and therefore irradiate their disks less, as well as have lower mass accretion rate, resulting in weaker accretion heating. The cooler nature of late-type stars manifests itself as fainter IR excess emission.

Zoom In Zoom Out Reset image size

We compare the shape of the 10 μm silicate features of the medians of Taurus, L1688, Chamaeleon I, and NGC 1333 in Figure 19. We do this in the exact same way as described in Furlan et al. (2009). The only difference between our Figure 19 and Figure 7 of Furlan et al. (2009) is the addition of the 10 μm silicate feature of NGC 1333 to the plot. Also plotted in Figure 19 is a "pristine," ISM-like profile that is represented by the average of the continuum-subtracted and normalized 10 μm silicate feature of LkCa 15 and GM Aur in Taurus (Sargent et al. 2009). This pristine profile is scaled to the mean 9.8 μm flux of the four regions. All four regions show evidence for grain growth and/or crystallization as they show a silicate profile wider than the "pristine" profile for ISM-like small amorphous silicate grains (sub-μm). NGC 1333 shows a long-wavelength wing similar to that of the other regions, which suggests similar grain sizes in the optically thin dust.

Zoom In Zoom Out Reset image size