A Spitzer survey of mid-infrared molecular emission from protoplanetary disks I: Detection rates

The young star+disk systems examined in this paper are drawn from several Spitzer-IRS short-high (SH) and long-high (LH) observing programs. Most of the Chamaeleon and Lupus sources were observed as part of a deep T Tauri star survey (PID 50641, PI. J. Carr, hereinafter GO-5) that was designed to obtain very high S/N spectra by maximizing redundancy and using dedicated background observations for each target. The subset of the GO-5 sample presented here was constructed by selecting disks spanning the known range in disk infrared colors and stellar X-ray luminosity, in particular those previously observed in the the cores to disks (c2d) Legacy program (Evans et al., 2003) that showed some hint of water emission (but which were not observed with sufficient redundancy and high enough S/N to produce firm detections). Another subset of targets from the high quality survey focused on the Taurus cluster, and will be the subject of a separate study (Carr et al., in prep). The sample covered here is summarized in Tables 1 and 2.

The non-Taurus ‘high-quality’ sample is supplemented by additional archival c2d spectra; and for comparison with the T Tauri stars that dominate the c2d selection, spectra from the Herbig Ae/Be star survey by J. Bouwman (PID 3470 Boersma et al., 2008) were extracted from the Spitzer archive and reduced using the same procedure as the remaining sample. This ensures that the sample includes a significant number of disks associated with stars of spectral types spanning from $M$ to $B$. The full sample includes disks from the young clusters in Perseus, Taurus, Chamaeleon, Lupus, Ophiuchus and Serpens, but is not complete in any strict sense. For instance, the sample represents a selection of the “best” available Spitzer spectra, and thus tends to be biased toward brighter, isolated, disks with low background emission (the c2d selection criteria are described in Evans et al., 2003). Nevertheless, given the size of the sample (75 disks), relative to the total number of class II disks brighter than $\sim$100 mJy at 8 $\mu$m in the major nearby star forming clouds surveyed by the c2d program ($\sim 120$, Evans et al., 2009), there is little room for strong biases, except for those imposed by the IRS sensitivity limit. A disk with a flat spectrum and a flux of 100 mJy at 8 $\mu$m roughly corresponds to a star with mass $\sim$0.1 M${}_{\odot }$ at the distance of Ophiuchus (125 pc), and 0.5 $M_{\odot }$ at the distance of Serpens (415 pc Dzib et al., 2010), using the evolutionary tracks of Siess et al. (2000).

Figure 1 shows the distribution of spectral indices $n_{13\mu m-31\mu m}$ of the sample, as defined in Kessler-Silacci et al. (2006) and Furlan et al. (2006). A few sources are labeled, for reference, including the transitional disks T Cha, LkHa 330, HD 135344B and SR 21. It can be seen that the sample spans the range from disks with slopes that are nearly photospheric, $n_{13\mu m-31\mu m}$$\sim$-2, to transitional disks with large inner holes and rising fluxes into the far-IR $n_{13\mu m-31\mu m}$$\sim$2. The majority of the disks center around $n_{13\mu m-31\mu m}$$\sim$-1 – 0, a distribution quite similar to that of the general disk populations of the included star forming regions (Furlan et al., 2009). The Herbig Ae/Be stars tend to be brighter than the T Tauri stars, but not enough to reflect their much higher luminosities – an indication that the Herbig Ae/Be stars are, on average, located at larger distances than the T Tauri stars. The sample also shows a deficit of transitional disks with low 31 $\mu$m flux. This may, to some extent, reflect an intrinsic property of transitional disks, although one canonical transitional disk that is not in our sample — GM Aur — only has a $31\mu$m brightness of $\sim 1$ Jy (Furlan et al., 2006).

The Spitzer-IRS spectra are reduced using IDL scripts optimized for deep integrations with high redundancy, i.e., up to 56 individual spectral frames for our data set, and dedicated background observations. The procedure is similar to that used by Carr & Najita (2008) for AA Tau, and the resulting spectra have been benchmarked to those results. The reduction begins with the droop frames. On-source spectra as well as background exposures are co-added. To avoid inappropriately weighting pixels at the edges of the entrance slit, the flat field is divided by low-order polynomial fits in both the dispersion and cross-dispersion directions. This produces a flat field that corrects the pixel-to-pixel response only. The background exposure is used to detect rogue pixels by flagging pixels with values that are more than 2$\sigma$ from the mean of that pixel. Next, the background frame is subtracted from the on-source frames (two nod positions) and the bad pixels are linearly interpolated in the dispersion direction. The ability to detect transient bad pixels by using a highly redundant (i.e. including many individual readouts) background observation taken essentially at the same time as the pointed observation produces excellent results. An example LH frame, before and after bad pixel cleaning, is shown in Figure 2.

Spectra are extracted from the 2D co-added spectral images using optimal extraction (Horne, 1986). A selection of 15 standard star observations of $\delta$ Dra and $\xi$ Dra were retrieved from the archive and reduced in the same manner, with the same flat field and extraction apertures, essentially producing a data base of spectral response functions (SRFs). Because the spectral response function depends on the location of the target in the slit, a data base of SRFs allows a selection of those that best match a given science observation. The best 6-12 standard stars were chosen by minimizing the noise at 10-11 $\mu$m (SH orders 19 and 20) and a flux calibrated spectrum is produced for each standard star observation. The flux calibrated spectra are carefully defringed using IRSFRINGE (Lahuis et al., 2007) ensuring that the real structure from the complex water emission spectrum is not affected. Setting the tolerance of a defringer too low will essentially result in the data being processed by a low-pass filter, removing high frequency structure, such as that produced by the densely packed water spectrum. The behavior of IRSFRINGE was tested by visually comparing the final product with spectra that were not defringed. In practice, most orders required no defringing at all, with most of the fringing pattern being removed by the SRFs. The orders are then combined using a weighted average with the blaze function as weights. No additional relative scaling of the orders is necessary, as the match is generally within the noise. Finally, the set of standard star divided spectra that contain the least amount of residual fringes are co-added to produce a final spectrum.

For the set of spectra that do not have dedicated background observations, such as the suite of c2d spectra, a background observation taken as close in time as possible ($<$6 months, typically 1-3 months) is used instead. Because no dedicated background observations were obtained for the c2d spectra, a few sources suffer from an incomplete background subtraction. However, given the spatial undersampling of Spitzer high resolution spectra, it is difficult to subtract the background without introducing strong artifacts.

Due to the large number of lines, almost all of which are blended with other lines of similar strength, the analysis of even a single Spitzer spectrum can be a daunting task. Further, for sources of somewhat lower quality (signal-to-noise $\sim$ 100), water emission close to the detection limit may be difficult to distinguish from noise and residual fringe patterns. Therefore, for the purposes of: 1) separating disks where water emission is clearly detected from non-detections and borderline cases, and 2) defining a simple excitation temperature parameter, a few strong lines were selected that could serve as tracers. Such tracer line complexes were selected among features that are relatively isolated in the spectrum, allowing a continuum to be fitted. Experiments were made in which a generic LTE model for water (see paper II) was correlated across a wide wavelength range of the spectrum. However, it was found that this did not produce results much different from the approach of using a few line complexes. The water line complexes at 15.17 and 17.22 $\mu$m are more isolated than most, and while they are both blends of several lines, the excitation temperature of the transitions contributing to each complex fall within a narrow range. The 15.17 $\mu$m complex traces rotational $J$ quantum numbers of 8-10 with excitation energies of 2300-2700 K, while the 17.22 $\mu$m complex traces $J$s of 9-12, corresponding to slightly higher energies of 2300-3000 K. Thus, these two complexes trace roughly the same gas, and they are expected to correlate.

A detection of water emission is defined as a detection of both 15.17 and 17.22 $\mu$m complexes at the 3.5$\sigma$ level. In a few cases, the 15.17 $\mu$m line was not formally detected, but other water lines were strong enough to still warrant a clear detection. Conversely, in some cases lines were formally detected, suggesting the presence of water emission that a visual inspection could not confirm. For instance, a lack of other water lines in the same spectral range would lead to the source being flagged as a non-detection. Due to the nature of the data set, such subjective analysis of a few borderline cases is unfortunately inevitable. We aimed to be conservative to ensure that there are no false positives in the sample, at the expense of rejecting a few real detections. Hence, all detection rates can be considered lower limits. Cases that were subjectively scrutinized but where we could not unambigously confirm detections of, in particular, water are flagged in Tables 3 and 4. Selected regions, including those containing the water tracer lines of the Spitzer-IRS spectra of T Tauri disks for which water is detected are shown in Figures 6 and 7, for the SH and LH modules, respectively. For comparison, spectral regions including the water tracer lines as well as the HCN and C${}_2$H${}_2$ Q branches for the Herbig Ae/Be disk sample are shown in Figure 8, although these spectra are all non-detections. In practice, the detection rate does not depend on the exact choice of tracer line complexes, except in a few borderline cases. Detections of OH are based on the doublets around 23.0, 27.6 and 30.5 $\mu$m, spanning upper level energies of $\sim$2400-4000 K. The organics HCN/C${}_2$H${}_2$ and CO${}_2$ are detected via their characteristic Q branches. For the latter species, it is noted that there may be a bias against their unambiguous detection in sources showing strong water emission, due to blending and confusion. In addition, these species are identified based on a single blended feature making it difficult to rule out the occasional spurious detection, especially if systematic errors or residual fringing exceed calculated errors. A few such spurious detections were eliminated with visual inspection.

The general low contrast of the molecular lines at R=600 relative to the strong continuum also affects the detection rate. Specifically, the fringe and flat field residuals, in addition to the photon statistics, scale with the continuum, so even for very high theoretical signal-to-noise ratios, there is a limit to how low the line-to-continuum ratio can be to allow the detection of lines. For our data, this “fidelity” limit, i.e. the line-to-continuum contrast where the data systematics become larger than the pure photon statistics, while difficult to quantify, seems to be better than 0.5%, with several 3.5$\sigma$ detections at the 1-2% level.

Figure 9 shows the relation between the continuum level and the flux level of the water tracer. The ability to detect a line of constant flux diminishes with increasing flux level, indicating the “fidelity” limit. On the other hand, a significant number of line detections are made in T Tauri disks with apparent brightnesses higher than those of many Herbig Ae/Be disks.

Using these criteria, water emission is detected in 22 disks, HCN in 25 disks, C${}_2$H${}_2$ in 17 disks, CO${}_2$ in 20 disks and OH in 18 disks, as summarized in Figure 10 and Table 3. Figure 10 also shows the detection rates of CO in the rovibrational fundamental band at 4.7 $\mu$m. The CO observations are described in Paper II.

The 15.17 and 17.22 $\mu$m water line complexes chosen to act as signposts for molecular emission trace temperatures in the middle of the range of the 10-36 $\mu$m water lines. However, the spectra also include line complexes from strong transitions with much lower excitation temperatures. Ratios of line fluxes between the high and low excitation line complexes are an important diagnostic of the nature of the molecular emission; indeed one of the first questions that we ask, based on the present sample, is whether the excitation conditions are similar for all disks, or whether there are significant differences. As discussed in Meijerink et al. (2009) for the particular case of mid-infrared water lines, transitions with lower excitation temperature trace larger radii of the disk. Hence, line ratios can be affected by both excitation conditions and the radial abundance distribution of the species in question. In any case, it will be necessary to identify lines with the widest possible baseline in excitation energies. A good low energy tracer is the ortho $7_{25}\to 6_{16}$ line at 29.85 $\mu$m with an excitation energy of 1100 K (see Figure 3). Spectra with higher resolving power will be able to separate lines tracing levels at least down to $J=5;K_a-K_c=3$ and $E\sim 900$K, but these are too blended with higher excitation lines at the Spitzer resolution to be clean tracers of specific conditions. Very high excitation tracers include lines at 13.3 $\mu$m (J=15-16, $E=4400-5200$K) and 29.5 $\mu$m (J=16-18, $E=4600-4800$K).

Table 3 summarizes the integrated fluxes of the 15.17, 17.22 and 29.85 $\mu$m line complexes. The line fluxes were calculated by fitting a Gaussian superposed on a linear continuum to the data. The selected line complexes appear to be spectrally unresolved, so the centers and widths of the Gaussian were kept constant for each line complex, corresponding to the theoretical center of the line complex and the nominal spectral resolving power of Spitzer-IRS. The wavelength range of the fit (determining where the continuum is constrained) is indicated by the models shown in Figures 6 and 7.

Figure 11 shows the distribution of detected water sources in a 3-line plot, matching the high excitation 15.17 and 17.22  $\mu$m complexes with the low excitation 29.85 $\mu$m complex. Since the low excitation line probes the disk surface at larger radii than the high excitation lines, the location of a source in this diagram could provide an indication of the radial water abundance structure of the disk surface: the smaller the 29.85/17.22 ratio, the sharper the cutoff of water emission beyond a certain radius. This was discussed in Meijerink et al. (2009), who presented generic non-LTE models of emission spectra, one for a constant abundance of water, and one in which the surface water abundance was lowered by 6 orders of magnitude beyond a radius of 0.7 AU, corresponding to the location of the mid-plane snow-line. The arrow in Figure 11 indicates the approximate range and direction of increasing “cold finger depletion” of the surface water abundance beyond the mid-plane snowline. It is seen that essentially all disks are located between these two extreme cases, with only a few disks exhibiting line ratios that are consistent with a high abundance of water throughout the disk surface. Meijerink et al. (2009) suggests a model, based on a few H${}_2$O Spitzer spectra from Salyk et al. (2008) and Carr & Najita (2008), in which the disk surface water vapor may be strongly depleted at radii larger than that of the midplane snow line due to vertical transport of water. Purely chemical effects may also contribute to generate structure in the radial surface abundance structure. In any case, the data, as shown in Figure 11, indicate that line ratios cluster for many sources, but significant variation does exist in the sample, with one possible interpretation being that this is due to structure of the radial emission profile of water.

The temperature structure of protoplanetary disks is determined by a detailed balance of heating and cooling processes. In the surface layers, collisional exchange between the gas and the dust is not sufficient to maintain a temperature equilibrium between the dust and the gas, leading to elevated gas temperatures. The degree to which the gas is superheated depends sensitively on the efficiency of line cooling. As a consequence, the finding that the mid-infrared wavelength range in typical protoplanetary disks is blanketed in molecular lines will have important consequences for models of disk structures.

Current disk models find water cooling rates of a few $\times {10}^{-5}{L}_{\odot }$ and individual cooling line luminosities from atomic species, as well as rotational lines of CO and H${}_2$, are ${10}^{-8}-{10}^{-5}{L}_{\odot }$ (Gorti & Hollenbach, 2008). However, we measure integrated water line luminosities in the 10-36 $\mu$m region that are 1-2 orders of magnitudes higher. Table 3 presents the 10-37 $\mu$m integrated H${}_2$O line luminosities for the disks in which H${}_2$O has been detected. The line luminosities were determined by scaling the fiducial model of Meijerink et al. (2009) to the SH and LH spectra, where different scaling values were allowed for the two spectral ranges. Typical values resulting from this range from $<5\times {10}^{-4}$ to almost ${10}^{-2}{L}_{\odot }$. Extrapolating to include lines outside the Spitzer range may increase these numbers by a factor 2 or more. Consequently, the detection of the mid-infrared molecular spectra from protoplanetary disks increases the total disk-averaged cooling rates by 1-2 orders of magnitude, relative to recent models, for gas temperatures below $\lesssim 2000$K; the infrared cooling budget is dominated by water in T Tauri stars. It should be stressed, however, that most of the water lines in the Spitzer range are likely formed in the innermost ($<$ a few AU) regions of the disk, the outer regions of which may be less affected by water as an exceptionally strong coolant. An immediate consequence of more efficient cooling seems to be that the surface layers of the disks will have a larger column of relatively cooler gas, facilitating the survival, and observability, of a rich chemistry in the disk surface.

Arguably the strongest observation, based on the current sample, is that the detection rate of molecular emission from Herbig Ae/Be stars with Spitzer-IRS is apparently very low. The upper limits for detections of the water tracers are given in Table 4. In fact, out of 25 Herbig Ae stars, there is not a single detection of water, HCN, C${}_2$H${}_2$ or OH, given our detection criteria; CO${}_2$, which is detected in HD101412, provides the sole exception. In comparison, the detection rate for disks around stars with spectral types later than F is $\sim$40%. Could this effect be due to some bias? First of all, as Figure 1 shows, the sample of Herbig stars are somewhat brighter than the T Tauri stars, by up to an order of magnitude. This indicates that the signal-to-noise ratios of the Herbig star spectra are unlikely to be significantly lower than those of the T Tauri stars. Figure 12 shows a comparison between the full IRS spectra of a Herbig star, a transitional disk and two classical T Tauri stars. Further, in Figure 13 a comparison of the line-to-continuum ratios of the H${}_2$O line tracers at 15.17 and 17.22 $\mu$m is given for the full sample. It is seen than the Herbig stars have 3.5$\sigma$ upper limits on their line-to-continuum ratios that are systematically smaller by a factor of 5-10 than the ratios in T Tauri stars where water is detected.

It is interesting to note that some Herbig Ae/Be disks do show tentative low-level emission features that may be due to water and OH at longer wavelengths, specifically in the LH module. Such tentative detections made using the 24.9-25.5 $\mu$m OH/H${}_2$O complex are displayed in Figure 14 and compared to the water spectrum from TW Cha. We describe them as tentative because they do not unambiguously match a water model over a wider range of wavelengths, although that could be explained by noise or systematics that vary with wavelength or spectral order. The 25 $\mu$m complex covers the $X{\pi }_{3/2}\to X{\pi }_{3/2}10.5$ and $X{\pi }_{1/2}\to X{\pi }_{1/2}9.5$ OH lines and water lines spanning excitation energies of 1500-2500 K, or somewhat cooler than those traced by the 15.17 and 17.22 $\mu$m complexes. If real, these detections still represent line-to-continuum ratios significantly lower than those of T Tauri disks, but suggest that, with a modest improvement in data quality, molecular emission features should also be detected in Herbig Ae/Be disks. Given that the tentative detections are at the same level of the current data systematics, it is not possible to further analyze them. Their presence, however, do suggest that the optically thick H${}_2$O lines at far-infrared wavelengths may be detected in Herbig Ae/Be disks by Herschel.

While the mid-infrared molecular line-to-continuum ratios are much smaller in Herbig Ae/Be stars relative to those in T Tauri stars, is it possible that the line fluxes are similar, but veiled by the stronger continuum fluxes of the Herbig stars? Most disks around Herbig Ae/Be stars, as well most transitional disks, do show strong lines from CO at the rovibrational band at 4.7 $\mu$m (Blake & Boogert, 2004; Salyk et al., 2009, and Paper II). New ground-based CO observations of the full sample are discussed in greater detail in Paper II.

A comparison of line-to-continuum ratios for CO M-band lines from Najita et al. (2003) and Blake & Boogert (2004) demonstrates lower line-to-continuum contrast, on average, for Herbig Ae/Be vs. T Tauri disks. The reason for this difference is still unclear, but if the water and other molecules behave similarly, this may explain the general lack of high contrast molecular line emission from Herbig Ae/Be stars. However, the difference for CO is only $\sim$ a factor of 2 — significantly smaller than the minimum difference required for water.

It is therefore tempting to conclude that the observed lower line-to-continuum ratios in Herbig Ae/Be disks are produced by lower molecular abundances, although other alternative explanations are possible. Differences in dust properties could affect line strengths, or the heating rate of the molecular layer may not scale linearly with stellar luminosity. Here, we discuss these possibilities.

Transitional disks are protoplanetary disks with a deficit of infrared emission due to a paucity of dust in the inner regions (Strom et al., 1989). This class of disks is thought to be in the process of clearing out their inner regions through various processes, possibly including planet formation. There are 5 transitional disks in this sample, LkHa 330, SR 21, T Cha, HD135344B and CoKu Tau/4. However, none show H${}_2$O, HCN, C${}_2$H${}_2$ or CO${}_2$ emission at Spitzer wavelengths detectable in our data. In contrast, LkHa 330, SR 21 and HD135344B show CO in emission in the fundamental rovibrational band at 4.7 $\mu$m (Pontoppidan et al., 2008; Salyk et al., 2009). Given the detection rates of molecular emission around regular disks, the current sample of transitional disks is too small to draw any firm statistical conclusions. HD135344B is an F star, so may be of sufficiently early type to lack molecular emission for the same reasons that the A and B stars lack it. CoKu Tau/4 is not a true disk in transition since dust and gas appears to have been cleared out by a stellar companion and not as a result of disk evolution (Ireland & Kraus, 2008). The chance that the remaining 3 transitional disks do not show H${}_2$O emission by coincidence, given the detection rate of the remaining disks, is only a few percent, but we do not consider this enough to conclude that no transitional disks will show molecular emission (especially if deeper spectra are acquired).

The absence of strong molecular emission from these systems is similar to the results of high signal-to-noise IRS spectroscopy of the transitional object TW Hya (Najita et al., 2010). The spectrum shows a striking lack of strong emission features of H${}_2$O, C${}_2$H${}_2$, and HCN in the 10-20 $\mu$m region, although weak emission is detected from other molecules (H${}_2$, OH, CO${}_2$, HCO${}^{+}$, and tentatively CH${}_3$). Najita et al. describe how the lack of strong molecular emission is consistent with the possibility that the inner disk has been cleared by an orbiting giant planet, although chemical and/or excitation effects may be responsible instead.

A lack of strong molecular emission from transitional disks is consistent with significant photodestruction due to insufficient shielding from dust in the optically thin inner region of the disks. In this scenario, CO is still seen in rovibrational transitions due to a combination of survival through self-shielding and the fact that CO has a mechanism for UV fluorescent excitation that is efficient also for cold, low density gas. The latter could be the case for SR 21, for which the CO gas originates from a ring at 6-7 AU, distances at which UV fluorescence of CO is likely operating. The CO gas from HD 135344B, on the other hand, originates mostly from small radii ($<1$AU Pontoppidan et al., 2008).

A significant sub-class of disks has been identified, consisting of sources showing strong 14.98 $\mu$m CO${}_2$ emission from the Q-branch of the fundamental bending mode, $v=(0,1,0)\to (0,0,0)$, but which show no other detectable molecular emission in the Spitzer wavelength range. Because the identification with CO${}_2$ is based on only one feature, it was confirmed that the emission is present in both the A and B nod positions and that it is not due to an artifact from the edge of order 13, which begins at 15.08 $\mu$m. Some oxygen-rich asymptotic giant branch (AGB) stars also show strong CO${}_2$ emission (Justtanont et al., 1998; Cami et al., 2000), but with much higher optical depth leading to strong features from combination bands at 13.87 and 16.18 $\mu$m. These (intrinsically weaker) bands are absent from protoplanetary disk spectra, suggesting relatively low optical depth of the bending mode. The data from the CO${}_2$ disks are shown in Fig. 15, and compared to a generic disk model from (Meijerink et al., 2009) generated using the radiative transfer code RADLite (Pontoppidan et al., 2009), assuming level populations in local thermodynamic equilibrium and a CO${}_2$ abundance of $6\times {10}^{-8}$ relative to H${}_2$. This abundance is 1-2 orders of magnitude lower than that inferred from chemical models (Willacy & Woods, 2009), but given that the critical density for the CO${}_2$ bending mode upper state is high (${10}^{10}-{10}^{12}c{m}^{-3}$, Castle et al., 2006) non-LTE calculations will probably result in higher abundances. Spectra at higher resolution should reveal strong R and P branches in the 14-16 $\mu$m range. Further chemical and radiative transfer modeling is required, not so much to explain the presence of the CO${}_2$ emission, but why that from water, HCN and C${}_2$H${}_2$ is absent. Are the abundances of these molecules really lower relative to CO${}_2$, or is the absence of emission due to small differences in the density and temperature structure of the molecular surface layer of these disks causing other molecules to be subthermally excited or shielded by dust?