EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 694 (February 2025) A&A, 694 (2025) A174 Full HTML
Open Access
Issue
A&A
Volume 694, February 2025
Article Number A174
Number of page(s) 15
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/202451629
Published online 11 February 2025

© The Authors 2025

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1 Introduction

Deuterium is expected to only be created during Big Bang nucleosynthesis (BBN). Afterward, the deuterium abundance in the universe is depleted as it is incorporated into stellar interiors and burned into 3He in the process of astration (Epstein et al. 1976).

There is no known nucleosynthetic process to replenish the supply of deuterium, and thus the deuterium abundance with respect to hydrogen ([D/H]) has been proposed as a robust measure of stellar chemical evolution; with more generations of star formation, there should be an increasing destruction of deuterium by astration (Sarkar 1996). The primordial [D/H] set by BBN has been accurately determined from UV observations of ultraviolet D I and H I absorption lines in high-redshift, metal-poor quasars to be 2.53 ± 0.04 × 10−5 (Cooke et al. 2014). This value of [D/H] is in good agreement with BBN model predictions based on cosmic microwave background observations of 2.58 ± 0.13 × 10−5 (Cyburt et al. 2016).

If [D/H] can be accurately measured in different locations through the Galactic disk, it should be an extremely valuable yet relatively simple probe of Galactic chemical evolution. Models of the formation of the Galactic disk generally predict that less star formation has occurred in the outer galaxy, and thus that [D/H] should increase with Galactocentric radius (RGC) (e.g. Chiappini et al. 2002; Romano et al. 2003, 2006; van de Voort et al. 2018). This picture is complicated by other effects, such as the infall of pristine and deuterium-rich gas which has experienced little astration (Prantzos 1996; Romano et al. 2006; van de Voort et al. 2018), and infall from metal-rich/deuterium-depleted gas originating in winds driven by starbursts (Tsujimoto et al. 2010; Tsujimoto 2011). To date, however, there have been few direct constraints on [D/H] outside a few kpc from the Sun with which to constrain [D/H] throughout the Galaxy.

Early attempts at measuring [D/H] as a function of RGC were made using radio-wavelength transitions of deuterated molecules such as DCN/HCN (Penzias et al. 1977; Penzias 1979), however, such surveys are now known to suffer from chemical fractionation effects which can alter the abundance of deuterated molecules by several orders of magnitude (Wootten 1987; Awad et al. 2014; Ceccarelli et al. 2014; Tielens 2013). To date, most measurements of [D/H] in our Galaxy come from UV absorption lines of D I and H I, however, their interpretation in terms of astration is debated (see Tosi 2010 for a review). Along sight lines with low total hydrogen column density within the local bubble (log(NH) < 19.2 cm−2), [D/H] is approximately constant at 1.5 ± 0.4 × 10−5, however, beyond the local bubble, variations of a factor of 4–5 in the derived [D/H] are seen between sight lines separated by only a few hundred parsec (Linsky et al. 2006). Linsky et al. (2006) have argued that this variation is in fact the result of variations in the depletion of the gas-phase deuterium onto dust grains. Since the C-D bond is stronger than the D-H bond, it is energetically favorable for deuterium to replace hydrogen in carbonaceous dust grains or polcyclic aromatic hydrocarbons (PAHs, Jura 1982; Tielens 1983; Draine 2006; Tielens 2008). Indeed, the C-D bond of aromatic and aliphatic carbon molecules has been detected in emission (Peeters et al. 2004, 2024). Linsky et al. (2006) show that [D/H] along the different sightlines is anticorrelated with the depletion of different refractory elements (such as Fe, Si). Regions with a high [D/H] and low levels of depletion may recently have experienced erosion or destruction of dust grains from supernova shocks, thus releasing D back into the gas phase. A tentative correlation between [D/H] and the H2 rotational temperature is also observed, which Linsky et al. (2006) suggest may be the result of clouds with warmer gas having a lower level of refractory depletion.

If depletion onto dust grains is significant, the highest [D/H] values observed are most representative of the total deuterium abundance in the local galactic disk. Taking dust depletion into account, a Bayesian analysis of the UV measurements suggest the corrected true [D/H] in the local Galactic disk is ≥2.0 ± 0.1 × 10−5 (Prodanović et al. 2010). This result is supported by recent updated measurements of [D/H] which have been improved with more precise H I measurements, although the correlation between refractory element depletion and [D/H] is somewhat weaker (Friedman et al. 2023). The depletion scenario has also been argued against by Hébrard (2010), who instead derive [D/H] from measuring [D/O] with UV absorption lines of D I and O I, and then scaling by the ISM oxygen abundance. They find that [D/O] is in fact homogenous along a similar spread of sightlines, and suggest [D/H] of the local galactic disk may be a much lower value of 7 ± 2 × 10−6 instead. Regardless, UV absorption lines become impractical with distances from the Sun >2.5 kpc due to blending of the D I/H I lines that occurs with a higher Hydrogen column density (Hoopes et al. 2003), increasing extinction in the Galactic plane, and the lack of bright background stars at large distances, and are thus also not well-suited to surveys over a wide range of RGC .

An alternative simple probe of [D/H] is provided by mid- to far-infrared rotational lines of H2 and HD, the simplest deuterated molecule. In the dense environments of molecular clouds, most deuterium should be in the form of HD, and the measurement of HD lines should provide direct and accurate means of obtaining the [D/H] ratio. As opposed to more complex deuterated species, chemical fractionation effects are not expected to be significant, though the chemical conversion of HD to atomic D can be important for abundance determination (Bertoldi et al. 1999). However, the HD molecule is difficult to observe in the mid-IR due to the low Einstein coefficients Aul of order 10−5 s−1 relative to other IR-active species such as CO2 (Aul ~ 1–100 s−1 ), high upper energy levels (Eup ≥ 1900 K) requiring strong excitation conditions, and its low abundance relative to H2 . Furthermore, rotational lines of HD are found in the mid to far-IR, and are not accessible from the ground. Only a handful of such HD detections have thus ever been made with the Infrared Space Observatory (ISO), Spitzer, and Herschel. All of these detections are in relatively nearby regions, in shocks found in protostellar outflows (<500 pc; Bertoldi et al. 1999; Neufeld et al. 2006a; Giannini et al. 2011; Yuan et al. 2012), photo dissociation regions (Wright et al. 1999), supernova remnants (<2 kpc Yuan et al. 2012), or protoplanetary disks (<160 pc; Bergin et al. 2013; McClure et al. 2016). A very tentative detection of a single HD transition with ISO made toward Sgr B2 in the Galactic center was claimed by Polehampton et al. (2002), which will not be discussed further here.

The James Webb Space Telescope (JWST) now provides a unique opportunity to survey the [D/H] abundance throughout the Galaxy using MIRI/MRS to observe rotational lines of HD and H2 with a single instrument. The sensitivity and spatial resolution of MIRI/MRS (0.3–1.1″) is greatly improved compared to previous mid-infrared facilities, allowing the H2 and HD column densities and excitation conditions to be spatially mapped in great detail. Moreover, the MIRI/MRS range (4.9–28.5 µm) covers a variety of shock-tracing lines from iron and sulfur that are accessible in the mid-IR, allowing a comparison of [D/H] and the HD column densities across a variety of different shock environments and providing a potential means to test the dust depletion hypothesis. Finally, the greatly improved sensitivity of JWST allows the faint HD lines to be detected out to much larger distances in the Galaxy than possible before, potentially allowing the [D/H] abundance across the Galaxy to be better constrained.

In this paper, we report the detections of HD lines in 5 nearby bright protostellar outflows and 5 distant high mass protostars at Galactocentric radii (RGC) from 4 to 11 kpc in the JWST Observations of Young protoStars (JOYS) program (Sect. 2). We compare the emission morphology and measure fluxes for H2, HD, and various lines tracing refractory elements in high- velocity jets ([Fe II], [S I]), use rotational diagrams to measure the H2 and HD column densities or upper limits, and determine [D/H] abundance toward our sources in Sect. 3. In Sect. 4, we compare the observed line fluxes and properties of the H2 and HD gas, and compare the [D/H] abundance between our sources and with Galactic chemical evolution models. We discuss the effects of the depletion of deuterium onto the dust grains, the possible variations in [D/H] across the Galaxy in the context of the Galactic chemical evolution models, and the outlook for future JWST studies of HD in Sect. 5. A summary and conclusions are provided in Sect. 6.

Table 1

Sources analyzed for JWST MIRI/MRS H2 and HD observations.

2 Observations

Our MIRI/MRS data are from the JOYS programs 1290 (PI: E. F. van Dishoeck) and 1257 (PI: T. Ray), which target embedded protostars with the pointings generally centered on the continuum position or blue-shifted outflow lobes. Here, we included observations of 5 high-mass star-forming regions, and 5 nearby and bright low-mass sources where the outflows are well resolved. The details of the observations for each source are listed in Table 1. We select low-mass targets with multiple pointings which are bright in H2 so that non-detections of HD can provide significant upper limits on [D/H], and all high-mass targets so that a wide range of Galactocentric radii are sampled, including IRAS 23385+6053 at 11 kpc and several targets in the inner Galaxy at ~4 kpc. We exclude the high-mass source G28-S, as no continuum source is detected even up to the longest MRS wavelength (28 µm).

The data reduction details for HH 211 are described in Caratti o Garatti et al. (2024). For all other sources in program 1290, the reduction details are as follows. The data were taken using the FASTR1 readout mode and a two-point dither pattern in the negative direction which was optimized for extended sources. All three gratings (A, B, C) were used, covering the full 4.9–28.5 µm wavelength range of MIRI-MRS. The number of pointings varies depending on the angular extent of the source and is listed in Table 1. Dedicated background observations were taken for all sources except G28 IRS 2 and G28 P1 to allow for a subtraction of the telescope background and detector artifacts. However, for IRAS 23385+6053, the background was not subtracted using the dedicated background since it also contains real astronomical continuum emission and PAH features (for details, see Beuther et al. 2023). Since the focus of this paper is on emission lines, the subtraction of telescope background does not significantly alter the results. An astrometric correction was performed for IRAS 23385+6053 using GAIA DR3 stars in the parallel imaged field (Beuther et al. 2023).

The data were processed through all three stages of the JWST calibration pipeline (Bushouse et al. 2024) using reference context jwst_1210.pmap of the JWST Calibration Reference Data System (CRDS; Greenfield & Miller 2016). Both the science and background data were first processed through the Detector1Pipeline using the default parameters. Next, the data were further processed through Spec2Pipeline, which included the subtraction of the dedicated backgrounds for G31 and IRAS 18089-1732 as well as applying the fringe flat for extended sources (Crouzet et al. in prep.) and the residual fringe correction on the detector level (Kavanagh et al. in prep.). A bad-pixel routine was applied to the resulting cal files using the Vortex Imaging Processing (VIP) package (Christiaens et al. 2023). Lastly, the final data cubes were constructed for each subband of each channel individually using the Spec3Pipeline. In this final step, the master background and outlier rejection steps were switched off.

3 Analysis

3.1 Integrated line intensity maps

Spatially extended emission from many H2, HD, and atomic lines is detected toward all of our sources (Table A.1 in Appendix A). The detected HD lines include transitions with much higher excitation energies than probed by previous mid- and far-IR observations, and range from R(4) (Eup = 1895 K) as high as the R(9) transition (Eup = 6656 K). To compare the morphology of the line emission, it is useful to inspect maps of the continuum subtracted and wavelength integrated line intensity. However, the MIRI/MRS point-spread-function (PSF) size varies significantly across the full wavelength range (0.3–1.1″), complicating the comparison. Furthermore, individual spaxels in the reconstructed MIRI/MRS data cubes can show amplitude rippling along the spectral axis from undersampling of the PSF at short wavelengths and residual fringing not removed by the pipeline. To enable an even comparison of the line intensities and mitigate these artifacts, we performed a Gaussian convolution of each data cube to a resolution of 1″, which is approximately the spatial resolution of MIRI/MRS at the wavelength of the longest analyzed line, [S I] at 25.25 µm.

We estimated the continuum emission underlying each emission line by taking the median intensity over the wavelength axis in a 1000 km s−1 wide range centered around the line’s rest wavelength. For each line, we then subtract the estimated continuum and integrate over the wavelength axis in a 400 km s−1 wide range centered on the line. This integration range is chosen to capture the full extent of the line because the velocity resolution of MIRI/MRS is only ~100–200 km s−1.

We show integrated line intensity maps of selected H2, HD, and atomic lines for IRAS 23385+6053 and HH 211 for selected lines in Figs. 1 and 2, the remaining maps and sources are shown in Appendix D on Zenodo. Where present, emission from HD is seen to be brightest toward bright shocks where the high-J H2 lines are detected.

In HH211, the locations where the HD emission peaks are strongly correlated with the [S I] and [Fe II] atomic lines tracing the high velocity inner jet and bow shock. Emission from [Fe I] is also bright in the inner jet where HD emission is detected, but not in the bow shock, likely reflecting differences in the ionization fraction of Fe. Similar correlations between the HD and atomic line emission are seen for the other nearby low-mass outflows, and for the high mass protostellar system IRAS 23385+6053 where multiple large scale outflows are being driven by embedded protostars, though over far larger scales of 1000s rather than 100s of au (Beuther et al. 2023; Gieser et al. 2023). HD emission is not clearly detected in the integrated line intensity maps in any other high-mass sources of the JOYS sample despite the similar sensitivities and distances.

thumbnail Fig. 1

Integrated line intensity maps for HH 211 of various lines and the continuum at 17 µm shown with a logarithmic stretch. The maps have been smoothed to a common resolution of 1″, shown by the white circle in the bottom-left. Apertures used for spectral extraction are shown by the green circles. An index for each aperture is provided in the top-left panel. The coordinates of the aperture centers can be found in Table B.1 of Appendix B.
thumbnail Fig. 2

Same as Fig. 1, but for IRAS 23385+6053.

3.2 Aperture spectral extraction

To quantitatively compare the line fluxes between different sources and outflow regions, we extracted spectra from the convolved cubes described in Sect. 3.1. From inspection of the integrated line intensity maps in the low-mass sources, the HD lines appear brightest in the outflow jet knots and bow shock, but there are also many outflow locations with bright H2 emission in the lower energy S(1) to S(3) lines (Eup ~ 1000–2500 K) where the HD emission is faint despite a large column density of H2 likely being present. To probe the origin of the HD, we thus extracted spectra toward positions representative of the HD- bright jets and bow-shocks, as well as the fainter positions in the outflow cavities. For the nearby low-mass sources, we extracted multiple apertures to compare conditions across the outflows, while we use only a single aperture for each high-mass source owing to their much greater distance. The locations and sizes of our extracted apertures are shown as green circles overlaid on the line intensity maps described in Sect. 3.1 and Appendix D on Zenodo, and we list their positions in Table B.1 of Appendix B.

For each position, we extracted a spectrum by taking the sum of the flux in a 1″ radius aperture. As all data cubes have been convolved to the same resolution of 1″, the aperture radius was kept constant with wavelength. We then applied a residual fringe correction to the 1D spectrum (Kavanagh et. al. in prep). To measure the flux of the H2, HD, and selected atomic lines, we fit each line with Gaussian profile and a first-order polynomial to represent the continuum. When a line was not detected, we took as an upper limit on the flux Fup = 3σRMSFWHM, where σRMS is measured in a small line-free region adjacent to the expected line, and FWHM is the full width at half-maximum of MIRI/MRS spectral response. This empirical approach in determining the upper limits incorporates both the sensitivity limits of MIRI/MRS and calibration uncertainties, such as residual high-frequency fringing remaining in the spectrum.

We show line detections and Gaussian fits for two representative apertures from HH 211 and IRAS 23385+6053 in Figs. 3 and E.1 of Appendix E on Zenodo, respectively. The fluxes for all lines and apertures are reported in Tables F.1 to F.3 in Appendix F on Zenodo. The lines are generally unblended, except for the HD R(6) line which can overlap with a bright transition of H2O in some cases, and the HD R(7) line which can become partially blended with CO2 Q and P-branch transitions. Bright OH transitions are often detected close to the HD R(6) line but do not overlap it.

thumbnail Fig. 3

Spectral line fits for H2, HD, and the [S I], [Fe I], and [Fe II] atomic lines for HH 211 aperture 8. The line flux (shaded blue region) is determined by the simultaneous fit of a Gaussian profile and a first order polynomial to represent the continuum. the rest-frame wavelength of each spectral line is shown by a dashed gray line, and non-detections are shown without a shaded region for the fit. Nearby molecular lines of are labeled.

3.3 Rotation diagram analysis

We used a rotation diagram analysis to measure the column density of H2 and HD, and thus determine the [D/H] abundance or upper limits for our targets. Under the assumption that the line emission is optically thin and in local thermodynamic equilibrium (LTE), a rotation diagram provides a simple way to determine the column density and kinetic temperature of the gas (Mangum & Shirley 2015). The extinction-corrected line intensities Iul (explained below) are used to compute the quantity for the y-axis of the diagram y=ln(Nu/𝑔u)=ln(4πIulAulhv𝑔u),$y = \ln \left( {{N_u}/{g_u}} \right) = \ln \left( {{{4\pi {I_{ul}}} \over {{A_{ul}}hv{g_u}}}} \right),$

where Nu is the column density of molecules in the upper state of the transition, gu the upper state degeneracy, Aul the Einstein A emission coefficient, and hν the photon energy. The x-axis of the diagram is the upper state energy x = Eu/k in units of K, and on such a diagram a linear fit of the form y = mx + b will yield the column density of the gas and temperature via: m = −1/T and b = ln(N/Q(T)), where N is the total column density of the molecule and Q(T) is the partition function.

To convert the wavelength-integrated line fluxes from Sect. 3.2 to the intensity needed for the rotation diagrams, we assumed that the H2 or HD emission fills the aperture. This assumption should be reasonable for the nearby low-mass sources, as modeling of the H2 excitation in shocks indicates that most should be spatially resolvable by JWST/MIRI for shock densities ~106 cm−3 (Kristensen et al. 2023, Figure 9). For the distant high-mass sources the emission solid angle may be much smaller than the aperture, but we nonetheless make the same assumption. In any case the assumed emission solid angle is the same for both H2 and HD and thus there is no impact on any column density ratios or the derived [D/H]. The observed line intensity Iul,obs is then Iul,obs = Ful,obsΩ, where Ful,obs is the sum of the integrated line flux in the extraction aperture with solid angle Ω equal to a circle with radius 1.″

The extinction correction needed for the line intensities can be derived using the H2 rotation diagrams and a wavelength dependent extinction law Aλ(AK). In Fig G.1 of Appendix G on Zenodo, we show the KP5 (Pontoppidan et al. 2024) extinction curve at AK values ranging from 1 to 10. This theoretical extinction curve constrained by observations is appropriate for dense molecular clouds, and contains absorption features from silicates and various ice species. It has been shown to perform well for extinction correction of H2 lines in JWST spectra (Narang et al. 2024). The J = 1–4 transitions are particularly useful as a measure of extinction, as the S(3) transition at ~9.66 µm is found within the broad 10 µm silicate feature, whereas the S(1), S(2), and S(4) transitions lie outside this feature and experience approximately the same extinction. Using the uncorrected intensities to produce a rotation diagram, the S(3) transition will thus appear relatively suppressed, as has been found in many previous H2 studies (e.g. Bertoldi et al. 1999; Neufeld et al. 1998). By fitting a single temperature component with a variable AK to the S(1) to S(4) transitions, the rotation diagram is effectively linearized. Specifically, a fit to the rotation diagram using the uncorrected intensities is done with the relation: yobs=ln(Nwarm exp(−x/Twarm)Q(Twarm )10−Aλ/2.5),${y_{{\rm{obs}}}} = \ln \left( {{N_{{\rm{warm }}}}{{\exp \left( { - x/{T_{{\rm{warm}}}}} \right)} \over {Q\left( {{T_{{\rm{warm }}}}} \right)}}{{10}^{ - {A_\lambda }/2.5}}} \right),$

where yobs=ln(4πIul,obsAulhv𝑔u).${y_{{\rm{obs}}}} = \ln \left( {{{4\pi {I_{ul,{\rm{obs}}}}} \over {{A_{ul}}hv{g_u}}}} \right).$

We performed this fit with equal weighting on all lines because of the importance of the relatively faint S(3) line in the determination of the extinction. We show the effects of this correction for HH 211 and IRAS 23385+6053 in Fig. G.2 of Appendix G on Zenodo.

After correction for extinction, our H2 rotation diagrams all show some amount of curvature in the higher upper energy S(5)- S(8) lines and are thus not well described by a single temperature component. Optically thick emission can produce curvature in a rotation diagram, but this should not be significant for H2 or HD rotational emission in protostellar outflows given the low column density of the gas and small Einstein coefficients. A more reasonable explanation is the presence of multiple temperature components in different shock environments in the aperture.

This also produces curvature in a rotation diagram, and has been noted in many previous studies of H2 emission from outflows (van Dishoeck 2004; Yuan & Neufeld 2011; Yuan et al. 2012). In addition to curvature, some of our H2 rotation diagrams show a zigzag pattern where the ortho (odd-J) transitions are suppressed relative to the para (even-J) transitions (See Fig. G.3 on Zenodo). This is characteristic of an ortho-to-para H2 column density ratio less than the typical value of 3 reached in local thermodynamic equilibrium (Neufeld et al. 2006b).

To fit the extinction corrected H2 rotation diagrams, we thus used a two-temperature model with an additional correction for the ortho-to-para correction term: y=ln(Nwarm exp(−x/Twarm )Q(Twarm )+Nhot exp(−x/Thot )Q(Thot ))+z(J),$y = \ln \left( {{N_{{\rm{warm }}}}{{\exp \left( { - x/{T_{{\rm{warm }}}}} \right)} \over {Q\left( {{T_{{\rm{warm }}}}} \right)}} + {N_{{\rm{hot }}}}{{\exp \left( { - x/{T_{{\rm{hot }}}}} \right)} \over {Q\left( {{T_{{\rm{hot }}}}} \right)}}} \right) + z(J),$

where Ntot = Nwarm + Nhot, for a warm component (Twarm = 250– 1000 K) and hot component (Thot = 1000–6000 K), and where z(J)={ ln(OPR3) odd J0, even J, $z(J) = \left\{ {\matrix{ {\ln \left( {{{{\rm{OPR}}} \over 3}} \right)} \hfill & {{\rm{ odd }}J} \hfill \cr {0,} \hfill & {{\rm{ even }}J} \hfill \cr } ,} \right.$

and OPR is the ortho-to-para ratio. Representative fits for IRAS 23385+6053 and HH 211 are shown in the top row of Fig. 4. We summarize the H2 fit results in Table C.1. We note that the error bars on the fit parameters include only the uncertainty propagated from the line flux measurements; we consider the effect of uncertainty in AK on deriving [D/H] in Sect. 3.4.

The majority of our H2 rotation diagrams are well described by two temperature components in LTE with OPR = 3. In some cases, OPR = 2–3 is found, but the values of the warm and hot components are still consistent with LTE. For a handful of positions associated with very bright bow-shocks, the rotation diagrams become nearly horizontal for the higher excitation transitions (S(4) to S(8)). Here, the temperature of the hot component is very poorly constrained and reaches the maximum allowed in the fit of 3000 K. This suggests a still hotter component of the H2 emission, which has been detected in rovibrational emission from other outflows with ISO-SWS and JWST/NIRspec (van Dishoeck et al. 1998; Giannini et al. 2008, 2015; Le Gouellec et al. 2024). In any case, the column density of the warm H2 component which we use in our estimate of [D/H] is determined almost entirely by the S(1) to S(4) lines, for which the temperature is well constrained.

The rotation diagram fits for HD are somewhat simpler, as we see little evidence of multiple temperature components, and HD does not have ortho and para forms. Since the HD is also only present in the outflow, it should experience the same extinction as H2, and we thus corrected the line intensities using the same value for AK, and performed a fit of a single temperature component to the rotation diagram. Representative fits for IRAS 23385+6053 and HH 211 are shown in the bottom row of Fig. 4. We note that the excitation energies of the HD R(4) to R(6) lines are similar to the H2 S(2) to S(4), so the single component fit to the HD transitions should trace gas with similar excitation conditions as the warm component in the H2 fits, which dominates Ntot . However, the much larger Einstein coefficients of the HD transitions than those of H2 result in a higher critical density, and deviations from LTE may become significant for the low densities in protostellar outflows. We discuss this further in Sect. 3.4, where we correct for non-LTE conditions in our determination of [D/H]. We note that we do not consider a temperature gradient as in previous works analyzing Spitzer detections of HD (Yuan & Neufeld 2011; Yuan et al. 2012), as the HD lines analyzed here are not sensitive to the cold (<300 K) gas that could be probed with the lower excitation HD lines detectable with Spitzer.

For our other JOYS targets where H2 is detected but no HD lines are observed, we can estimate upper limits on the HD column density. We converted each HD line upper limit to a column density assuming the temperature of the warm component of H2 in the rotation diagram fits, and applied the same extinction correction to the line intensities as was derived for H2 . We took the lowest column density upper limit of the three transitions as our overall 3σ upper limit for each observation. We note that despite the much lower sensitivity in channel 4 of MIRI/MRS, the most constraining upper limit is typically from the 23.04 µm R(4) line of HD, as this transition lies at the lowest upper energy level (see Table A.1). We summarize the results of our HD rotation diagram fits in Table C.1. Where we have detections of HD, 3–5 transitions are detected, which provides much more precise constraints on the column density and temperature of HD in protostellar outflows than possible with previous facilities. In general, we find the excitation temperatures of the single component for HD to be intermediate between the temperatures of the warm and hot H2 components, typically 500–1000 K, and the column density of HD to be ~5 orders of magnitude lower than that of the warm H2 component. We can thus rule out that enhancement of HD by orders of magnitude relative to H2 is occurring, as is often the case for other complex deuterated species (Ceccarelli et al. 2014).

thumbnail Fig. 4

Rotational diagrams for the detected H2 (left panels) and HD (right panels) transitions in IRAS 23385 + 6053 aperture 1 and HH 211 aperture 8. Corrections for extinction and an ortho-to-para ratio <3 have been applied (see text and App. G on Zenodo).

3.4 [D/H] abundance and upper limits

As a starting point to determine [D/H] abundances, we can assume LTE holds, and compare the column density from the warm component of H2 with that of HD traced the same approximate ranges in upper energy of the transitions. For our protostellar outflows, we assumed that the outflows are almost entirely molecular H2, except for a core of atomic H associated with the jet, as has been found for HH 211 (Caratti o Garatti et al. 2024). The [D/H] abundance or upper limit is then simply [D/H]=12NHD/NH2$[{\rm{D}}/{\rm{H}}] = {1 \over 2}{N_{{\rm{HD}}}}/{N_{{{\rm{H}}_2}}}$. However, at densities <1 × 106 cm−3 that can occur in protostellar outflows, non-LTE deviations in the level populations for the J > 4 HD lines are significant (Bertoldi et al. 1999). Prior estimates of the H2 from a non-LTE analysis of a handful of detected HD transitions toward the L1448 outflow suggest that densities <1 × 106 cm−3 are indeed present in some regions of the flow (Giannini et al. 2011). Bertoldi et al. (1999) show that not accounting for the non-LTE excitation can result in an overestimate of the HD column density. For their analysis of HD in the Orion OMC-1 outflow, they found that non-LTE excitation increases the inferred column density by a factor of ~1.47 relative to an LTE model, though the exact correction factor is sensitive to the shock density, temperature, and dissociation fraction.

Another important effect on the derived deuterium abundance considered by Bertoldi et al. (1999) is the chemical conversion of HD to D via: HD + H + ∆H0 ⇌ D + H2, where ∆H0 = 418 K is the enthalpy difference. This can be significant in the warmest and densest parts of C-type shocks present in some regions of protostellar outflows, though HD reforms in the cooler post-shock regions. For the Orion OMC-1 outflow, Bertoldi et al. (1999) found that the average effect over the entire shock region was to reduce the HD abundance relative to H2 by a factor of 1.67. When the effect of chemical conversion is combined with the non-LTE excitation, their true [D/H] abundance is a factor of 2.45 higher than what would be inferred from a purely LTE analysis. Interestingly, Bertoldi et al. (1999) find that this result is not particularly sensitive to how dissociative the shock is. In a less dissociative shock with a smaller fraction of atomic H, chemical conversion of HD to D is weaker but the correction for non-LTE effects is larger, and the combined correction factor is approximately constant. To account for the combined effects of chemical conversion and non-LTE excitation of HD, we thus followed Bertoldi et al. (1999) and multiplied our LTE [D/H] estimates by a factor of 2.45.

An estimate for the uncertainty for [D/H] should also include the effect of extinction. However, the uncertainty in the column densities of warm H2 and HD from Sect. 3.3 includes only the propagated errors on the line fluxes, as AK was held fixed to the value determined from the S(1) to S(4) transitions of H2 . To estimate the effect of AK uncertainty on the derived [D/H], we reran our two component fits with AK up to two magnitudes larger or smaller than the best fit value. As an example, we show this procedure for IRAS 23385+6053 in Fig G.4 in Appendix G on Zenodo. The resulting column densities for warm H2 and HD are shown as a function of AK in the upper panel, while the derived [D/H] is shown in the lower panel. Since both the HD and H2 column densities are strongly correlated with any change in AK, the overall effect on [D/H] is small. We took as our uncertainty in AK the result from the fit of a single temperature component to the S(1) to S(4) transition (green shaded region in lower panel of Fig. G.4), and converted this to an uncertainty in [D/H] (blue shaded region in lower panel of Fig. G.4). We included the uncertainty of AK it in our final estimate of [D/H] by adding it in quadrature with the uncertainty propagated from the fits to the rotation diagrams.

4 Results

We summarize our inferred deuterium abundances or upper limits in Table 2, and also list the Galactocentric radius RGC to each source. For apertures and sources with a non-detection of HD and low column densities and/or low temperatures of the warm H2 component, the upper limits on the HD column density result are low as the HD lines are not expected to be bright. This results in corresponding upper limits on [D/H] which are well above the primordial abundance of 2.58 ± 0.13 × 10−5 (Cyburt et al. 2016), and thus not constraining.

We also include in Table 2 literature [D/H] measurements that have been derived from rotational lines of HD and H2 observed in a protostellar outflow and photodissociation region (PDR) with ISO (Bertoldi et al. 1999; Wright et al. 1999). We note that similar rotation-line-based [D/H] measurements are available from Spitzer observations of two nearby proto- stellar outflows and a supernova remnant (Yuan et al. 2012); however, the Yuan et al. (2012) [D/H] values have systematic uncertainties of a factor of a few depending on whether spectroscopic or imaging data is used to estimate line fluxes, and thus are not particularly constraining. A detection of DCN in a cloud ∼0.1 pc from the Galactic center was reported by Lubowich et al. (2000), and through comparison with an HCN line and chemical fractionation modeling [D/H]= 1.7 ± 0.3 × 10−6 was claimed; however, this result may also vary by an order of magnitude or more given the systematic uncertainty associated with fractionation models and their strong temperature dependence. Finally, we also include the [D/H]= 2.0 ± 0.1 × 10−5 in the local Galactic disk derived through UV absorption lines of D I and H I (Prodanović et al. 2010), and the [D/H]=7.0 ± 2.0 × 10−6 derived from D I and O I lines and scaling by the ISM oxygen abundance (Hébrard 2010).

Fig. 5 presents our sources with detected HD lines and a [D/H] measurement (including the upward correction factor of 2.45 for non-LTE effects and chemical conversion, see Sect. 3.4), or a [D/H] upper limit below the primordial abundance. Significant variations in [D/H] abundance of a factor ∼4 are seen between detected sources. The outflows from low-mass protostars (green symbols in Fig. 5) originate from star-forming regions at similar positions in the Galactic disk. These regions should have formed from gas that has experienced similar levels of astration and thus with nearly the same [D/H] abundance (discussed further in Sect. 5). This should certainly be true when comparing different positions from the same outflow (for example, HH 211 and L1448-mm); however, a similar magnitude of variation as between different sources is seen instead. We note that in HH 211 [D/H] is robustly lower in the bow-shock positions (apertures 1–3) than the jet positions (apertures 4– 6), hinting at a link between the shock environment and [D/H]. Overall, it is thus unlikely that the [D/H] variations seen are the result of differences in the degree of astration for the low-mass sources. However, it is crucial to note that the [D/H] derived from observations, whether from rotational H2 and HD lines or H I and D I absorption lines, is only sensitive to the gasphase component. Galactic chemical evolution (GCE) models make predictions for the total [D/H], including both molecular and atomic forms and any deuterium sequestered in dust grains or deuterated PAHs.

We note that the factor of ∼4 abundance variation is similar to that seen by Linsky et al. (2006) and Prodanović et al. (2010) when comparing their UV absorption [D/H] measurements along different sight lines in the local Galactic disk, which they attribute to variations in the degree of deuterium depletion onto dust grains. It is thus plausible that the same depletion effect could explain the variations we see between different outflows and outflow positions, which we investigate further in the following Section.

Table 2

[D/H] abundances from mid-IR rotational lines.

thumbnail Fig. 5

Comparison of [D/H] abundance (Table 2) in this work. A correction factor of 2.54 has been applied to our [D/H] measurements (see text). Orange symbols show [D/H] in apertures for high-mass sources, while green symbols show [D/H] in apertures from low-mass sources. The dotted gray lines separate different sources. The blue shaded region shows the primordial total [D/H] abundance of Cyburt et al. (2016). [D/H] upper limits above the primordial abundance are omitted.

4.1 Line flux and rotation diagram correlations

As noted in Sect. 3.1, the integrated line intensity maps show that the morphology of the HD emission is clearly associated with the jets and bow-shocks in our sources. The HD emission is particularly bright toward outflow positions that are also bright in [Fe I] and [S I]. iron and sulfur are both known to be depleted onto dust grains (Jenkins 2009), though the exact reservoir of the latter is not well-understood (Anderson et al. 2013; Kama et al. 2019). A possible explanation for the association between HD, [Fe I], and [S I] is that deuterium is significantly depleted onto dust grains, and is subsequently released along with Fe and S in strong shocks where the dust grains are destroyed. Free atomic D can be readily converted to HD via D + H2 → HD + H.

To perform a more quantitative comparison, we can investigate correlations between the extracted line fluxes and the rotation diagram fit parameters in all of our apertures. We first apply the same extinction correction derived from the H2 rotation diagrams to all of our line fluxes using the method described in Sect. 3.3. We note that since the emission is extended and should fill the extracted apertures, the line flux and intensity are different only by a constant factor of the aperture solid angle. We show in the top row of Fig. 6 correlations between the HD R(6) line and 3 other lines associated with the outflow jets or bowshocks: [S I], [Fe I], [Fe II]. A particularly strong correlation (ρ = 0.95) is seen between the HD R(6) line and [S I], while the correlation with [Fe I] and [Fe II] exhibit much more scatter (ρ = 0.36 and ρ = 0.57 respectively). This is possibly due to variations in the ionization conditions between different shocks that affect the intensity of the [Fe I] and [Fe II] lines, but which are not as significant for [S I]. Namely, the ionization potential of the [Fe I] line is only 7.9 eV while for [S I] itis 10.36 eV. In the bottom row of Fig. 6, we show correlations between the HD R(6) line and the H2 S(1), S(4), and S(7) lines. Strong correlations are seen particularly with the higher upper-energy S(7) line (Eup ∼ 7200 K), suggesting the HD R(6) emission is associated with the denser and higher temperature components of the outflow.

While the association of the HD emission with dense and shocked regions of the outflows is clear, stronger evidence of D depletion and release in shocks would be a correlation between [D/H] and the gas phase abundance of refractory species or some measure of the degree of dust grain destruction. Determining the degree of dust destruction in the outflow shocks is nontrivial, however. For the case of HH211, modeling of the gas-phase abundance of Fe toward bright H2 knots in the inner jet and bow shock has been performed by Caratti o Garatti et al. (2024) using the models of Hollenbach & McKee (1989). Their results show that the gas-phase iron abundance in all shock knots are depleted with respect to the solar value, and the knots in the inner jet are less depleted than those in the bow shock. (Caratti o Garatti et al. 2015, Fig. 13). This is the opposite of the expectation that the gas-phase iron abundance should be higher further along the flow as more grain destruction occurs and more Fe is released (Nisini et al. 2005; Podio et al. 2006). However, it is consistent with the increased [D/H] in the inner jet (apertures 4–6) relative to the bow-shocks in HH 211 (apertures 1–3), suggesting depletion of deuterium into dust grains is indeed a plausible mechanism. We note that two important caveats in the modeling of Caratti o Garatti et al. (2024) are the uncertainty in the shock velocity and density, which makes the relative Fe depletion somewhat uncertain, and the assumption that sulfur is undepleted in the line ratios used, which affects the absolute gasphase abundance. To confirm the link between dust depletion and [D/H] variations in the outflow, further modeling of the shock conditions and degree of refractory depletion is needed.

thumbnail Fig. 6

Correlations between the extinction corrected flux of the HD R(6) line measured in all apertures for our sample with the extinction corrected flux of various atomic and H2 lines. The line flux is related to the line intensity by a constant factor of the aperture solid angle (see text). Top row: correlation with the [S I], [Fe I], and [Fe II] line flux. Bottom row: correlation with the H2 S(1), S(4), and S(7) line flux.

4.2 [D/H] variations and Galactic chemical evolution

Galactic chemical evolution models make predictions for the total [D/H] abundance throughout the Galactic disk that can in principle be compared with observations. This is useful for both placing constraints on the GCE models themselves and assessing the potential effects of deuterium depletion on the observed [D/H]. We thus plot in Fig. 7 as a function of RGC the average [D/H] of our low-mass sources, the [D/H] measurement or upper limits for our high-mass sources, [D/H] from previous ISO measurements, and [D/H] measurements based on UV absorption lines. We overlay the primordial deuterium abundance of Cyburt et al. (2016), as well as total [D/H] gradients from GCE models selected from Romano et al. (2006) and van de Voort et al. (2018).

If we interpret our [D/H] measurements derived from HD and H2 as representative of total [D/H], their values are generally much more similar to the low estimates of the local disk [D/H] of Hébrard & Moos (2003) and Hébrard (2010) than the high estimates of Linsky et al. (2006) and Prodanović et al. (2010). As noted in Sect. 1, the debate on which value represents the present Galactic disk [D/H] abundance within a few kiloparsec of the Sun remains unsettled, and the significance of the inferred depletion onto dust grains and the uncertainties in abundances and GCE effects introduced by using [D/O] as a proxy are both debated (Tosi 2010). However, the molecular line measurements and the low local disk [D/H] are typically a factor of ∼2–4 below the predictions of the GCE models. The model of Romano et al. (2006) is based on the GCE code of Chiappini et al. (2001), and simultaneously satisfies existing chemical abundance constraints while being marginally able to reproduce the depletion-corrected [D/H] of Prodanović et al. (2010). Romano et al. (2006) have explored whether GCE models can reproduce the low values of [D/H] in the local Galactic disk of <1.0 × 10−5 from Hébrard (2010) by increasing the efficiency of star formation in the solar neighborhood, and thus the amount of astration. However, the increased star formation also produces more enrichment of the local Galactic disk with metals, and a resulting metallicity distribution for G and K type stars in the solar neighbourhood with much larger [Fe/H] than is observed. Similar problems in fitting this metallicity distribution occur if the amount of infall onto the Galactic disk of deuterium-rich and metal-poor gas is instead decreased.

It is also interesting to compare the [D/H] predictions from GCE models which are not aimed at specifically reproducing the abundance constraints from our Galaxy. In the context of the evolution of the [D/H] abundance in galaxies with cosmic time, van de Voort et al. (2018) have performed GCE simulations based on cosmological zoom-in simulations for galaxies forming over a wide range of redshifts and metallicity. The selected model shown in Fig. 7 shows the total [D/H] gradient for a typical milky way-like galaxy. The high [D/H] across the disk is much more consistent with the local disk [D/H] of Prodanović et al. (2010) that includes a correction for dust depletion.

GCE models generally predict an increasing total [D/H] with galactocentric radius owing to the increased star formation, and thus astration toward the Galactic center. Our high-mass protostar targets with detections or constraining upper limits on [D/H] are located in the inner and outer Galaxy, and may therefore be useful for measuring this predicted gradient, with the caveat of possible dust depletion. We detect HD lines only in the local low- mass protostars and in IRAS 23385+6053 at RGC=11 kpc, but not in any of our inner Galaxy targets. Taken at face value, the upper limit for IRAS 18089-1732 in particular would be consistent with a [D/H] gradient; however, we cannot rule out the possibility that deuterium is more depleted onto dust in the inner Galaxy than in the local disk or outer Galaxy.

thumbnail Fig. 7

[D/H] as a function of Galactocentric radius based on our HD and H2 rotational line measurements within JOYS (squares), other rotational line values from ISO (circles) (Wright et al. 1999; Bertoldi et al. 1999), and local Galactic disk values from UV absorption lines (triangles) (Prodanović et al. 2010; Hébrard 2010). Our measurements and those for the Orion OMC outflow include an upward correction factor of 2.54 to [D/H] for chemical conversion of HD to D. We also show Galactic chemical evolution models of our Galaxy from (black dashed line, reproduced from solid line in Fig. 5 of Romano et al. 2006) and from a Milky Way type Galaxy in a cosmological zoom-in simulation from (gray line, reproduced from model “m12i” in Fig. 5 of van de Voort et al. 2018).

5 Discussion

In summary, our results show large variations in [D/H] between our local low-mass sources, and the derived values of [D/H] are generally a factor of ∼2–4 lower with respect to the predictions of Galactic chemical evolution models. Several of our low-mass sources also show a [D/H] at nearly the primordial value (Fig. 5). If depletion is the underlying cause for the variations in [D/H], this would suggest that the total [D/H] is underestimated in many of our objects when derived from the molecular lines alone. Variable dust depletion would resolve the tension with the various GCE models shown in Fig. 7, which are generally unable to produce a low [D/H] and satisfy other abundance constraints. However, we do not see any clear trends in [D/H] abundance in comparison with the fluxes of the [S I] and [Fe I] lines, which could trace refractory material in dust grains released in shocks. One possible reason is that deuterium is thought to be depleted by replacing H-C bonds in carbonaceous dust grains or deuterated PAHs, with the latter experiencing different shock conditions, and thus the atomic [S I] and [Fe I] lines are not tracing the solid reservoir into which deuterium is depleted.

An additional source of uncertainty in our [D/H] measurements is the value of the chemical conversion and non-LTE correction factors derived by Bertoldi et al. (1999). While they show that the combination of these two correction factors is approximately constant across a range of shock conditions, more detailed modeling of individual spectra and shocks is warranted to see what fraction of our variations in [D/H] can be attributed to differences in the excitation conditions and chemical conversion.

To make progress on measuring total [D/H] with JWST, further work is needed on both the observations and modeling. On the observational front, deeper integrations would provide much better upper limits on the HD column density in quiescent regions of outflows where dust destruction is unlikely to occur. Observations of additional distant targets in the inner and outer Galaxy are needed, especially close to the Galactic center, if GCE models are to be constrained. In addition to pro- tostellar outflows, rotational HD lines have also been detected in the IC 443C supernova remnant (Yuan et al. 2012). Observations of HD emission in supernova remnants may provide a useful probe of [D/H] in higher velocity shock environments where grain destruction may be more complete. Observations with MIRI/MRS cover many rotational lines of H2 and HD, but rovibrational transitions of both are found in the near-infrared accessible with NIRSpec. The additional rovibrational transitions would be useful for determining the shock conditions. Rovibrational transitions of HD have already been detected in PDRs (Peeters et al. 2024) with NIRspec. Hints of Infrared emission bands from deuterated PAHs have been detected with JWST (Peeters et al. 2004, 2024), and their interpretation will be important for understanding their role as a reservoir for deuterium depletion.

In terms of modeling, it is particularly important to understand the shock conditions and degree of dust destruction. Shock models for H2 (e.g. Kristensen et al. 2023) and forbidden atomic line emission (e.g. Hollenbach & McKee 1989) can be used to constrain the shock density, and velocity, and excitation conditions. (Caratti o Garatti et al. 2015, 2024). These are important inputs for determining the degree of non-LTE excitation and chemical conversion of HD that takes place in the shock positions.

6 Summary and conclusions

Using JWST/MIRI observations from the JOYS program (PIDs: 1290, 1257) we have presented an analysis of rotational H2 and HD lines in outflows from 5 high-mass and 5 low-mass protostars with the goal of constraining the [D/H] abundance. Our results and conclusions can be summarized as follows:

  • H2 lines are detected toward all selected sources, but HD lines are only detected toward the brightest positions tracing the outflow shock knots in the low-mass sources, and in the one outer Galaxy (RGC = 11 kpc) high-mass protostar IRAS 23385+6053. Notably, no HD emission is detected in the inner Galaxy high-mass protostars (RGC ∼4 kpc) where a higher degree of deuterium astration is expected.

  • In resolved low-mass outflows, the morphology and flux of the HD emission is strongly correlated with high-excitation H2 transitions, [S I], and [Fe I] line emission, showing a clear association with bright shock knots within the jet and bowshocks.

  • The column density or upper limits of H2 and HD toward these sources has been obtained through an LTE rotation diagram analysis. Additional corrections to the HD column density accounting for non-LTE behavior and chemical conversion have been applied following Bertoldi et al. (1999). JWST/MIRI can now provide much more precise constraints on the H2 and HD column density than previous facilities, thanks to the high sensitivity and large number of H2 and HD lines in the mid-IR, which allow the excitation conditions to be more robustly constrained than previously possible.

  • The H2 rotation diagrams typically shows two components: a warm (∼400–900 K) and high column density (∼1020 – 1022 cm−2) component associated with the wide- angle outflows in the integrated line intensity maps, and a hot (∼1000–2000 K) and much lower column density component (∼1017 – 1019 cm−2) associated with bright knots and bowshocks. In most cases the H2 shows and ortho-to-para ratio of ∼3. For a handful of positions, an ortho-to-para ratio <2.5 and evidence of enhanced excitation for the higher energy H2 transitions is seen.

  • The HD rotation diagrams are well-described by a single temperature component, and typically show column densities ∼105 orders magnitude smaller than the warm H2 component, and excitation temperatures intermediate between the warm and hot H2 (∼700–1000 K.)

  • Our derived [D/H] abundances show a factor of ∼4 variation between low-mass protostar sources and between different positions within the same outflows. However, the gas compared between low-mass protostars and within a given outflow is expected to have experienced the same degree of astration. Furthermore, our abundances are a factor of ∼2–4 below the predictions of GCE models. This suggests that we are not probing the total [D/H] in the observed protostars.

  • The observed variations may be caused by the significant depletion of D into carbonaceous dust grains and varying degrees of grain destruction. The enhanced [D/H] and gasphase abundance of Fe in the inner jet relative to the bowshocks in HH 211 is consistent with this scenario.

  • Most of our measurements of [D/H] or similar literature measurements using rotational lines suggest [D/H] ≲ 1.0 × 10−5 over Galactocentric radii from 4.6 to 11 kpc, which is very difficult to reconcile with existing Galactic chemical evolution models while still satisfying other chemical abundance constraints. Depletion of D into dust and its release in shocks would be consistent with the variations in [D/H] between low-mass sources and the typically low values of [D/H] with respect to Galactic chemical evolution models. Further modeling of the effects of the grain destruction in the shocks are needed.

If the total [D/H] abundance can be accurately estimated, the detection of HD in IRAS 23385+6053 at a galactocentric distance of 11 kpc, and the upper limits in the inner Galaxy from IRAS 18089-1732 present an exciting opportunity. Additional observations of HD and H2 in protostellar or supernova shocks when combined with models for the effects of depletion and chemical conversion could for the first time sample [D/H] across the Galactic disk, especially closer to the Galactic center and the outer Galaxy where little is known. This would provide powerful new constraints on models for the formation and chemical evolution of our Galaxy.

Data availability

Appendices D–G are available on Zenodo: https://zenodo.org/records/14592550

Acknowledgements

We are grateful to John H. Black for useful discussions on modeling of H2 and HD emission, and to the anonymous referee for their comments on the paper. This work is based on observations made with the NASA/ESA/CSA James Webb Space Telescope. The data were obtained from the Mikulski Archive for Space Telescopes at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-03127 for JWST. These observations are associated with program 1290. Astrochemistry in Leiden is supported by funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No. 101019751 MOLDISK), by the Netherlands Research School for Astronomy (NOVA), and by grant TOP-1 614.001.751 from the Dutch Research Council (NWO). A.C.G. acknowledges support from PRIN-MUR 2022 20228JPA3A “The path to star and planet formation in the JWST era (PATH)” funded by NextGeneration EU and by INAF-GoG 2022 “NIR-dark Accretion Outbursts in Massive Young stellar objects (NAOMY)” and Large Grant INAF 2022 “YSOs Outflows, Disks and Accretion: toward a global framework for the evolution of planet forming systems (YODA)”. P.J.K. acknowledge support from the Science Foundation Ireland/Irish Research Council Pathway program under grant No. 21/PATH-S/9360.

Appendix A Analyzed Spectral Lines

We report the H2, HD, and forbidden atomic lines analyzed in this work in Table A.1.

Table A.1

Analyzed H2, HD, and forbidden atomic lines.

Appendix B Aperture Locations

We describe the position of the apertures used for spectral extraction toward each source in Table B.1.

Table B.1

Position of our extracted aperture centers in each source. Each aperture is 1″ in radius.

Appendix C Rotation Diagram Fit Results

In table C.1, we provide the results of the rotation diagram fits for H2 and HD to each aperture described in Section 3.3.

Table C.1

Rotation Diagram Fit Results.

References

  1. Anderson, D. E., Bergin, E. A., Maret, S., & Wakelam, V. 2013, ApJ, 779, 141 [Google Scholar]
  2. Awad, Z., Viti, S., Bayet, E., & Caselli, P. 2014, MNRAS, 443, 275 [NASA ADS] [CrossRef] [Google Scholar]
  3. Bergin, E. A., Cleeves, L. I., Crockett, N. R., Favre, C., & Neill, J. L. 2013, in Astronomical Society of the Pacific Conference Series, 476, New Trends in Radio Astronomy in the ALMA Era: The 30th Anniversary of Nobeyama Radio Observatory, eds. R. Kawabe, N. Kuno, & S. Yamamoto, 185 [NASA ADS] [Google Scholar]
  4. Bertoldi, F., Timmermann, R., Rosenthal, D., Drapatz, S., & Wright, C. M. 1999, A&A, 346, 267 [NASA ADS] [Google Scholar]
  5. Beuther, H., van Dishoeck, E. F., Tychoniec, L., et al. 2023, A&A, 673, A121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Bushouse, H., Eisenhamer, J., Dencheva, N., et al. 2024, https://doi.org/18.5281/zenodo.18878758 [Google Scholar]
  7. Caratti o Garatti, A., Stecklum, B., Linz, H., Garcia Lopez, R., & Sanna, A. 2015, A&A, 573, A82 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Caratti o Garatti, A., Ray, T. P., Kavanagh, P. J., et al. 2024, A&A, 691, A134 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Ceccarelli, C., Caselli, P., Bockelée-Morvan, D., et al. 2014, in Protostars and Planets VI, eds. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 859 [Google Scholar]
  10. Chiappini, C., Matteucci, F., & Romano, D. 2001, ApJ, 554, 1044 [NASA ADS] [CrossRef] [Google Scholar]
  11. Chiappini, C., Renda, A., & Matteucci, F. 2002, A&A, 395, 789 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Christiaens, V., Gonzalez, C., Farkas, R., et al. 2023, J. Open Source Softw., 8, 4774 [NASA ADS] [CrossRef] [Google Scholar]
  13. Cooke, R. J., Pettini, M., Jorgenson, R. A., Murphy, M. T., & Steidel, C. C. 2014, ApJ, 781, 31 [NASA ADS] [CrossRef] [Google Scholar]
  14. Cyburt, R. H., Fields, B. D., Olive, K. A., & Yeh, T.-H. 2016, Rev. Mod. Phys., 88, 015004 [NASA ADS] [CrossRef] [Google Scholar]
  15. Draine, B. T. 2006, in Astronomical Society of the Pacific Conference Series, 348, Astrophysics in the Far Ultraviolet: Five Years of Discovery with FUSE, eds. G. Sonneborn, H. W. Moos, & B. G. Andersson, 58 [Google Scholar]
  16. Epstein, R. I., Lattimer, J. M., & Schramm, D. N. 1976, Nature, 263, 198 [NASA ADS] [CrossRef] [Google Scholar]
  17. Friedman, S. D., Chayer, P., Jenkins, E. B., et al. 2023, ApJ, 946, 34 [NASA ADS] [CrossRef] [Google Scholar]
  18. Giannini, T., Calzoletti, L., Nisini, B., et al. 2008, A&A, 481, 123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Giannini, T., Nisini, B., Neufeld, D., et al. 2011, in IAU Symposium, 280, The Molecular Universe, eds. J. Cernicharo, & R. Bachiller, 329 [Google Scholar]
  20. Giannini, T., Antoniucci, S., Nisini, B., Bacciotti, F., & Podio, L. 2015, ApJ, 814, 52 [NASA ADS] [CrossRef] [Google Scholar]
  21. Gieser, C., Beuther, H., van Dishoeck, E. F., et al. 2023, A&A, 679, A108 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Greenfield, P., & Miller, T. 2016, Astron. Comput., 16, 41 [NASA ADS] [Google Scholar]
  23. Hébrard, G. 2010, in Light Elements in the Universe, 268, eds. C. Charbonnel, M. Tosi, F. Primas, & C. Chiappini, 59 [Google Scholar]
  24. Hébrard, G., & Moos, H. W. 2003, ApJ, 599, 297 [Google Scholar]
  25. Hollenbach, D., & McKee, C. F. 1989, ApJ, 342, 306 [Google Scholar]
  26. Hoopes, C. G., Sembach, K. R., Hébrard, G., Moos, H. W., & Knauth, D. C. 2003, ApJ, 586, 1094 [NASA ADS] [CrossRef] [Google Scholar]
  27. Jenkins, E. B. 2009, ApJ, 700, 1299 [Google Scholar]
  28. Jura, M. 1982, in NASA Conference Publication, 2238, ed. Y. Kondo, 54 [NASA ADS] [Google Scholar]
  29. Kama, M., Shorttle, O., Jermyn, A. S., et al. 2019, ApJ, 885, 114 [Google Scholar]
  30. Kristensen, L. E., Godard, B., Guillard, P., Gusdorf, A., & Pineau des Forêts, G. 2023, A&A, 675, A86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Le Gouellec, V. J. M., Lew, B. W. P., Greene, T. P., et al. 2024, ApJ, [arXiv:2418.11895] [Google Scholar]
  32. Linsky, J. L., Draine, B. T., Moos, H. W., et al. 2006, ApJ, 647, 1106 [Google Scholar]
  33. Lubowich, D. A., Pasachoff, J. M., Balonek, T. J., et al. 2000, Nature, 405, 1025 [NASA ADS] [CrossRef] [Google Scholar]
  34. Mangum, J. G., & Shirley, Y. L. 2015, PASP, 127, 266 [Google Scholar]
  35. McClure, M. K., Bergin, E. A., Cleeves, L. I., et al. 2016, ApJ, 831, 167 [Google Scholar]
  36. Narang, M., Manoj, P., Tyagi, H., et al. 2024, ApJ, 962, L16 [NASA ADS] [CrossRef] [Google Scholar]
  37. Neufeld, D. A., Melnick, G. J., & Harwit, M. 1998, ApJ, 506, L75 [NASA ADS] [CrossRef] [Google Scholar]
  38. Neufeld, D. A., Green, J. D., Hollenbach, D. J., et al. 2006a, ApJ, 647, L33 [NASA ADS] [CrossRef] [Google Scholar]
  39. Neufeld, D. A., Melnick, G. J., Sonnentrucker, P., et al. 2006b, ApJ, 649, 816 [NASA ADS] [CrossRef] [Google Scholar]
  40. Nisini, B., Bacciotti, F., Giannini, T., et al. 2005, A&A, 441, 159 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Peeters, E., Allamandola, L. J., Bauschlicher, C. W., J., et al. 2004, ApJ, 604, 252 [NASA ADS] [CrossRef] [Google Scholar]
  42. Peeters, E., Habart, E., Berné, O., et al. 2024, A&A, 685, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Penzias, A. A. 1979, ApJ, 228, 430 [NASA ADS] [CrossRef] [Google Scholar]
  44. Penzias, A. A., Wannier, P. G., Wilson, R. H., & Linke, R. A. 1977, ApJ, 211, 108 [NASA ADS] [CrossRef] [Google Scholar]
  45. Podio, L., Bacciotti, F., Nisini, B., et al. 2006, A&A, 456, 189 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  46. Polehampton, E. T., Baluteau, J. P., Ceccarelli, C., Swinyard, B. M., & Caux, E. 2002, A&A, 388, L44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Pontoppidan, K. M., Evans, N., Bergner, J., & Yang, Y.-L. 2024, RNAAS, 8, 68 [NASA ADS] [Google Scholar]
  48. Prantzos, N. 1996, A&A, 310, 106 [NASA ADS] [Google Scholar]
  49. Prodanović, T., Steigman, G., & Fields, B. D. 2010, MNRAS, 406, 1108 [NASA ADS] [Google Scholar]
  50. Romano, D., Tosi, M., Matteucci, F., & Chiappini, C. 2003, MNRAS, 346, 295 [NASA ADS] [CrossRef] [Google Scholar]
  51. Romano, D., Tosi, M., Chiappini, C., & Matteucci, F. 2006, MNRAS, 369, 295 [NASA ADS] [CrossRef] [Google Scholar]
  52. Sarkar, S. 1996, Rep. Progr. Phys., 59, 1493 [CrossRef] [Google Scholar]
  53. Tielens, A. G. G. M. 1983, A&A, 119, 177 [NASA ADS] [Google Scholar]
  54. Tielens, A. G. G. M. 2008, ARA&A, 46, 289 [NASA ADS] [CrossRef] [Google Scholar]
  55. Tielens, A. G. G. M. 2013, Rev. Mod. Phys., 85, 1021 [NASA ADS] [CrossRef] [Google Scholar]
  56. Tosi, M. 2010, in Light Elements in the Universe, 268, eds. C. Charbonnel, M. Tosi, F. Primas, & C. Chiappini, 153 [NASA ADS] [Google Scholar]
  57. Tsujimoto, T. 2011, MNRAS, 410, 2540 [NASA ADS] [CrossRef] [Google Scholar]
  58. Tsujimoto, T., Bland-Hawthorn, J., & Freeman, K. C. 2010, PASJ, 62, 447 [CrossRef] [Google Scholar]
  59. van de Voort, F., Quataert, E., Faucher-Giguère, C.-A., et al. 2018, MNRAS, 477, 80 [CrossRef] [Google Scholar]
  60. van Dishoeck, E. F. 2004, ARA&A, 42, 119 [NASA ADS] [CrossRef] [Google Scholar]
  61. van Dishoeck, E. F., Wright, C. M., Cernicharo, J., et al. 1998, ApJ, 502, L173 [NASA ADS] [CrossRef] [Google Scholar]
  62. Wootten, A. 1987, in Astrochemistry, 120, eds. M. S. Vardya, & S. P. Tarafdar, 311 [NASA ADS] [CrossRef] [Google Scholar]
  63. Wright, C. M., van Dishoeck, E. F., Cox, P., Sidher, S. D., & Kessler, M. F. 1999, ApJ, 515, L29 [NASA ADS] [CrossRef] [Google Scholar]
  64. Yuan, Y., & Neufeld, D. A. 2011, ApJ, 726, 76 [NASA ADS] [CrossRef] [Google Scholar]
  65. Yuan, Y., Neufeld, D. A., Sonnentrucker, P., Melnick, G. J., & Watson, D. M. 2012, ApJ, 753, 126 [NASA ADS] [CrossRef] [Google Scholar]

All Tables

Table 1

Sources analyzed for JWST MIRI/MRS H2 and HD observations.

In the text
Table 2

[D/H] abundances from mid-IR rotational lines.

In the text
Table A.1

Analyzed H2, HD, and forbidden atomic lines.

In the text
Table B.1

Position of our extracted aperture centers in each source. Each aperture is 1″ in radius.

In the text
Table C.1

Rotation Diagram Fit Results.

In the text

All Figures

thumbnail Fig. 1

Integrated line intensity maps for HH 211 of various lines and the continuum at 17 µm shown with a logarithmic stretch. The maps have been smoothed to a common resolution of 1″, shown by the white circle in the bottom-left. Apertures used for spectral extraction are shown by the green circles. An index for each aperture is provided in the top-left panel. The coordinates of the aperture centers can be found in Table B.1 of Appendix B.
In the text
thumbnail Fig. 2

Same as Fig. 1, but for IRAS 23385+6053.
In the text
thumbnail Fig. 3

Spectral line fits for H2, HD, and the [S I], [Fe I], and [Fe II] atomic lines for HH 211 aperture 8. The line flux (shaded blue region) is determined by the simultaneous fit of a Gaussian profile and a first order polynomial to represent the continuum. the rest-frame wavelength of each spectral line is shown by a dashed gray line, and non-detections are shown without a shaded region for the fit. Nearby molecular lines of are labeled.
In the text
thumbnail Fig. 4

Rotational diagrams for the detected H2 (left panels) and HD (right panels) transitions in IRAS 23385 + 6053 aperture 1 and HH 211 aperture 8. Corrections for extinction and an ortho-to-para ratio <3 have been applied (see text and App. G on Zenodo).
In the text
thumbnail Fig. 5

Comparison of [D/H] abundance (Table 2) in this work. A correction factor of 2.54 has been applied to our [D/H] measurements (see text). Orange symbols show [D/H] in apertures for high-mass sources, while green symbols show [D/H] in apertures from low-mass sources. The dotted gray lines separate different sources. The blue shaded region shows the primordial total [D/H] abundance of Cyburt et al. (2016). [D/H] upper limits above the primordial abundance are omitted.
In the text
thumbnail Fig. 6

Correlations between the extinction corrected flux of the HD R(6) line measured in all apertures for our sample with the extinction corrected flux of various atomic and H2 lines. The line flux is related to the line intensity by a constant factor of the aperture solid angle (see text). Top row: correlation with the [S I], [Fe I], and [Fe II] line flux. Bottom row: correlation with the H2 S(1), S(4), and S(7) line flux.
In the text
thumbnail Fig. 7

[D/H] as a function of Galactocentric radius based on our HD and H2 rotational line measurements within JOYS (squares), other rotational line values from ISO (circles) (Wright et al. 1999; Bertoldi et al. 1999), and local Galactic disk values from UV absorption lines (triangles) (Prodanović et al. 2010; Hébrard 2010). Our measurements and those for the Orion OMC outflow include an upward correction factor of 2.54 to [D/H] for chemical conversion of HD to D. We also show Galactic chemical evolution models of our Galaxy from (black dashed line, reproduced from solid line in Fig. 5 of Romano et al. 2006) and from a Milky Way type Galaxy in a cosmological zoom-in simulation from (gray line, reproduced from model “m12i” in Fig. 5 of van de Voort et al. 2018).
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.