© ESO, 2015

1. Introduction

The physical and chemical conditions present in low- and high-mass star-forming regions differ significantly. Massive star-forming regions are found to have higher UV radiation fields and levels of turbulence than their low-mass counterparts (see Stäuber et al. 2007; Herpin et al. 2012). The temperatures, feedback mechanisms, magnetic fields, accretion rates, and outflow forces are different between low- and high-mass young stellar objects (for more details see Bontemps et al. 1996; Behrend & Maeder 2001; Beuther et al. 2002, 2007; Palla & Stahler 1993; Cesaroni et al. 2005; Zinnecker & Yorke 2007).

However, many studies have shown that high-mass YSOs behave in certain aspects as scaled-up versions of their low-mass counterparts (van der Tak et al. 1999, 2000; Shepherd 2003; Johnston et al. 2012; San José-García et al. 2013, 2015; Karska et al. 2014a). In addition, the lifetime of the embedded phase of high-mass YSOs (0.07−0.4 Myr, Mottram et al. 2011) is comparable to that of low-mass YSOs (0.15 Myr for Class 0, 0.5 Myr for Class 0+I, Dunham et al. 2014), even if massive objects evolve more in this period. The line luminosity of molecules like CO, HCO⁺, and OH scales with bolometric luminosity and envelope mass, as well as the degree of turbulence in the warmer inner regions of protostellar envelopes (San José-García et al. 2013; van der Tak et al. 2013; Wampfler et al. 2013; Benz et al. 2015). Moreover, the kinematics of outflows and envelopes seem to be linked independently of the mass of the forming star (San José-García et al. 2013). Therefore, a further characterisation of the physical conditions and dynamics of these components will help to identify the differences and similarities between low- and high-mass YSOs and better understand the fundamental processes in the formation of stars.

Water is unquestionably a key molecule for studying the energetics and dynamical properties of protostellar environments (van Dishoeck et al. 2011). In particular, the analysis of the velocity-resolved water data provided by the Heterodyne Instrument for the Far-Infrared (HIFI; de Graauw et al. 2010) on board of Herschel Space Observatory (Pilbratt et al. 2010) allows us to characterise the emission from molecular outflows, which play a crucial role in the formation of young stars and in the feedback on their surroundings (Kristensen et al. 2012; van der Tak et al. 2013; Mottram et al. 2014). Given that the bulk of the water data on extragalactic sources out to the highest redshifts and the data provided by the other Herschel instruments (the Photodetector Array Camera and Spectrometer, PACS, Poglitsch et al. 2010; and the the Spectral and Photometric Imaging Receiver, SPIRE, Griffin et al. 2010) on galactic sources are spectrally unresolved, it is important to quantify the different components that make up the observed lines.

Outflows remove angular momentum effectively, which is necessary for the formation of a disk and mass accretion onto the forming star (see review by Lada & Kylafis 1999). The power agent of these structures (either jets or winds from the star/disk system) triggers not only the formation of the outflows, but also extreme and complex physical and chemical processes across the protostellar environment. In particular, different types of shocks take place in the outflow cavity wall at the interface between the cavity and the envelope. Non-dissociative outflow-cavity shocks are localised along the outflow cavity wall (Visser et al. 2012; Kristensen et al. 2012; Mottram et al. 2014). On the other hand, dissociative shocks take place either along the jet, revealed through perturbations known as extremely high velocity components (EHV; Bachiller et al. 1990; Tafalla et al. 2010), or at the base of the outflow cavity wall where the jet or wind impacts directly (Kristensen et al. 2012, 2013; Mottram et al. 2014). These shocks are also called spot shocks. Therefore, shocks and turbulent motions injected into the cavity wall propagating within this physical structure are products linked to the activity of the molecular outflow (Arce et al. 2007). The dynamical nature of these two phenomena (turbulence and shocks) is different, and they also differ from that characterising the entraining gas in classical outflows. To comprehensibly interpret the molecular emission of H₂O and ¹²CO within the outflow-cavity system, it is important to investigate whether the dynamical properties of low-mass objects can be extrapolated to more massive YSOs.

This was one of the goals of the “Water In Star-forming regions with Herschel” key programme (WISH; van Dishoeck et al. 2011), which observed several water lines, as well as high-J CO and isotopologue transitions, for a large sample of YSOs covering early evolutionary stages over a wide range of luminosities. An extensive study has been performed on all HIFI water observations for low-mass protostars (Mottram et al. 2014) and for low-lying water transitions within the high-mass sub-sample (van der Tak et al. 2013). The line profiles of the water transitions were analysed and decomposed into multiple velocity components, which are associated to different physical structures of the protostellar system. These studies investigated trends with luminosity, mass, and evolution and explored the dynamical and excitation conditions probed by these lines. In addition, observations with PACS reveal that ¹²CO J> 20 transitions originate mostly in shocks for both low- and high-mass YSOs (Herczeg et al. 2012; Manoj et al. 2013; Green et al. 2013; Karska et al. 2013, 2014a,b). The excitation of warm CO is also similar across the luminosity range, but rotational temperatures in high-mass objects are higher than in their low-mass counterparts in the case of H₂O, due to their higher densities (Karska et al. 2014a).

In order to link these two studies, this paper focuses on the analysis of the excited water lines across the entire WISH sample of YSOs. The ground-state water transitions of high-mass sources show absorption features from foreground clouds which complicates the extraction of velocity information from these lines (van der Tak et al. 2013), a reason why these lines are not considered in this study. Results from the line profile, line luminosity and excitation condition analysis are connected from low- to high-mass YSOs and interpreted together with those obtained from high-J¹²CO observations (J ≥ 10). In addition, the obtained results may help to interpret and understand those of extragalactic sources. The aim is to better constrain the dynamical properties of molecular outflows across a wide range of luminosities and complement the study presented in San José-García et al. (2013) based on the analysis of high-J CO and isotopologue transitions for the same sample of YSOs.

We start by introducing the selected sample, the studied H₂O and ¹²CO observations and the reduction and decomposition methods in Sect. 2. The results from the water line profile and line luminosity analysis are presented in Sect. 3 and compared to those obtained for CO. In this section, the excitation conditions are also derived from the line intensity ratios. The interpretation of these results are discussed in Sect. 4, and extrapolated to extragalactic sources. Finally, in Sect. 5 we summarise the main conclusions of this work.

2. Observations

2.1. Sample

The sample of 51 YSOs is drawn from the WISH survey and is composed of 26 low-mass, six intermediate-mass and 19 high-mass YSOs. In addition, the intermediate-mass object OMC-2 FIR 4 (Kama et al. 2013) taken from the “Chemical HErschel Surveys of Star forming regions” key programme (CHESS; Ceccarelli et al. 2010) is added to enlarge the number of sources of this sub-group. The intrinsic properties of each source such as coordinates, source velocity (ν_LSR), bolometric luminosity, distance (d), and envelope mass (M_env) can be found in Mottram et al. (2014), Wampfler et al. (2013) and van der Tak et al. (2013) for the low-, intermediate- and high-mass YSOs respectively.

The sample covers a wide interval of luminosity and each sub-group of YSOs contains a mix of different evolutionary stages: both low- and intermediate-mass objects range from Class 0 to Class I; and the high-mass YSOs from mid-IR-quiet/mid-IR-bright massive young stellar objects (MYSOs) to ultra-compact H ii regions (UC H ii). The focus of this paper is to analyse different physical properties across the luminosity range; trends within the low-mass sample are discussed in Mottram et al. (2014); the intermediate- and high-mass samples are too small to search for trends within their several evolutionary stages.

2.2. Water observations

The targeted para-H₂O 2₀₂–1₁₁ (988 GHz) and 2₁₁–2₀₂ (752 GHz) lines and the ortho-H₂O 3₁₂–3₀₃ (1097 GHz) line were observed with the HIFI instrument. The upper energy level (E_u), frequency, HIFI-band, beam efficiency (η_MB), beam size (θ) and spectral resolution of each water transition are given in Table 1. The beam efficiencies of the different HIFI-bands have been recently updated¹ and in general the values decrease by 15−20% for the band considered here. The presented observations have not been corrected with the new η_MB parameters because the analysis in this paper was completed before the new numbers were available and for consistency with our previous CO study. For completeness, the new beam correction factor of each HIFI-band are listed in Table 1.

The H₂O 2₀₂–1₁₁ line was observed for the entire WISH sample; the 2₁₁–2₀₂ line for 24 out of the 26 studied low-mass protostars and for all intermediate- and high-mass YSOs; and the 3₁₂–3₀₃ transition was observed for only 10 low-mass protostars, two out of six intermediate-mass sources and all high-mass YSOs. These water lines are detected for all observed intermediate- and high-mass YSOs and for 75% of the observed low-mass protostars (see Mottram et al. 2014).

The data were observed simultaneously by the Wide Band Spectrometer (WBS) and the High Resolution Spectrometer (HRS), in both vertical (V) and horizontal (H) polarisation (more details in Roelfsema et al. 2012). We present the WBS data because the baseline subtraction for the HRS data becomes less reliable due to the width of the water lines, which is comparable to the bandwidth of the HRS setting. Single pointing observations were performed for all targets in dual-beam-switch (DBS) mode with a chop throw of 3′. Contamination from emission by the off-position has only been identified in the H₂O 2₀₂–1₁₁ spectrum of the low-mass protostar BHR71 (further information in Mottram et al. 2014). The allocated observation numbers for each source and line, designated with the initial obsIDs, are indicated in Table A.2 of Mottram et al. (2014) for the low-mass protostars, and in Table A.1 of this manuscript for the intermediate- and high-mass YSOs.

2.3. Additional ¹²CO observations

Complementing the water observations, ¹²CO J = 10–9 and 16−15 spectra are included in this study to extend the comparison to other components of the protostellar system traced by this molecule and set a reference for abundance studies. The J = 10–9 transition was observed as part of the WISH key programme for the entire low- and intermediate-mass sample and for the high-mass object W3-IRS5 (see San José-García et al. 2013). The J = 10–9 spectrum was obtained for AFGL 2591 from the work of Kaźmierczak-Barthel et al. (2014). For the other high-mass sources, ¹²CO J = 3–2 spectra are used as a proxy (San José-García et al. 2013).

¹²CO J = 16–15 observations of 13 low-mass Class 0 protostars were observed within the OT2_lkrist01_2 Herschel programme (Kristensen et al., in prep.). Finally, this transition was obtained for three high-mass YSOs: W3-IRS5 (OT2_swampfle_2 Herschel programme; Wampfler et al. 2014), and for AFGL 2591 (Kaźmierczak-Barthel et al. 2014) and NGC6334-I (Zernickel et al. 2012) as part of the CHESS key programme (Ceccarelli et al. 2010).

2.4. Reduction of the H₂O data

The calibration process of the water observations was performed in the Herschel Interactive Processing Environment (HIPE²; Ott 2010) using version 8.2 or higher. The intensity was first converted to the antenna temperature ${\mathit{T}}_{\mathrm{A}}^{\mathrm{\ast }}$ scale and velocity calibrated with a ν_LSR precision of a few m s^-1. Further reduction was performed with the GILDAS-CLASS³ package. The spectra observed in the H and V polarisations were averaged together to improve the signal-to-noise ratio (S/N) and the intensity scale converted to main-beam brightness temperature scale, T_MB, using the specific beam efficiencies for each band (Roelfsema et al. 2012). Finally, a constant or linear baseline was subtracted.

All data were then resampled to 0.27 km s^-1 in order to compare the results among the water lines and to those from the high-J CO lines (San José-García et al. 2013) in a systematic manner. The rms noise of the spectra at that resolution, the maximum peak brightness temperature, ${\mathit{T}}_{\mathrm{MB}}^{\mathrm{peak}}$, the integrated intensity, W = ^∫T_MBdν, and the full width at zero intensity, FWZI, are presented in Tables A.2 to A.4. The latter parameter was measured as in Mottram et al. (2014): first by resampling all spectra to 3 km s^-1 to improve the S/N, then re-calculating the rms and finally considering the “zero intensity” where the intensity of the spectrum drops below 1σ of that rms. The velocity range constrained by the FWZI is used to define the limits over which the integrated intensity of the line is calculated.

Since the spectra have not been corrected with the recently released beam efficiencies of the different HIFI-bands¹ (Sect. 2.2), the results presented in Tables A.2 to A.4, as well as those shown in Figs. 3 and 9, should be divided by the correction factor indicated in Table 1. The line profiles do not change if the new η_MB values are applied and the variation of the line ratios is of the order of 1%.

Finally, the C¹⁸O J = 10–9 emission line is detected in the line wing of the H₂O 3₁₂–3₀₃ spectrum for five low-mass protostars (NGC 1333 IRAS2A and IRAS4B, Ser SMM1, GSS30 and Elias 29) and four high-mass YSOs (G5.89-0.39, W3-IRS3, NGC 6334-I and W51N-e1). Therefore, a Gaussian profile with the same FWHM, position and amplitude as those constrained in San José-García et al. (2013) was used to remove the contribution of C¹⁸O J = 10–9 line from the reduced H₂O 3₁₂–3₀₃ spectrum for each of these sources. The data of these specific YSOs are then analysed and plotted after subtracting the C¹⁸O line.

2.5. Decomposition method

As shown by Kristensen et al. (2010, 2012), van der Tak et al. (2013) and Mottram et al. (2014), the water line profiles are complex and can be decomposed into multiple velocity components. The purpose of decomposing the line profile is to disentangle the different regions probed within the protostar, which are characterised by specific physical conditions and kinematics.

Generally, these velocity components can be well reproduced by Gaussian-like profiles; other types of profiles do not give improved fits (Mottram et al. 2014). Depending on the water transition and luminosity of the source, the number of components needed to fit the profile varies. For most of the low-mass protostars, the spectra can be decomposed into a maximum of four different Gaussian components, as shown in Kristensen et al. (2010, 2012, 2013) and Mottram et al. (2014). In order to determine the number of velocity components of the water lines for the intermediate- and high-mass YSOs, these spectra were initially fit with one Gaussian profile using the idl function mpfitfun. Then, the residual from this fit was analysed and since it was larger than 3 sigma rms for all lines, an extra Gaussian component was added to the decomposition method to improve the fitting. The procedure is the same but now considering two independent Gaussian profiles. A self-absorption feature at the source velocity was detected in the H₂O 2₀₂–1₁₁ line for 9 out of 19 high-mass objects, so for those objects an extra Gaussian component in absorption was added in the decomposition method. In some high-mass sources this component is weaker or non-detected in the other studied transitions and it can be negligible (for an example, see the DR21(OH) observations). Therefore, the number of components is determined by the spectrum itself and its S/N and not by the assumed method.

As in Mottram et al. (2014), we force the FWHM and central position of each component to be the same for all H₂O transitions of a given sources. This procedure is adopted because, as for the water observations of the low-mass protostars, the width of the line profiles does not change significantly between the different observed transitions of a source (see also Fig. 2 of Kristensen et al. 2013, for several low-mass protostars), suggesting that in each case the emission from the excited water lines comes from the same parcels of gas. While these two parameters are constrained simultaneously by all available spectra of a given YSO, the intensities of each Gaussian component are free parameters that can vary for each transition. The resulting fits were examined visually as a sanity check. The values of the FWHM, T_peak, ν_peak, and integrated intensity of each Gaussian component are summarised in Tables A.5 to A.7.

2.6. Association with physical components

The multiple velocity components needed to reproduce the H₂O line profiles can be related to physical components in protostellar systems (Kristensen et al. 2011, 2012; van der Tak et al. 2013; Mottram et al. 2014).

Quiescent inner envelope gas produces a Gaussian profile in emission with the smallest FWHM centred at the source velocity (see Sect. 3.2.2 of Mottram et al. 2014). In previous studies these velocity components were called narrow components. The cold outer protostellar envelope can cause a self-absorption, which is more common in ground-state H₂O lines and in objects with massive and cold envelopes (e.g. van der Tak et al. 2013; Mottram et al. 2013).

The chemical and physical conditions present in shocks increase the abundance of water molecules by sputtering from the grain mantles (Codella et al. 2010; Van Loo et al. 2013; Neufeld et al. 2014) and/or by the action of the high-temperature water formation route in the warm post-shock gas (van Dishoeck et al. 2013; Suutarinen et al. 2014). The line profile resulting from shocked water emission depends on the nature and kinematical properties of the shocks generating it, which translate into velocity components with different features (see Table 3 and Sect. 3.2 of Mottram et al. 2014).

The emission from non-dissociative shocks in layers along the outflow cavity wall produces velocity components with the largest FWHM (>20 km s^-1) and are generally centred near the source velocity (Kristensen et al. 2010, 2013; van Kempen et al. 2010; Nisini et al. 2010; Suutarinen et al. 2014; Santangelo et al. 2014). However, these broad Gaussian-like profiles, named cavity shock components (Mottram et al. 2014) or simply broad components, should be differentiated from the broad velocity component identified in low- and mid-J (J< 10) CO spectra, even if shape and width are similar. The reason is that the water emission originates in shocks in the cavity while the CO emission comes from cooler material deeper in the wall and closer to the quiescent envelope (Raga et al. 1995; Yıldız et al. 2013).

In contrast, spot shocks occur in small localised regions and are associated to hotter and more energetic dissociative shocks (Mottram et al. 2014). This emission may originate in extremely-high velocity (EHV) gas along the jets (Bachiller et al. 1990; Tafalla et al. 2010; Kristensen et al. 2011) or at the base of the outflow cavity (previously referred to as either the medium or the offset component; Kristensen et al. 2013). These Gaussian profiles show smaller FWHMs than those measured for the cavity shock component and are in general more offset from the source velocity. A more detailed characterisation and discussion can be found in Mottram et al. (2014), van der Tak et al. (2013) and Kristensen et al. (in prep.).

The contribution of the cavity shock and envelope components with respect to the total integrated intensity of the water lines for the low-, intermediate- and high-mass YSOs are presented in Table 2 together with the analogous contribution from the entrained outflowing material (broad) and envelope gas (narrow) components for the ¹²CO J = 10–9 line (ratios derived from San José-García et al. 2013). The values for the low-mass Class 0 and Class I protostars were calculated by Mottram et al. (2014) for different water transitions as well as the fraction corresponding to the spot shock component (see Table 4 of that manuscript).

In the following, the different velocity components of the water and CO lines are distinguished according to the terminology based on the probable physical origin of the emission and the width of the profile.

3. Results

The basic properties of the spectra and their decomposition are introduced in Sect. 3.1. In Sect. 3.2 the results of the line profile decomposition are compared to those obtained for the high-J CO observations described in San José-García et al. (2013). The water line luminosity properties are also compared to those of CO in Sect. 3.3. In Sect. 3.4 integrated intensity ratios calculated for different pairs of water transitions are presented and these line ratios are further studied as across the line profiles in Sect. 3.5. Finally, the excitation conditions of the water lines are analysed with the non-LTE radiative transfer code radex in Sect. 3.6.

3.1. Water line profile characterisation

The observed H₂O 2₀₂–1₁₁, 2₁₁–2₀₂, and 3₁₂–3₀₃ spectra for the intermediate- and high-mass YSOs are presented in Figs. A.1−A.3, respectively. The Gaussian profile fitting the broad (cavity shock) component of each water transition and source is indicated with a pink line. The spectra of the low-mass protostars are shown in Appendix A of Mottram et al. (2014).

In order to easily compare all the data, normalised averaged spectra of the H₂O 2₀₂–1₁₁, 2₁₁–2₀₂ and 3₁₂–3₀₃ transitions are computed for the low-mass Class 0 and Class I protostars, the intermediate-mass sources and high-mass YSOs, as shown in the three first panels of Fig. 1. These spectra are calculated for each transition by shifting each spectrum to zero velocity, normalising it to its peak intensity and averaging it together with the observations of the corresponding sub-group of objects. The presence of self-absorption features for some of these sources, which are stronger for certain water transitions, prevents the normalised averaged spectra to reach unity for several of these lines. Independently of this issue, the normalised averaged spectra obtained for the three H₂O transitions are similar for each sub-type of YSOs. Only the H₂O 2₀₂–1₁₁ high-mass averaged spectrum shows a slightly different profile with respect to the other two water lines because a larger number of YSOs show strong and deep self-absorption features. Except for the Class I protostars, the averaged spectrum for a given transition seems to be broader for the low-mass objects, but at the base of the spectra the widths are similar, independent of the luminosity.

In the right-hand panel of Fig. 1 the three H₂O transitions (2₀₂–1₁₁ in red, 2₁₁–2₀₂ in blue and 3₁₂–3₀₃ in purple) are over-plotted for four different sources: a low-mass Class 0 (NGC 1333 IRAS4B), a low-mass Class I (GSS 30 IRS), an intermediate-mass objects (NGC 2071) and a high-mass YSO (W33A). For each source, the shapes of the three water line profiles are similar but scaled-up in intensity, a result that is confirmed from the visual inspection of the water line profiles of all YSOs. In particular, the line wings are very similar. This indicates that the three water transitions are probing the same dynamical properties in each source.

Moving to the outcomes from the line decomposition explained in Sect. 2.5, the analysis suggests that the quiescent envelope and cavity shock components are the only two physical components consistently present in the H₂O 2₀₂–1₁₁, 2₁₁–2₀₂ and 3₁₂–3₀₃ spectra of all low-, intermediate- and high-mass YSOs. The spot shock components are not detected for the excited water lines presented here towards high-mass YSOs and six out of seven intermediate-mass objects, though they have been seen in absorption against the outflow in some ground-state H₂O lines for some high-mass sources (for more information see van der Tak et al. 2013).

As shown in Table 2, for a given sub-sample of YSOs the averaged contribution of the cavity shock component with respect to the total integrated intensity of the line is the same for the three water transitions. This fraction seems to decrease from low- to high-mass, but no statistically significant trend with L_bol can be claimed because the specific contribution of the cavity shock emission varies from source to source. The remaining emission comes from the envelope in the case of the low-mass Class I, intermediate- and high-mass YSOs and from spot shock components for low-mass Class 0 protostars (Mottram et al. 2014). This picture is consistent with the average spectra presented in Fig. 1, where the envelope component of the water lines is more prominent for the high-mass sources than for their low-mass Class 0 counterparts. In addition, this narrow component associated with the envelope is also less prominent in the H₂O 3₁₂–3₀₃ transition regardless of the YSO mass, as expected since the envelope is presumably composed of cool quiescent gas.

Independently of these numbers, in this paper we focus on characterising the physical conditions causing the line-wing emission in the water line profiles, i.e., the broader velocity component associated to the shock emission along the outflow cavity.

3.2. Comparison of the H₂O and ¹²CO line profiles

The fourth panel in Fig. 1 includes the normalised averaged ¹²CO J = 10–9 spectrum of each sub-sample of YSOs. The procedure followed to obtain these spectra is the same as that used for the water data. These averaged profiles are clearly narrower than those of water, especially for the low-mass Class 0 protostars, and the width of the spectra seems to increase from low- to high-mass. Therefore, just from a basic visual inspection of the water and the high-J CO normalised averaged spectra we can point to differences in the shape of the line profiles of these two molecules and a different trend in the width from low- to high-mass.

To consistently compare the dynamical conditions of the entrained outflowing material traced by CO and the shocked gas along the outflow cavity, the line-wing emission of these two species is analysed using two parameters: the FWHM and the FWZI (see Sects. 2.4 and 2.5). Both variables are used because FWHM characterises the average extent of emission from the source velocity while FWZI characterises the fastest material. For simplification, the FWHM of the Gaussian profile reproducing the cavity shock component is differentiated from the other velocity components by using the sub-script b to indicate that this is the broader velocity component obtained from the line decomposition.

The top panel of Fig. 2 shows the FWZI for the H₂O 2₀₂–1₁₁ transition as a function of bolometric luminosity (in Fig. C.2 the FWZI of the other water lines are also plotted versus L_bol and envelope mass). Similarly, the constrained FWHM_b from the cavity shock component (same for all three lines) versus L_bol and M_env are presented in Fig. C.1. The FWZI values vary from 15 to 189 km s^-1, while the FWHM_b range from 13 to 52 km s^-1. Due to the large scatter and dispersion of the data points, no trend or correlation with luminosity can be claimed in either case. The smaller FWZI and FWHM_b values are those of the low-mass Class I protostars, consistently lying at the bottom of these figures. In addition, the low- and intermediate-mass YSOs which show EHV components are marked with squares in Fig. 2 to indicate that their FWZI was calculated including the spot shock emission and to investigate if there is any particular trend for these objects. The spectra of the marked high-mass object do not have EHV components but their line profiles are characterised with broad and triangular shapes. More information about these specific sources can be found in Appendix B.

As indicated in Sect. 2.3, the ¹²CO J = 10–9 transition was not observed for most of the high-mass YSOs. However, those sources for which both ¹²CO J = 10–9 and J = 3–2 transitions were available, the values of the constrained FWHM_b and also FWZI are similar within the uncertainty (see San José-García et al. 2013). Therefore, the J = 3–2 transition is used as a proxy of the J = 10–9 spectra for the study of the kinematical structure of the outflowing gas. The FWZI and FWHM_b for the ¹²CO observations (middle panels of Figs. 2 and C.1 respectively) are spread across a smaller velocity range than that for water. In addition, the FWZI shows less scatter than the FWHM_b. There is a statistically significant trend of larger FWZI for more luminous sources (4.7σ⁴) with a Pearson correlation coefficient r = 0.72, which is also seen to a lesser extent for the FWHM_b.

Table 3 presents the mean FWHM_b and FWZI values for H₂O and ¹²CO and the averaged rms in a 0.27 km s^-1 bin, σ_rms. In the case of the high-mass YSOs, the derived values are not affected by the higher σ_rms in those data (as left panels of Fig. 1 already show) since the actual signal to noise, S/N, on the water spectra themselves given by the peak intensity relative to the rms are higher (averaged S/N value of ~60) than those of their low-mass counterparts (averaged S/N of ~20).

Without considering the low-mass Class I protostars, which are more evolved and therefore have weaker, less powerful outflows (Mottram et al. 2014), the average FWHM_b and FWZI values derived for H₂O are similar from low- to high-mass (Table 3). A decrease of the mean FWZI values with increasing luminosity is only hinted at for the H₂O 2₀₂–1₁₁ transition. Combining the results from both FWHM_b and FWZI we conclude that the extent of the water line emission is similar for the entire sample. In contrast, and as suggested by the middle panel of Fig. 2, the averaged values of both FWZI and FWHM_b for the ¹²CO observations seem to increase with luminosity.

The dispersion observed for the FWZI and FWHM_b in both H₂O and CO could be related to the intrinsic properties of the source, such as its inclination, evolutionary stage, clustering, etc. In order to minimise possible effects caused by these inherent characteristics, the ratio of the FWZI derived for the ¹²CO observations and the FWZI of the water lines is plotted versus the bolometric luminosity in the bottom panel of Fig. 2. The same procedure is followed for FWHM_b of the ¹²CO and water spectra (see bottom panels of Fig. C.1).

Independently of the use of FWZI or FWHM_b, a correlation between these ratios and L_bol is measured for each of the three water lines with statistical significance between 3.3σ and 5.0σ (Pearson correlation coefficients between 0.50 and 0.75). While ¹²CO J = 3–2 and 10–9 may trace different layers in the outflow (Yıldız et al. 2013; Santangelo et al. 2013), the FWZI, i.e., the maximum offset velocity (ν_max), of the CO lines does not change with the J transition (Kristensen et al., in prep.). Therefore, the use of ¹²CO J = 3–2 as a proxy for ¹²CO J = 10–9 for the high-mass sources will not affect this trend.

The average ratios of FWZI and FWHM_b are also given in Table 3 and generally increase with increasing L_bol. The FWZI ratio indicates that the velocity of the material traced by the wings of ¹²CO is larger than that of water for the high-mass YSOs. The FWHM_b values are basically the same for ¹²CO and H₂O, independently of whether these molecules are probing different regions and physical conditions within the outflow. However, this is not the case for the low-mass protostars, which show larger line-wings for the water observations than for the ¹²CO observations, consistent with Kristensen et al. (2012). This trend is further analysed and interpreted in Sect. 4.1.

3.3. Line luminosity study

The integrated intensity of the water emission lines, W = ^∫T_MBdν, is calculated as described in Sect. 2.4 and presented, together with its uncertainty, in Tables A.2 to A.4. For those H₂O 2₀₂–1₁₁ spectra with an absorption at the source velocity, the integrated intensity was determined by masking the absorption feature and using the area traced by the two Gaussian emission profiles fitting the spectrum.

Since the studied sample of YSOs covers a wide range of luminosities and distances (from 0.13 to 7.9 kpc), the integrated intensity is converted to line luminosity, L_H2O. This parameter takes into account the distance and can be compared to values obtained in extragalactic studies. These quantities have been calculated using Eq. (2) of Wu et al. (2005), assuming a Gaussian beam with size according to Table 1 and point source objects. The uncertainties in the line luminosity are calculated from the rms of the spectrum and assuming a distance uncertainty of ~20%.

Figure 3 presents the logarithm of the line luminosity for each water transition, log(L_H2O), versus the logarithm of the bolometric luminosity, log(L_bol). A strong correlation between the logarithm of these quantities is measured (solid black line) for each of the three lines. The Pearson correlation coefficient, r, is larger than 0.9 and the trends extend over more than six orders of magnitude in both axes. Similarly, a strong correlation is also measured between the logarithm of L_H2O and the logarithm of the envelope mass, M_env (see Fig. C.4). The parameters of the fit and correlation coefficient for all water lines and for both L_bol and M_env are presented in Table 4.

The slope of the correlations are close to unity within the uncertainty. As there is no saturation of the line luminosity, the emission must be optically thin or at least effectively thin (i.e. optically thick but sub-thermally excited, see Linke et al. 1977). In either case, the intensity scales as the product of the beam-averaged column density and average volume density in the emitting gas. Assuming that the volume and column density in the regions probed by water in low- and high-mass protostars and the source geometry are similar, a correlation between log(L_H2O) and log(L_bol) with slope close to unity suggests that the emitting region size increases proportionally with L_bol (see Sects. 3.6 and 4.3 for further discussion). The 3₁₂–3₀₃ transition shows the smallest number of detections for the low-mass protostars, which could explain the difference in the value of the slope with respect the other two water lines. This linear correlation of log(L_H2O) and log(L_bol) should not be extrapolated to the ground-state water transitions due to the presence of broad and deep absorption features which complicate the comparison of the water line luminosity, particularly in the case of massive YSOs (see van der Tak et al. 2013).

The same correlation between the ¹²CO 10–9 line luminosity and L_bol was measured by San José-García et al. (2013, their Fig. 6). The calculated slope of the linear fit is the same within the uncertainty to those derived for the water lines (see Table 4 of San José-García et al. 2013, for comparison). As also noted by Mottram et al. (2014), the low-mass Class I sources have lower water line luminosities than the Class 0 sources, due to decreases in both FWHM and FWZI (see Table 2), as well as peak intensities, causing a clear separation between these two classes of low-mass protostars. However, in the case of CO, the separation between the low-mass Class 0 and Class I was small. Both ¹²CO J = 10–9 and especially water are more sensitive to fast moving material in the outflows, which differentiate sources at different evolutionary stages, while the CO isotopologues are mostly probing the bulk of quiescent envelope material. This confirms the conclusions of Mottram et al. (2014) that the difference in water emission between Class 0 and Class I sources is related to the decreasing strength of the outflow and not to the removal of the envelope.

The integrated intensity emission coming from the cavity shocks, i.e., the area of the Gaussian profile fitting the line-wings, is also converted to line luminosity, L_{broad H2O}. Figure C.5 shows the logarithm of this quantity as a function of log(L_bol). As for the total line luminosity, a strong correlation which extends across the entire luminosity range is measured but with slightly smaller values of the slope (averaged slope of 0.78±0.04). Figure 3 includes for each water transition the relation between the logarithm of L_{broad H2O} and L_bol (dashed black line) for comparison with the total line luminosity, L_H2O, (solid black line). The parameters from the log-log correlations of the broad component are also presented in Table 4. Consistent with the larger envelope contribution for high-mass sources (see Table 2) the differences between the fits for L_H2O and L_{broad H2O} are larger for all three lines for the high-mass sub-sample.

The fact that we found the same strong correlation between log(L_bol) and log(L_H2O) or log(L_{broad H2O}), i.e., similar values for the slope close to unity, confirms that the total water line luminosity is dominated by the broad component (outflow cavity shocks).

Finally, we found that the line width of the cavity shocked gas seems independent of the source luminosity (see previous section), i.e., the kinetic energy of the outflow increases little from low- to high-mass YSOs. Thus, this suggests that the physical processes originating and powering the outflows (the dynamical conditions) are similar from low- and high-mass and it is the mass of the outflow, i.e., the total cavity shocked water gas, what increases at higher L_bol through a larger emitting area.

3.4. Integrated intensity ratios

The integrated intensity ratios of two different water transitions across the luminosity range can be used to analyse the excitation and physical conditions of protostellar systems. Mottram et al. (2014) show that there is little variation in line ratio between outflow-related components for the low-mass sample. They showed that the lines lie in the optically thick but effectively thin regime. This means that the water lines presented in that study are sub-thermally excited due to their high critical densities compared to the actual gas density of the emitting region. This condition lowers the chance of collisional de-excitation, so photons will scatter and eventually escape as if the lines were optically thin. We will extend this study to intermediate- and high-mass YSOs, first considering the total integrated intensity ratios and then focusing on just the cavity shock component.

The observations of the different water transitions have different beam sizes which must be corrected by the corresponding beam-size ratio, θ₁/θ₂. However, the correction depends on the size and shape of the emitting region within the protostellar system covered by each beam. If the emission comes from a point source, the correction factor is (θ₁/θ₂)²; if the emitting area is conical, meaning that it fills the beam in one axis and is point-like in the other, then the correction applied should be (θ₁/θ₂); and if the emitting region covers both axes the factor used is 1 (Tafalla et al. 2010). For the considered water lines, the maximum scaling correction would be 2.1, corresponding to a point-source scaling of the 2₁₁–2₀₂ (752 GHz) and 3₁₂–3₀₃ (1097 GHz) beams.

The emission from the excited water transitions for the low-mass protostars mostly comes from the outflow cavity, which extend beyond a single Herschel beam (see Nisini et al. 2010; Santangelo et al. 2012, 2014; Mottram et al. 2014). In this case, the emitting area could be taken as conical or even assumed to fill both beams. So, either the beam correction factor is (θ₁/θ₂) or 1 (or something in between for other, more complex geometries). For the intermediate- and high-mass YSOs it is more complicated to define the area covered by the beams as the larger distances mean that the outflows are often unresolved in the Herschel beam. For this reason, we consider the three scenarios introduced in the previous paragraph and the entire range of intensity ratios (applying the three possible beam correction factors) are evaluated together with the results obtained from their low-mass counterparts. The θ₁/θ₂ value for each pair of water lines are indicated in Table 5.

The ratios of the total integrated intensity of H₂O 2₁₁–2₀₂/ 2₀₂–1₁₁; 3₁₂–3₀₃/2₀₂–1₁₁; and 3₁₂–3₀₃/2₁₁–2₀₂ are shown as a function of bolometric luminosity in the top-, middle- and bottom-left panels of Fig. 4, respectively. The value of these ratios without applying the beam corrections are indicated by filled dots. The dashed lighter line connected to each dot represents the value of the ratio after applying the larger beam size correction assumed for each sub-sample of objects: θ₁/θ₂ linear correction for the low-mass protostars and (θ₁/θ₂)² for the intermediate- and high-mass YSOs. Therefore, the lighter vertical lines illustrate the range of values that the integrated intensity ratio could take for each source depending on the shape of the emission region.

Figure 4 also shows the analytical ratios of the integrated intensity derived for two extreme scenarios: when both water lines are optically thin (purple dotted-dashed line) and when they are optically thick (orange dashed line). In both cases local thermodynamical equilibrium (LTE), an excitation temperature, T_exc, of 300 K and the same beam size for the compared transitions are assumed. The specific numbers corresponding to each regime for each pair of water lines are presented in Table 5. We take T_exc ≈ 300 K because this is the temperature found for the cavity shock component by (Mottram et al. 2014). Also, Karska et al. (2013, and references cited) observe with Herschel PACS a warm component for CO with an excitation temperature around 300 K and water is found to be spatially associated with this component in the outflow cavity (Santangelo et al. 2013). Note that for high densities, the excitation temperature approaches the kinetic temperature.

The ratio of the 2₁₁–2₀₂ (752 GHz) and 2₀₂–1₁₁ (988 GHz) transitions (top-left panel of Fig. 4) shows a significant increase with L_bol. The trend of the non-beam-corrected ratios (dots) is indicated by the dotted black line and confirmed by a Pearson correlation coefficient, r, equal to 0.83. This increase with luminosity is also seen when the maximum beam correction factor for each sub-group of YSOs (delineated by the dash) is applied (r = 0.85). Therefore, the result is consistent across the luminosity range and regardless of emitting region shape. Furthermore, the intensity ratio of these H₂O transitions is above 1 only for those sources with L_bol> 10³ L_⊙.

The 3₁₂–3₀₃ (1097 GHz) to 2₀₂–1₁₁ (988 GHz) intensity ratio also seems to increase slightly with luminosity (middle-left panel Fig. 4). However, no trend can be claimed on a statistical basis because the significance of the Pearson correlation coefficient is less than 2.5σ (r< 0.5 and N = 32). These ratios lie below 1 for all sources, independent of whether the beam correction factor has been applied (which is small for this combination of lines). The intensity ratio of the 3₁₂–3₀₃ (1097 GHz) to 2₁₁–2₀₂ (752 GHz) (bottom-left panel) is generally smaller than 1 for the three sub-sample of YSOs and it seems to decrease with L_bol. Once again, a trend cannot be claimed because the significance of this correlation is less than 2.3σ (r< −0.5 and N = 31).

Comparison of the left and right panels of Fig. 4 show that similar results are obtained if the integrated intensity ratio corresponding to only the broad component is used for the same pair of water transitions as a function of L_bol. A trend with luminosity is also only observed for the 2₁₁–2₀₂/2₀₂–1₁₁ intensity ratio (with a smaller correlation coefficient r = 0.66). The increase of the 2₁₁–2₀₂/2₀₂–1₁₁ ratio with L_bol suggests that the emission from the shocked gas along the outflow cavity becomes more excited, i.e., warmer and/or denser, for more massive objects. The 3₁₂–3₀₃/2₀₂–1₁₁ ratio seems to increase with L_bol while the 3₁₂–3₀₃/2₁₁–2₀₂ ratio decreases with luminosity, both of them tendencies not statistically significant.

The mean and standard deviation of the intensity ratios for the broad component are presented in Table 5, together with the beam size ratios and the estimated intensity ratios for the optically thin and thick solutions calculated assuming LTE and T_exc = 300 K. These calculated values ratify the results and tends from low- to high-mass described in the previous paragraph.

We can rule out the optically thin LTE solution for all sources and lines, except when the largest beam size correction is applied to the massive YSOs in the case of the 2₁₁–2₀₂/2₀₂–1₁₁ integrated intensity ratio. The optically thick LTE solution for the 3₁₂–3₀₃/2₀₂–1₁₁ ratios can be also excluded for the entire sample of objects. In order to further constrain the excitation conditions, the non-LTE radiative transfer code radex will be used in Sect. 3.6.

3.5. Intensity ratios versus velocity for H₂O and ¹²CO

Here we investigate whether the excitation conditions change with velocity. Figure 5 shows the calculated line ratios of the H₂O 2₁₁–2₀₂ (752 GHz; top-panel) and 3₁₂–3₀₃ (1097 GHz; bottom-panel) transitions over the 2₀₂–1₁₁ (988 GHz) line as a function of velocity with respect to the source velocity. The grey shaded areas correspond to the averaged line ratio and standard deviation of the low-mass Class 0 protostars presented in Mottram et al. (2014). The green and red lines are associated with the intensity ratios of the intermediate- and high-mass YSOs respectively.

To obtain these ratios, the spectra of each set of water transitions for a given source are resampled to 1 km s^-1 bins. Then, the studied line (either the 752 or 1097 GHz transition) is divided by the spectrum of the 988 GHz line within the range of velocities over which the intensities of both lines are above 3σ of the re-binned data. Finally, the red and blue wings are averaged together and the obtained intensity ratio as a function of velocity is averaged over all sources that compose each sub-sample of YSOs. Only velocity ranges well offset from the quiescent envelope emission are plotted, typically | ν− ν_lsr | > 3−5 km s^-1.

As found for the low-mass Class 0 protostars, the line intensity ratios seem to be constant as a function of velocity for both intermediate- and high-mass YSOs. This result confirms the fact that the shape of the three water line profiles is similar for a given source and just scales in intensity (Fig. 1, right). For the high-mass objects, the ratio increases slightly for offset velocities >25 km s^-1. However, the significance of the increase in the ratio is of the order of the uncertainty of those channels. The 2₁₁–2₀₂/2₀₂–1₁₁ ratio is above 1 for this sub-sample of objects at all velocities, which is consistent with the result obtained in the previous sub-section for this set of water lines. A constant line ratio as a function of velocity rules out the optically thick LTE solution for all lines as this would require high column densities even at high velocities, and thus a flatter line profile than observed, or high opacity across the entire profile.

A different outcome is obtained if the CO line intensity is compared to that of water as a function of velocity. Figure 6 shows the ¹²CO J = 10–9 and J = 16–15 over the H₂O 2₀₂–1₁₁ (988 GHz) line ratio as a function of velocity offset (top and bottom panel respectively). In this case, we observe two different behaviours for the two sub-samples of YSOs. For the low-mass Class 0 protostars the ¹²CO J = 10–9/H₂O 2₀₂–1₁₁ line ratio decreases with velocity, while the ¹²CO J = 16–15/H₂O 2₀₂–1₁₁ ratio is invariant with velocity.

On the other hand, for the high-mass YSOs the line ratio of both CO transitions with H₂O 2₀₂–1₁₁ is constant with velocity across the studied velocity interval. This agrees with the ratio of the FWZI and FWHM_b for these two molecules being slightly greater and around unity respectively (see red crosses in Figs. 2 and C.1). Therefore, in massive objects the emission of these CO and water transitions may originate in regions of the outflow cavity wall not too distant from each other and characterised by similar excitation and dynamical conditions. For the low-mass Class 0 sources, the ¹²CO 16–15 line seems to trace the same material as H₂O but the 10–9 transition does not (Kristensen et al., in prep.).

3.6. Excitation conditions

In order to characterise the excitation conditions responsible for the broad cavity shock water emission, we compared the line fluxes presented above, together with the fluxes for the H₂O 1₁₁–0₀₀ transition from M^cCoey et al. (in prep.) and van der Tak et al. (2013), for the intermediate and high-mass YSOs to a grid of non-LTE radex calculations (van der Tak et al. 2007). We use the same modifications and parameters as adopted by Mottram et al. (2014) for the low-mass cavity shock components and we assume 1D slab geometry, an ortho-to-para ratio of H₂ of 3 and a circularly emitting region at the distance of the source. The area of this emitting region in the plane of the sky is derived by comparing the modelled and observed H₂O 2₀₂–1₁₁ integrated intensity and by determining the fraction of the beam for which the modelled and integrated intensity values match.

We find a good fit with post-shock volume densities n_H2 = 10⁵−4 × 10⁶ cm^-3 and water column densities N_H2O = 10¹⁷–10¹⁸ cm^-2 for 5/7 of the intermediate-mass sources. The results for NGC 7129 are shown as an example in the top panels of Fig. 7. These results are similar to lower-density solution obtained for the low-mass sources by Mottram et al. (2014) but with slightly larger emitting regions sizes, equivalent to radii of 300−500 AU if circular. All studied water transitions are optically thick but effectively thin in these conditions. In one of the remaining two intermediate-mass sources (L1641 S3MMS1) we cannot constrain the fit very well, while for the other source (NGC 2071) a good fit cannot be found within the grid (${\mathit{\chi }}_{\mathrm{best}}^{\mathrm{2}}\mathit{>}\mathrm{100}$), with the emission in the H₂O 2₁₁–2₀₂ line underproduced and emission in the 1₁₁–0₀₀ line overproduced by the models. For the high-mass YSOs, good fits can only be found for 6 of the 19 sources, with the best fits again giving similar densities and column densities and emitting regions ranging from 800 to 6000 AU. However, even in these cases the best-fit models consistently underproduce the H₂O 2₁₁–2₀₂ transition. To conclude, if we assume a circular geometry for all sources, the radius of this emitting region is ~10–200 AU for low-mass protostars, ~200–600 AU for intermediate-mass YSOs and ~600–6000 AU for non-pumped high-mass sources.

In order to find a good solution for these high-mass sources, we first ran additional grids of models for the poorly fit source G10.47+0.03 where the kinetic temperature was increased from the default value of 300 K up to 1500 K. However, the radex models continued underproducing the H₂O 2₁₁–2₀₂ transition, independently of the adopted temperature, and no better solutions were found. This suggests that the excitation of these observed water lines (with E_u in range up to ~250 K) is not sensitive to temperature, as also found for their low-mass counterparts (cf. Fig. 11 of Mottram et al. 2014). Although a realistic shock environment has a range of temperatures this would provide similar line ratios within the uncertainty. This scenario is confirmed by the fact that the intensity ratios do not change with velocity (see Fig. 5), even if the line profile is composed of cavity and spot shock emission, i.e., gas which is at different temperatures (see Fig. 9 of Mottram et al. 2014). In addition, this conclusion is consistent with the results of Karska et al. (2014a), who find similar rotational temperatures for warm CO emission observed with PACS for both low- and high-mass sources.

Since the observed H₂O 2₁₁–2₀₂/2₀₂–1₁₁ line ratio for most high-mass YSOs is larger than unity, an alternative scenario is explored for one of these massive YSOs (G10.47+0.03). In this case, we study the effect of pumping by an infrared radiation field with the same shape as the source SED from van der Tak et al. (2013), and scaled assuming that the source of infrared radiation is at varying distances from the water emitting region. Such an approach was found by Mottram et al. (2014) to be inconsistent with the observations of low-mass protostars, but González-Alfonso et al. (2014) show that this is needed to reproduce H₂O observations towards some galaxies. For the high-mass YSO G10.47+0.03, we find that pumping by infrared photons can indeed provide a better fit to the data (compare middle and lower panels of Fig. 7). This result is consistent with previous studies based on observations of vibrational excited molecules such as HC₃N and HCN (Wyrowski et al. 1999; Rolffs et al. 2011). We explored a range of distances between the source of the infrared radiation and the H₂O emitting region, from 10 to 1000 AU, and find the best fit for a distance of ~100 AU. The best fit for this strong radiation field has both high densities (~10⁹–10¹¹ cm^-3) and column densities (~10¹⁸ cm^-2) for emitting regions with sizes of order 500 AU. All other high-mass sources were then tested over a more limited grid of distances. For 10 out of 19 sources, inclusion of an infrared radiation field improved the fit to the line ratios, while for the remaining 9 YSOs the fit was worse. This latter group contains 3/5 of the mid-IR quiet HMPOs and 3/4 of the hot molecular cores within the WISH sample, but also 3/5 of the UCHii regions. More detailed investigation is needed to ascertain why certain sources seem to shield water from infrared radiation better than others. The implications and reasonableness of these results will be discussed further in Sect. 4.2.

4. Discussion

4.1. Disentangling the dynamical properties of H₂O and CO emission

The data presented in Sect. 3 suggest that H₂O appears in a different physical component from some of the CO lines. It is now well established that the broad component seen in low- (J< 5) and mid-J (5 ≤ J< 10) ¹²CO lines probes the colder entrained outflow gas, whereas only the high-J (J ≥ 10) ¹²CO transitions trace the currently shocked gas also seen in water (van Kempen et al. 2009, 2010; Nisini et al. 2010; Santangelo et al. 2012; Tafalla et al. 2013; Goicoechea et al. 2012; Yıldız et al. 2013; San José-García et al. 2013). The cartoon presented in Fig. 8 illustrates a simplified version of a possible framework for the origin of ¹²CO (yellow, orange and red areas) and water (blue lines and opaque region) within the outflow cavity wall for a low- (LM) and high-mass (HM) YSO. Turbulent motions, likely driven by the outflow (San José-García et al. 2013, 2015) are represented with spiral symbols and the entrainment with swirls.

The yellow areas are related to colder layers of the outflow cavity wall directly attached to the envelope where the broad velocity component of low-J¹²CO transitions originates. The orange regions are associated with warmer entrained outflowing material within the outflow cavity wall and the red areas indicate material that is currently mixing and undergoing shocks along the outflow cavity. Finally, the velocity of the layers forming the outflow cavity wall and the temperature of the gas decrease with distance from the cavity and the central source (Raga & Cantó 1997). The shocked gas along the layers closest to the outflow cavity should have higher average velocities than that originating near the envelope.

In this context, Fig. 2 and Fig. C.1 show that both the values of FWZI and FWHM_b behave differently with L_bol for H₂O and high-J¹²CO. In the case of the low-mass protostars, the FWHM_b of the water lines (which probes the shocked gas along the outflow cavity) is larger than that of ¹²CO J = 10–9 (which traces the warm entrained outflowing material). In addition, the decrease of the ¹²CO J = 10–9/H₂O 2₀₂–1₁₁ ratio as a function of velocity (Fig. 6) suggests that these two transitions are kinematically distinct, and originate in different layers of the outflow cavity wall. On the other hand, the ratio of the CO 16–15/2₀₂–1₁₁ line intensity is constant with velocity, so water and ¹²CO 16–15 emission go together along the outflow cavity wall.

For the high-mass YSOs, the FWHM_b for high-J¹²CO transitions is comparable to that of H₂O (see bottom panel of Fig. C.1). Furthermore, the line intensity ratio of these two species is constant as a function of velocity for both the ¹²CO J = 10–9 and 16–15 transitions. This implies that the emission of the ¹²CO J = 10–9 transition already comes from a region with dynamical properties that are similar to those of water. This is in contrast to what is seen for their low-mass counterparts where ¹²CO J = 10–9 and water do not go together (Yıldız et al. 2013) suggesting that the water emission arises in warmer gas as probed by higher-J CO lines, e.g. ¹²CO J = 16–15 (E_up/k_B = 750 K), consistent with what is seen at off-source shock positions (e.g. Santangelo et al. 2013). In the following we will discuss two scenarios to explain these results.

The first scenario is based on the fact that high-mass objects have stronger UV fields than those measured for their low-mass counterparts (Stäuber et al. 2007). The strong UV radiation would principally destroy water molecules in the warmer layers of the outflow cavity wall closer to the cavity because colder entrained material will be shielded. Therefore, water may be effectively photo-dissociated (and photo-desorbed) deeper into the outflow cavity wall for high-mass objects than for low-mass protostars. CO molecules are not as easily photo-dissociated as H₂O, but still the destruction of CO will be more efficient in the layers nearest to the outflow cavity. This means that for massive YSOs, the water molecules present in layers of the outflow cavity wall closer to the wall, which move at higher velocities, are destroyed while the water molecules located in regions near to the entrained layer, which move at slower velocities, may survive.

Due to the weaker UV fields, the photo-dissociation of water molecules along layers closer to the outflow cavity is less efficient for low-mass protostars. Therefore, a higher percentage of fast shocked gas emission along the regions closest to the cavity may survive, emission which contributes to the broadening of the line-wings (quantitative models in Visser et al. 2012).

An alternative scenario is based on the fact that the protostellar environments of high-mass YSOs are more turbulent than those of their low-mass counterparts (see Herpin et al. 2012; San José-García et al. 2013, 2015, and bottom panel of Fig. 8). Outflow systems inject turbulent motions into the structure forming the outflow cavity wall and these motions could propagate to other physical components of the protostellar system such as the entrained outflowing material. In this case, the material from both components may have a high level of mixing and the separation between the entrained and shocked material within the outflow cavity wall would not be as clearly delineated as for the low-mass protostars. Since the turbulent motions in low-mass protostars are weaker, the different layers constituting the cavity wall are better defined. Thus, the kinematical properties of each J transition are better differentiated and the emission of CO J = 10–9 material is separated from that originating in shocked water gas along the cavity. Therefore, the dynamical properties of these two structures would be more closely linked for the high-mass YSOs than for their less massive counterparts, being harder to disentangle. This could explain why the calculated ratio of the FWZI or FWHM_b for the ¹²CO J = 10–9 and H₂O profiles is around 1 for the high-mass sub-sample but is lower than 1 for their less massive counterparts.

Finally, we can rule out that high-mass YSOs are analogous to low-mass Class I protostars, for which the line profiles are narrower due to decreased wind velocities (cf. Kristensen et al. 2012; Mottram et al. 2014). The high-mass sources considered here are closer to the low-mass Class 0 objects in terms of evolutionary stage because the ratios of the envelope mass to gravitational stellar mass, M_env/M_g, are >20 (San José-García et al. 2015). Massive YSOs are therefore still deeply embedded and follow the relation with luminosity of the Class 0 rather than to the Class I low-mass sources.

4.2. Excitation condition across the luminosity range

As discussed in Sect. 3.6, the excitation conditions of water in low- and intermediate-mass sources seem to be similar. However, for high-mass YSOs it is difficult to reproduce the observed line ratios, particularly those involving the H₂O 2₁₁–2₀₂ line (752 GHz), with a similar model. Since water excitation is not sensitive to temperature (Mottram et al. 2014), the other plausible scenario is that for the majority of the high-mass sources the lines may be pumped by infrared radiation, a process probably not present in low- and intermediate-mass objects due to lower luminosities. However, the best-fit excitation conditions obtained for the high-mass source G10.47+0.03 with the inclusion of a radiation field require very high densities, enough that water emission should be in LTE, over considerable regions. Note that from the grid without radiative pumping we can rule out that high densities alone can reproduce the observed water line ratios. The regions around these high-mass YSOs undoubtably have a complex geometry and in several cases, including G10.47+0.03, are known to harbour multiple sources within the Herschel beam. Therefore, the assumption of a slab geometry, while instructive, may not accurately constrain the range of excitation conditions present. In addition, we have only performed a limited test using a scaled version of the SED as the radiation field, while the “true” radiation field irradiating the water emitting gas in these regions may include shorter-wavelength photons. Furthermore, we assume all the gas is irradiated by the same radiation field, when the geometry of the water emitting region could lead to variation in the field experienced by different parts of the gas.

With all these caveats, we therefore can conclude two things. Firstly, the excitation of water in most high-mass sources is different to low- and intermediate-mass objects. Secondly, while radiative pumping can be ruled out in low-mass protostars, it is a plausible explanation for high-mass sources given that higher kinetic temperatures do not lead to strong enough emission in the 2₁₁–2₀₂ transition. In addition, a similar gas temperature for low- and high-mass YSOs is consistent with the results obtained from PACS observations, which show a comparable CO rotation temperature for low- and high-mass objects (Karska et al. 2014a).

4.3. Correlations with bolometric luminosity

Figure 3 illustrates the relation between the luminosities of the water lines and the bolometric luminosity of the sources from low- to high-mass. The inferred slope on the log-log scale ranges from 0.76 to 0.94. San José-García et al. (2013) found a similar relation for the various high-J CO isotopolog lines, with slopes around unity.

What is the origin of these strong correlations? As shown in Benz et al. (2015, Appendix D), the luminosity of any optically thick line originating from the protostellar envelope scales with the radius, r_thick, where the line becomes optically thick (i.e. the τ = 1 photosphere) and the temperature at that radius: ${\mathit{L}}_{\mathrm{line}}\mathrm{\propto }{\mathit{r}}_{\mathrm{thick}}^{\mathrm{2}}{\mathit{B}}_{\mathit{\nu }}\mathrm{\left(}\mathit{T}\mathrm{\right(}{\mathit{r}}_{\mathrm{thick}}\mathrm{\left)}\mathrm{\right)}\mathit{.}$(1)It can be shown that r_thick scales roughly as ${\mathit{r}}_{\mathrm{thick}}\mathrm{\propto }{\mathit{L}}_{\mathrm{bol}}^{\mathrm{0.25}\mathrm{-}\mathrm{0.4}}$. For a power-law density structure, $\mathit{T}\mathrm{\propto }{\mathit{L}}_{\mathrm{bol}}^{\mathrm{0.2}}\hspace{0.17em}{\mathit{r}}^{-0.4}$ in the outer envelope. Inserting these relations into the above equation one can show that ${\mathit{L}}_{\mathrm{line}}\mathrm{\propto }{\mathit{L}}_{\mathrm{bol}}^{\mathit{\alpha }}$ with the exponent α in the 0.60–0.82 range. This is entirely consistent with the slope of 0.76 observed for species such as HCO⁺ through its J = 6–5 line (Benz et al. 2015).

This analysis should hold for any molecule or line that is optically thick and whose emission arises primarily from the envelope. The slope depends on where the line becomes optically thick, so it will vary with different molecules and transitions. These arguments could thus explain the strong relation found, for example, for the ¹³CO J = 10–9 lines by San José-García et al. (2013), which are dominated by narrow envelope emission. The somewhat steeper relation of ~1.0 found in that case could be due to additional heating for high-mass sources due to UV radiation, enhancing the temperature.

For ¹²CO J = 10–9, a significant fraction of the emission comes from the outflow rather than the envelope. Table 2 shows that the fraction of envelope emission is rather constant at around 50% across the luminosity range, so that a similar relation could still hold, albeit with extra scatter due to the variation in outflow contribution for each source.

The situation is in principle different for the H₂O lines, which are dominated by outflow material and which are effectively optically thin (see Mottram et al. 2014, and Table 2). In this case, the line luminosity should scale with the total mass of water. The (swept-up) outflow mass as measured from low-J CO lines has been found observationally to scale as ${\mathit{M}}_{\mathrm{outflow}}\mathrm{\propto }{\mathit{L}}_{\mathrm{bol}}^{\mathrm{0.56}}$ (Wu et al. 2005). If we assume that the same relation holds for the shocked gas where H₂O emission originates and if the gas temperature still has a weak dependence on L_bol, one could imagine that various factors conspire to give a slope around 0.75. As argued in Sects. 3.3 and 3.6, the emitting area increases from low- to high-mass YSOs. Table 2 also shows that the fraction of water emission from the envelope increases strongly from low- to high-mass sources, from 0% to close to 40%. This envelope water emission is optically thick and thus the above analysis holds, which would introduce an extra slope on top of the outflow relation.

4.4. From Galactic to extragalactic sources

San José-García et al. (2013) found that the correlation between the logarithms of the ¹²CO J = 10–9 line luminosity and L_bol for the WISH sample held when they extended their analysis to include a sample of extragalactic sources. Similarly, we investigate if the significant correlation measured between the logarithm of the water luminosities, L_H2O, and L_bol continues in the regime of extragalactic sources.

Figure 9 shows the values of L_H2O for a series of AGN and other extragalactic objects presented in Yang et al. (2013). A correlation between the logarithm of these parameters, which extends over more than 12 orders of magnitude in both axis, is measured for the three transitions. The value of the slope is the same for each H₂O line and is equal to 1, which in the case of the 3₁₂–3₀₃ transition is notably larger than the slope obtained considering only the WISH sample. Therefore, if the extragalactic sources are added, the proportionality between the logarithms of the line luminosity and L_bol is emphasised for all water lines and the slope becomes unity within the uncertainties, as for the ¹²CO J = 10–9 analysis of San José-García et al. (2013). As previously discussed in Sect. 4.3 of that paper, this result suggests that both high-J CO and excited water transitions can be used as a tracer of dense gas, similar to HCN (Gao & Solomon 2004a,b) and that the power-law index does not become sub-linear for higher excitation transitions. Indeed, Liu et al. (2015) also show that for ¹²CO J = 10–9 this relation continues to be linear even up to L_bol = 10¹³L_⊙. Regardless of the underlying physical cause, the fact that the line and bolometric luminosities are strongly correlated across such a wide range of L_bol means that water and mid-J CO are good and relatively clean tracers of star formation.

The observed line intensity ratios of the cavity shock component are determined for the three water lines with respect to the 2₀₂–1₁₁ transition (see Fig. 10). The results are presented in the context of the extragalactic framework by comparing the excitation conditions of low-, intermediate- and high-mass YSOs to those investigated in the extragalactic sources discussed in Yang et al. (2013).

Surprisingly, the 752 GHz/988 GHz line intensity ratio only lies above 1 for the high-mass YSOs. This peculiar result, already highlighted within our studied sample of YSOs, also stands out when the line ratios are compared with those measured for extragalactic sources. However, as discussed in Sect. 4.2, the modelling used in this paper is quite simple and is only for lines with E_u ≲ 250 K. Better characterisation of the excitation conditions of more H₂O lines than modelled here is therefore required to conclusively interpret this result.

To summarise, the fact that the slope of the correlation between L_H2O and L_bol (Fig. 9) gets closer to unity, that the significance of the relation improves, and that the dashed line drops below the low-mass Class 0 protostars might suggest a plausible scenario in which the H₂O emission from extragalactic sources comes from more evolved and likely massive sources, and therefore indirectly measures the star formation potential. However, a more fundamental question is whether the extragalactic emission comes from star formation, PDRs, AGN or some combination thereof. Once there is better understanding of where extragalactic water emission originates, further progress can be made on comparing to Galactic observations such as those presented in this paper and on explaining the measured difference between the dashed and solid lines.

5. Conclusions

From the analysis of excited H₂O and high-J¹²CO observations across the entire WISH sample of YSOs we derive the following conclusions:

• The shape of the H₂O 2₀₂–1₁₁, 2₁₁–2₀₂ and 3₁₂–3₀₃ line profiles for a given source are similar but scaled in intensity. In addition, their averaged line intensity ratios are constant with increasing offset from the source velocity. These results suggest that the dynamical properties probed by the three water transitions are similar for a given source.

• On average, more than 60% of the total integrated H₂O line intensity comes from the cavity shock component for all YSOs. The calculated contribution fractions can be used to interpret and compare velocity resolved and unresolved observations of water.

• FWZI and FWHM_b change little between low- and high-mass YSOs for the H₂O spectra, while the FWZIs of the ¹²CO J = 10–9 and 3–2 transitions increase slightly with L_bol. For both molecules, there is significant scatter in the values of FWZI and FWHM_b over the studied luminosity interval, particularly in the case of water (FWZI ranges from 15 to 189 km s^-1).

• This scatter, related to the intrinsic properties of the sources and the varying line-of-sight orientation of the source, is reduced when the ratio of the ¹²CO and H₂O FWZI or FWHM_b were considered. The averaged ratio of their FWZI increases with the source bolometric luminosity ranging from 0.5 to 1.9 and from 0.6 to 0.9 in the case of the FWHM_b. This suggests that for low-mass protostars the entrained material traced by the wings of the ¹²CO J = 10–9 (and lower-J lines) is kinematically different to that traced by water in the cavity shocked gas. In the case of the high-mass YSOs, the CO J = 10–9 emitting region seems to coexist with that of water, so both may come from the same or from closer layers of the outflow cavity wall (Fig. 8).

• The line intensity ratio of ¹²CO J = 10–9 to water decreases with increasing velocity for low-mass YSOs but is approximately constant for high-mass YSOs (Fig. 6, top), confirming the previous conclusion. However, the line intensity ratio of ¹²CO J = 16–15 and H₂O remains constant for both low- and high-mass YSOs as a function of velocity (Fig. 6, bottom). This suggests that higher-J¹²CO transitions (J> 10) trace deeper warmer layers of the outflow cavity characterised by dynamical properties likely linked to those of the shocked water emission, independently of the mass of the object.

• Two scenarios were proposed to explain the previous result. The first scenario assumes that the increased UV radiation in high-mass YSOs leads to water being more effectively photo-dissociated in those layers of the cavity wall closer to the outflow cavity than for the low-mass protostars. An alternative scenario was also discussed, where the stronger turbulence in high-mass star-forming regions, injected by the outflows, propagate into layers of the outflow cavity wall further from the outflow cavity and cause a larger mixing of the material forming those layers than in low-mass protostars. Therefore, the high level of turbulence helps the interaction of warm entrained and cavity shocked material.

• The water excitation for the high-mass YSOs is somewhat different from that of their low- and intermediate-mass counterparts. Since the excitation of the observed H₂O lines is not sensitive to temperature for any of the three sub-types of objects, pumping from an infrared radiation field is the most likely explanation for how the results are reproduced for the majority of the massive YSOs.

• As for CO, the logarithm of the H₂O line luminosity (as well as the emission from the cavity shock component) strongly correlates with the logarithm of the bolometric luminosity and envelope mass of the object with the slope of the correlations being unity within the uncertainty. Although not yet fully understood, the relation seems to imply that as L_bol increases, whatever drives the outflow (wind/jet from the source) can heat material further from the source and/or deeper into the cavity wall, and that the shock cooling structure (and thus n and source-averaged N) is universal. This then means that the emission is always dominated by material at similar conditions but more luminous sources keep more gas at those conditions.

• The correlation between the total H₂O line luminosity can be extended to extragalactic sources. In addition, the values of the line intensity ratios of the extragalactic objects are similar to those calculated for the low-mass protostars, but these ratios differ from those measured for the intermediate- and high-mass YSOs. Thus, both high-J CO and H₂O trace star formation.

Overall, despite some differences in excitation, the properties of excited water line emission are similar between low- and high-mass YSOs, suggesting that a common physical mechanism is at work in the outflows of sources across a wide range of luminosity. High spatial resolution and fidelity imaging of tracers probing both the entrained and currently shocked outflow material will help us to draw a more comprehensive picture of the mechanisms driving the outflows of all YSOs as well and the physical conditions of the material constituting their outflow cavity wall, something ALMA and JWST will enable us to do in the coming years for the first time.

---

Acknowledgments

The authors are grateful to the anonymous referee for the constructive and thorough reviews, to Sylvie Cabrit for the useful feedbacks that helped to improve the manuscript and to Mihkel Kama for providing data and suggestions to the paper. We would like to thank Umut Yıldız for helping on the data reduction as well as the WISH team for many inspiring discussions. This work is supported by the European Community’s Seventh Framework Programme FP7/2007−2013 under grant agreement 238258 (LASSIE) and by the Space Research Organisation Netherlands (SRON). Astrochemistry in Leiden is supported by the Netherlands Research School for Astronomy (NOVA), by a Spinoza grant, A-ERC grant 291141 CHEMPLAN and grant 614.001.008 from the Netherlands Organisation for Scientific Research (NWO). HIFI has been designed and built by a consortium of institutes and university departments from across Europe, Canada and the United States under the leadership of SRON Netherlands Institute for Space Research, Groningen, The Netherlands and with major contributions from Germany, France and the US. Consortium members are: Canada: CSA, U.Waterloo; France: CESR, LAB, LERMA, IRAM; Germany: KOSMA, MPIfR, MPS; Ireland, NUI Maynooth; Italy: ASI, IFSI-INAF, Osservatorio Astrofisico di Arcetri- INAF; Netherlands: SRON, TUD; Poland: CAMK, CBK; Spain: Observatorio Astronómico Nacional (IGN), Centro de Astrobiología (CSIC-INTA). Sweden: Chalmers University of Technology − MC2, RSS & GARD; Onsala Space Observatory; Swedish National Space Board, Stockholm University − Stockholm Observatory; Switzerland: ETH Zurich, FHNW; USA: Caltech, JPL, NHSC.

References

Appendix A: Spectra of the excited water lines

The H₂O 2₀₂–1₁₁, 2₁₁–2₀₂ and 3₁₂–3₀₃ spectra for the intermediate- and high-mass YSOs and the observation number identification are presented in this section (see Figs. A.1 to A.3 and Table A.1). The basic properties derived from these data, such as the rms of the spectra, ${\mathit{T}}_{\mathrm{MB}}^{\mathrm{peak}}$, integrated intensity and FWZI are summarised from Table A.2 to A.4. In addition, the results from the Gaussian decomposition explained in Sect. 2.5 are shown in Tables A.5 to A.7.

Appendix B: Specific sources

The low-mass Class 0 protostars indicated in Fig. 2 are characterised for showing bullet emission in their water profiles. These already studied sources are: L1448-MM, NGC 1333 IRAS2A, BHR 71, Ser SMM1 and L1157. In the case of the low-mass Class I object IRAS 12496, a spot shock component significantly offset from the source velocity is identified in the 988 GHz water line in emission but in absorption in the other ground-based transitions. More information about these objects and their observations in Mottram et al. (2014).

The excited water spectra of the intermediate-mass YSO NGC 2071 show an offset component at around 35 km s^-1 from the source velocity, which could be considered as a spot shock

component, in particular an EHV component. This source is known for being formed by several YSOs, which line profile could be then composed by the emission of several sources and molecular outflows (van Kempen et al. 2012). Therefore, interferometric data is needed disentangle and better understand the emission.

Finally, the line profiles of the high-mass YSO G5.89-0.39 are complex and composed by velocity components with non-Gaussian profiles. As for NGC 2071, this region is crowded with several protostars and energetic outflows (Su et al. 2012) with some of the highest velocities measured for these structures (Choi et al. 1993) Therefore, the interpretation of the emission should go together with extra single-dish and interferometric observations.

Appendix C: Additional figures

All Tables

All Figures