© ESO, 2017

1. Introduction

The general, cartoon picture of how stars form has been agreed for some time: a dense core within a molecular cloud becomes gravitationally unstable, causing material to fall inwards towards the centre; a protostar forms and launches a bi-polar molecular outflow; over time the outflow and infall combine to remove the envelope, eventually starving the protostar, which then slowly settles to the main sequence (e.g. Shu et al. 1987). However, a more detailed understanding is still required, particularly on infall and outflow, in order to quantitatively track the conversion of matter into stars and accurately predict the evolution and outcome of the star-formation process for individual sources, stellar clusters and even whole galaxies.

The first step is improved quantification of the basic physical properties (e.g. L_bol, M_env) and evolutionary state of low-mass protostars, on which considerable progress has been made. Improvements in detectors and telescopes have lead to full-wavelength coverage from optical to radio wavelengths at better sensitivity and resolution, while dedicated very long baseline interferometry (VLBI) campaigns in the radio are providing much more accurate distances for nearby star-forming regions (e.g. Loinard 2013, for a recent review).

A framework for defining the evolutionary status of protostars has also been developed, dividing protostellar sources into one of five categories (Class 0, Class I, Flat, Class II and Class III) using various ways of quantifying the shift in the spectral energy distribution (SED) to shorter wavelengths as the source evolves: the infrared spectral index (α_IR, e.g. Lada & Wilking 1984; Lada 1987; Greene et al. 1994); the submillimetre (λ> 350  μm) to bolometric luminosity ratio (L_submm/L_bol used as a proxy for M_env/L_bol, e.g. André et al. 1993); and bolometric temperature (T_bol, e.g. Myers & Ladd 1993; Chen et al. 1995). For this latter measure, which is the intensity-weighted peak of the SED, these classifications are defined as: Class 0 (T_bol< 70 K), Class I (70 ≤ T_bol< 650 K), Class II (650 ≤ T_bol< 2800 K) and Class III (T_bol ≥ 2800 K). Flat-SED sources have T_bol values in the 350−950 K range with a mean around 650 K (Evans et al. 2009).

The Spitzer Space Telescope (Gallagher et al. 2003) and more recently the Herschel Space Observatory (Pilbratt et al. 2010) have allowed the full potential of this evolutionary framework to be exploited in constraining how the properties of protostars change as the source evolves through large-area, high spatial resolution, uniform photometric surveys of many nearby star-forming regions (e.g. Evans et al. 2003, 2009; André et al. 2010; Rebull et al. 2010; Megeath et al. 2012; Dunham et al. 2014; Furlan et al. 2016). Furthermore, the statistics available from such large surveys have enabled estimates of the relative lifetimes of the different Classes to be obtained, showing in particular that the combined Class 0 and I phases, where the majority of the protostellar mass is accreted and the final properties of the star and its accompanying disk are imprinted, last approximately 0.4−0.7 Myr (Dunham et al. 2015; Heiderman & Evans 2015; Carney et al. 2016).

For a 1 M_⊙ star, such lifetimes imply typical time-averaged mass-accretion rates onto the protostar of approximately 10^-6M_⊙ yr^-1. Since not all material in the core will end up on the star, the infall rate in the envelope must presumably be higher than this by at least a factor of 2 or 3. Searches to quantify the infall in protostars have presented candidates using molecular line observations (e.g. Gregersen et al. 1997; Mardones et al. 1997) based on the doppler-shift of infalling material causing asymmetries in the line profile (Myers et al. 2000). However, confirming and quantifying infall in protostellar envelopes remains extremely challenging, limiting our understanding of the rate at which, and route by which, material reaches the disk and protostar, as well as how this changes with time and depends on the mass of the core/star.

Bipolar molecular outflows also play an important role in the evolution and outcome of star formation, as they remove mass from and inject energy into the envelope and surrounding material. However, the driving mechanism for protostellar outflows is still uncertain (e.g. Arce et al. 2007; Frank et al. 2014). A decrease in the driving force was measured between Class 0 and I sources, in addition to relations with L_bol and M_env, by Bontemps et al. (1996) using ground-based observations of CO. They attributed the decrease in outflow driving force with Class to a decrease in the accretion/infall rate as the source evolves. However, their study only included ten Class 0 sources, as few were known at the time.

Recent observations of H₂O and highly-excited CO using the Heterodyne Instrument for the Far-Infrared (HIFI; de Graauw et al. 2010) and Photodetector Array Camera and Spectrometer (PACS; Poglitsch et al. 2010) with Herschel have shown that these primarily trace active shocks related to the outflow and/or warm disk winds heated by ambipolar diffusion, rather than the entrained outflow as is accessible with ground-based CO observations (Nisini et al. 2010; Kristensen et al. 2013; Tafalla et al. 2013; Santangelo et al. 2013, 2014; Mottram et al. 2014; Yvart et al. 2016). The line-width and intensity in these tracers decreases between Class 0 and I while the excitation conditions (T, N, n) remain the same (Mottram et al. 2014; Manoj et al. 2013; Karska et al. 2013, 2014a; Green et al. 2013a; Matuszak et al. 2015). However, these studies have typically considered relatively small samples (N ≲ 30) of bright, well-known sources and so the statistical significance of trends with evolution and other source parameters has, in some cases, been low.

Two of the main surveys studying nearby Class 0/I protostars with Herschel were the “Water in star-forming regions with Herschel” (WISH) guaranteed time key program (van Dishoeck et al. 2011), which observed 29 Class 0/I protostars with HIFI and PACS plus ground-based follow-up, and the “Dust, Ice, and Gas in Time” (DIGIT) Herschel key program (Green et al. 2013a, 2016), which observed a further 13 Class 0/I protostars, primarily with full-scan PACS spectroscopy. Both the WISH and DIGIT surveys selected their samples to target well known, archetypal sources, ensuring success in detecting water, CO and other species and the availability of complementary data. As a result, these samples favoured luminous sources with particularly prominent and extended outflows, which may not be representative of the general population of protostars. In addition, both programs together only included a total of 42 low-mass sources split between Classes 0 and I, limiting the statistical significance of trends with evolution that might otherwise have been expected, for example between integrated intensity in water emission and T_bol.

The motivation of the William Herschel Line Legacy (WILL) survey was therefore to further explore the physics (primarily infall and outflow) and chemistry of water, CO and other complementary species in Class 0/I protostars in nearby low-mass star forming regions using a combination of Herschel and ground-based observations, building on WISH and DIGIT. The aim was to increase the number of Class 0/I protostars observed, thus improving the statistical significance of the existing correlations found by for example Kristensen et al. (2012), and allowing shallower correlations to be tested, as well as improving the sampling of fainter and colder sources.

This paper is structured as follows. Section 2 discusses the selection of the WILL sample, the basic physical properties of the sources and evaluates the properties of the combined WISH+DIGIT+WILL sample. Section 3 gives the details and basic results of both the Herschel observations and a complementary ground-based follow-up campaign. More detailed results and analysis are then presented thematically, centred around outflows (Sect. 4) and envelope emission (Sect. 5), followed by a discussion on the variation of water with evolution (Sect. 6). Finally, we summarise our main conclusions in Sect. 7.

2. Sample

2.1. Selection

The starting point for selecting a flux-limited sample of low-mass protostars was the catalogue of Class 0/I protostars identified as part of photometric surveys with the Spitzer Space Telescope of the closest major star-forming clouds that make up the Gould Belt (Gould 1879). In particular, these were drawn from the Spitzer c2d (Evans et al. 2009), Spitzer Gould Belt (Dunham et al. 2015) and Taurus Spitzer (Rebull et al. 2010) surveys.

The initial catalogue was compiled from individual cloud catalogues for the Perseus, Taurus, Ophiuchus, Scorpius (also known as Ophiuchus North), Corona Australis and Chameleon star-forming regions (for more details, see Jørgensen et al. 2007; Rebull et al. 2007, 2010; Padgett et al. 2008; Jørgensen et al. 2008; Hatchell et al. 2012; Peterson et al. 2011; Alcalá et al. 2008). At the time of selection in 2011, the Herschel Gould Belt (André et al. 2010) survey had also produced catalogues of protostellar candidates in the Aquila Rift region (Maury et al. 2011), so these were also considered in an attempt to extend the coverage of the WILL survey to particularly young (cold) embedded young stellar objects (YSOs).

From this master catalogue of protostars in major star-forming regions within 500 pc, the following criteria were used to select the final WILL sample:

• (i)

  infrared slope (2−24  μm) α_IR> 0.3 or non-detection;

• (ii)

  T_bol< 350 K;

• (iii)

  L_bol> 0.4 L_⊙ for Class 0 (T_bol< 70 K), L_bol ≥ 1 L_⊙ for Class I (70 ≤ T_bol< 350 K);

• (iv)

  δ< 35°.

The distinction between Class I and II sources is normally made at T_bol = 650 K (Chen et al. 1995), however Evans et al. (2009) found that Flat SED sources cover the range 350−950 K with a mean around 650 K and therefore likely consist of more evolved Class I or younger Class II sources. An upper limit of 350 K was therefore imposed in order to exclude more evolved Class I sources from the sample. Water emission is typically weaker for Class I sources than Class 0s and is generally higher for more luminous sources (e.g. Kristensen et al. 2012), so a higher L_bol cut was used for Class I sources in an attempt to ensure detections. Criteria i−iii were therefore designed to ensure that the sample includes only young, deeply embedded protostars that are bright enough to be detected in H₂O and related species based on the experience of the WISH and DIGIT surveys. Criterion iv ensures that all WILL sources can be observed with ALMA to allow high spectral and spatial resolution ground-based interferometric follow-up of interesting sources.

Unfortunately, edge-on disks, reddened background sources and evolved asymptotic giant-branch (AGB) stars all have the potential to present similar infrared colours and thus contaminate any sample selected purely based on continuum properties. As first highlighted by van Kempen et al. (2009) for a sample of sources in Ophiuchus, molecular emission tracing dense gas can help to break this degeneracy. More specifically, the high critical density of HCO⁺J = 4−3 or J = 3−2 means it will not be strong in foreground cloud material, while the rarity of C¹⁸O similarly means that the J = 3−2 transition is only bright and concentrated in protostellar sources. In addition, more evolved disk sources will not present strong emission in single-dish HCO⁺ spectra due to beam-dilution. Such data, particularly for HCO⁺, have been collected and used to remove contaminants in a number of Gould Belt samples by Heiderman et al. (2010), Heiderman & Evans (2015) and Carney et al. (2016), which have some overlap with the initial candidate sample. Therefore, following the cuts detailed above, non-detection in HCO⁺J = 4−3 or 3−2 was used, where data were available, to exclude contaminant sources.

Most of the sources observed by the WISH and DIGIT surveys also conform to the above criteria, so any initial candidates within 5′′ of a WISH or DIGIT source were also excluded to avoid repeat observations. However, two sources, PER 03 and PER 11, have enough overlap with the WISH observations of L1448-MM (offset by 7.7′′) and NGC 1333-IRAS4B (offset by 6.4′′), respectively, particularly in the H₂O 1₁₀ −1₀₁ (557 GHz) ground-state line obtained in a 39′′ beam, that they are removed from the WILL sample as presented here. Finally, source TAU 05 was removed as it is the young and active Class II source DG Tau B, which has an edge-on disk (Podio et al. 2013).

2.2. Properties and evaluation

The properties of the final sample of 49 sources that make up the WILL sample are presented in Table 1. For simplicity, we give each a name based on the region and a number ordered by right ascension, but many are already well known and therefore the table also gives details of common names used by previous studies for the same sources.

The following distances are used for the various regions covered by our sample: 235 pc for Perseus (Hirota et al. 2008), 140 pc for Taurus (Kenyon et al. 2008), 125 pc for Ophiuchus and Scorpius (de Geus et al. 1989), 130 pc for Corona Australis (Knude & Høg 1998), 150 pc for Chameleon I and 178 pc for Chameleon II (Whittet et al. 1997). For Aquila, W40 and Serpens South, Ortiz-León et al. (2017) recently found that these regions, as well as Serpens Main, are at a common distance of 436 pc.

The determination of the source properties and evolutionary classification is discussed in detail in Appendix A. To summarise briefly, the SED for each source is constructed from the near-IR to (sub-)mm and used to calculate L_bol, L_submm/L_bol, T_bol and α_IR. M_env is obtained from sub-mm or mm photometry assuming that the dust is optically thin, while ν_LSR is calculated from molecular line observations. Finally, the classification of each source is reached by considering the spatial and spectral properties of both the gas and dust associated with each source (see Appendix A.7 for more details).

The sample comprises 23 Class 0, 14 Class I, 8 Class II and 4 uncertain, potentially pre-stellar sources. In the case of this last group of sources, all in W40, they are faint or not detected at <160  μm, show few detections in PACS and have no indications of outflow activity, but the presence of the W40 PDR, detected in some of the HIFI and ground-based lines, leaves some ambiguity. These and other cases of note are discussed in more detail in Appendix C.

Figure 1 shows the L_bol, T_bol and M_env distribution of the WILL sample, along with the WISH and DIGIT samples for comparison. The properties of the WISH sample are taken from Kristensen et al. (2012) while those for the DIGIT sample are taken from Green et al. (2013a) and Lindberg et al. (2014). These are corrected to the distances for the various regions discussed above where needed. It should be noted that M_env values are not available for the DIGIT sample, leading to the difference in the number of sources between the upper-left and upper-right panels.

The probability (p) that a given value of the Pearson correlation coefficient (ρ) for sample size n represents a real correlation (i.e. the likelihood that a two-tailed test can reject the null-hypothesis that the two variables are uncorrelated with ρ = 0) can be expressed in terms of the standard deviation of a normal distribution, σ, as: $\mathit{p}\mathrm{=}\mathrm{\left|}\mathit{\rho }\mathrm{\right|}\sqrt{\mathit{n}\mathrm{-}\mathrm{1}}\mathit{\sigma ,}$(1)following (Marseille et al. 2010). We consider p = 3σ (i.e. 99.7%) to be the threshold for statistical significance. Thus, for a sample size of 30, values of | ρ | > 0.56 indicate real, statistically significant correlations while for a sample size of 50, this is true for | ρ | > 0.43. While one might expect correlations between some of the observed properties of embedded protostars due to the related nature of their different components (e.g. envelope, outflow and driving source), such tests are a simple way of ascertaining whether or not the data are able to support such links. As mentioned above, the extension of the sample of sources studied in spectral lines with PACS and HIFI enabled by the WILL survey and presented here allows us to study these more completely for the first time.

The evolutionary tracks between L_bol and M_env shown in the top-right panel of Fig. 1 are taken from Duarte-Cabral et al. (2013). They assume an exponential decrease of M_env and a core-to-star formation efficiency of 50%, such that the net accretion rate is given by: $\mathit{Ṁ}\mathrm{acc}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{=}\mathrm{0.5}\hspace{0.17em}\frac{{\mathit{M}}_{\mathrm{env}}\mathrm{\left(}\mathit{t}\mathrm{\right)}}{\mathit{\tau }}\mathit{,}$(2)where τ is the e-folding time, which is assumed to be 3 × 10⁵ yr.

The WILL sample doubles the number of low-mass YSOs observed, which have slightly lower values of L_bol and M_env, as well as lower T_bol for Class 0 sources, than the WISH and DIGIT samples. Comparing to Spitzer Gould Belt (SGB) sources with T_bol ≤ 350 K, taken from Dunham et al. (2015), it can be seen in Fig. 1 that the combined WILL+DIGIT+WISH sample is representative of the overall Class 0/I population and contains most sources above ~1 L_⊙. Below this luminosity, the sample rapidly becomes incomplete, and thus the combined sample is still biased towards higher mean L_bol compared with the general distribution, but the addition of the WILL sources shifts the completeness limit approximately a factor of three lower. In terms of T_bol, the sample is biased towards lower values, but judging from upper-left panel of Fig. 1, the higher T_bol sources in the SGB data are primarily those below our L_bol limit, that is, the mean L_bol decreases as T_bol increases for SGB sources. The differences between the values of Dunham et al. (2015) and those given here for individual sources are likely due to our inclusion of far-IR data in these determinations.

It is worth mentioning a couple of caveats. Firstly, the sample of Class 0 sources is dominated by sources in the Perseus molecular cloud, while the Class I sources are drawn from a number of regions that vary in the concentration and activity of their star formation (e.g. Taurus vs. Ophiuchus). There may well be regional differences due to environmental effects, which we cannot test due to the overall small sample size for a given region. Secondly, by excluding older Class I and flat-spectrum sources, we introduce a bias towards younger Class I sources, so the properties of an average Class I source may well be slightly different from those determined with this sample. However, in general for the part of parameter space that WILL, WISH and DIGIT are designed to probe, the addition of the WILL survey leaves the combined sample broadly complete.

3. Observations and results

The primary observations for the WILL survey were taken with Herschel¹ using the Heterodyne Instrument for the Far-Infrared (HIFI, de Graauw et al. 2010) and Photodetector Array Camera and Spectrometer (PACS, Poglitsch et al. 2010) detectors between the 31st October 2012 and 27th March 2013. The observing modes, observational properties, data reduction and detection statistics are described for each instrument separately in Sects. 3.1 and 3.2. Complementary spectroscopic maps obtained through follow-up observations of the sample with ground-based facilities are then described in Sect. 3.3.

3.1. HIFI

3.1.1. Observational details

HIFI was a set of seven single-pixel dual-sideband heterodyne receivers that combined to cover the frequency ranges 480−1250 GHz and 1410−1910 GHz with a sideband ratio of approximately unity. Spectra were simultaneously observed in two polarisations, H and V, which pointed at slightly different positions on the sky (~6.5′′ apart at 557 GHz decreasing to ~2.8′′ at 1153 GHz), with two spectrometers simultaneously providing both wideband (WBS, 4 GHz bandwidth at 1.1 MHz resolution) and high-resolution (HRS, typically 230 MHz bandwidth at 250 kHz resolution) frequency coverage.

The HIFI component of the WILL Herschel observations consists of single pointed spectra at four frequency settings, principally targeting the H₂O 1₁₀ −1₀₁, 1₁₁ −0₀₀ and 2₀₂ −1₁₁ transitions at 557, 1113 and 988 GHz respectively and the ¹²CO J = 10−9 transition at 1152 GHz, which also includes the H₂O 3₁₂ −2₂₁ transition. All observations were carried out in dual-beam-switch mode with a nod of 3′ using fast chopping. The specific central frequencies of the settings were chosen to maximise the number of observable H₂O, CO and H$\begin{array}{c}\mathrm{18}\\ \mathrm{2}\end{array}$O transitions, the details of which are given in Table 2 along with the corresponding instrumental properties, spectral and spatial resolution, and observing time. The main difference compared to the WISH HIFI observations of low-mass sources (see Kristensen et al. 2012; Mottram et al. 2014) was that the frequency of the WILL observations for the H₂O 1₁₀ −1₀₁ and 1₁₁ −0₀₀ settings was set so that the corresponding H$\begin{array}{c}\mathrm{18}\\ \mathrm{2}\end{array}$O transition was observed simultaneously, and longer observing times were used for the H₂O 1₁₀ −1₀₁ setting. The observation ID numbers for all WILL HIFI observations are given in Table B.1.

Initial data reduction was conducted using the Herschel Interactive Processing Environment (hipe v. 10.0, Ott 2010). After initial spectrum formation, any instrumental standing waves were removed. Next, a low-order (≤2) polynomial baseline was subtracted from each sub-band. The fit to the baseline was then used to calculate the continuum level, compensating for the dual-sideband nature of the HIFI detectors (the initial continuum level is the combination of emission from both the upper and lower sideband, which we assume to be equal). Following this the WBS sub-bands were stitched into a continuous spectrum and all data were converted to the T_MB scale using the latest beam efficiencies (see Table 2). Finally, for ease of analysis, all data were converted to FITS format and resampled to 0.3 km s^-1 spectral resolution on the same velocity grid using bespoke python routines.

Few differences have been found in line-shape or gain between the H and V polarisations (e.g. Kristensen et al. 2012; Yıldız et al. 2013; Mottram et al. 2014), so after visual inspection the two polarisations were co-added to improve signal-to-noise. The velocity calibration is better than 100 kHz, while the pointing uncertainty is better than 2′′ and the intensity calibration uncertainty is ≲10% (Mottram et al. 2014).

3.1.2. Results

Figures 2 and 3 present the observed HIFI ortho-H₂O 1₁₀ −1₀₁ (557 GHz) ground-state transition and ¹²CO J = 10−9, respectively, for all WILL sources. The water spectra are complex, containing multiple components, some absorption, which is usually narrow, and emission up to ±~100 km s^-1 from the source velocity, similar to other Herschel HIFI observations of water towards Class 0/I sources (e.g. Kristensen et al. 2012). ¹²CO J = 10−9 typically shows two gaussian emission components with a lower total velocity extent than H₂O. Strong, narrow absorption in ¹²CO J = 10−9 for W40 sources 01, 03 and 06 (see Fig. 3) indicates that contamination in at least one of the reference positions affects these spectra and also likely affects most of the H₂O transitions for these sources as well. The narrow yet bright nature of the ¹²CO J = 10−9 seen in six sources (OPH 01, W40 01 and W40 03−06, see Fig. 3), combined with the narrow and low-intensity nature of the H₂O emission, suggests that they are related to photon-dominated regions (PDRs, cf. for example CO observations of the Orion Bar PDR, Hogerheijde et al. 1995; Jansen et al. 1996; Nagy et al. 2013).

The detection statistics for all transitions are given in the last column of Table 2, excluding W40 sources 01, 03 and 06 due to the contamination of these spectra. The H₂O 1₁₀ −1₀₁ transition is detected towards 39/46 sources in total, including 33/36 confirmed Class 0/I sources (not detected in CHA 02, PER 04 and W40 07, see Fig. 2), while ¹²CO J = 10−9 is detected towards 40/46 sources in total including 32/36 Class 0/Is (not detected in CHA 02, PER 07, PER 15 and W40 07). H$\begin{array}{c}\mathrm{18}\\ \mathrm{2}\end{array}$O 1₁₀ −1₀₁ is only detected towards the source with the strongest H₂O emission, SERS 02, while C¹⁸O J = 9−8 is only detected towards four sources (PER 02, SERS 02, W40 04 and  05).

A more detailed analysis of the kinematics of the HIFI lines is presented and discussed in San José-García (2015), including the results of Gaussian decomposition of the lines using the methods outlined for the WISH sample by Mottram et al. (2014) and San José-García et al. (2013) for H₂O and CO, respectively. In summary, the minimum number of Gaussian components is found that results in no residuals above 3σ, with these components then categorised between the envelope and C or J-type outflow-related shocks depending on their width and offset from the source velocity. A global fit is used for the H₂O transitions with the component peak velocity and line-widths constrained by all lines and the intensity allowed to vary between transitions because the lines all have a consistent shape. The different CO transitions are fit independently as their line profile shapes vary between different transitions.

3.2. PACS

3.2.1. Observational details

PACS consisted of four detectors, two photoconductor arrays with 16 × 25 pixels for integral field unit (IFU) spectroscopy and two bolometer arrays with 16 × 32 and 32 × 64 pixels for broad-band imaging photometry. In IFU spectroscopy mode, observations were taken simultaneously in the red 1st order grating (102−210  μm) and one of the 2nd or 3rd order blue gratings (51−73  μm or 71−105  μm) over 5 × 5 spatial pixels (spaxels), which covered a 47′′ × 47′′ field of view. For details of the 70, 100 and 160  μm PACS (and 250, 350 and 500  μm SPIRE, Griffin et al. 2010) photometric maps used to determine the continuum flux densities for the SEDs (discussed in Sect. A.1) see André et al. (2010).

WILL PACS observations were carried out using the IFU in line-scan mode where deep observations were obtained for targeted wavelength regions (bandwidth Δλ/λ = 0.01) around selected transitions. Two wavelength settings were used, each including observations in both the blue and red gratings, as summarised in Table 3. The principle transitions within these regions are from H₂O, OH, [O i], CO and [C ii], the properties of which are given in Table A.3. While WILL targeted the key lines observed by WISH, some of the wavelength ranges were shifted slightly in order to allow for better baseline subtraction or additional line detections (e.g. those around 82 and 90  μm were shifted to slightly longer wavelengths) while others were omitted to save time (e.g. around CO J = 14–13). The velocity resolution of PACS ranges from 75 km s^-1 at the shortest wavelength to 300 km s^-1, with only [O i] sometimes showing velocity resolved line profiles in a few sources. All observations used a chopping/nodding observing mode with off-positions within 6′ of the target coordinates. The obsids for WILL PACS observations are given in Table B.1. For one source, TAU 08, PACS data were not obtained because the coolant on Herschel ran out before they could be successfully observed.

Data reduction was performed with hipe v.10 with Calibration Tree 45, including spectral flat-fielding (see Herczeg et al. 2012; Green et al. 2013a, for more details). The flux density was normalised to the telescopic background and calibrated using observations of Neptune, resulting in an overall calibration uncertainty in flux densities of approximately 20% (Karska et al. 2014b). 1D spectra were obtained by summing over a number of spaxels chosen after inspection of the 2D spectral maps (Karska et al. 2013), with only the central spaxel used for point-like emission multiplied by the wavelength-dependent instrumental correction factors to account for the PSF (see PACS Observers Manual²).

3.2.2. Results

An overview of the PACS spectra for all sources is shown in Fig. 4, while an overview of the detection of all transitions is given in Table A.4. An extensive analysis of the PACS data for WILL sources in the Perseus molecular cloud was published in Karska et al. (2014b), while a global study of PACS spectroscopy towards all WILL, DIGIT and WISH sources will be presented in Karska et al. (in prep.). Line flux densities were extracted from the PACS data as described in Karska et al. (2013).

The detection statistics for the main transitions are also given in Table A.4. The most frequently detected line is [O i], which is detected in 42 out of 48 sources. Those sources not showing [O i] detections (AQU 03−06, PER 13 and W40 06) are generally not detected in other PACS lines. These sources have weak and/or narrow lines, where detected, in the HIFI observations (cf. Fig. 2). There are 30 sources detected in at least one PACS water transition, while 27 are detected in at least one OH line and 32 in at least one CO line, with a detection more likely in the lower-energy transitions.

3.3. Ground-based follow-up

Follow-up ground-based observations were conducted towards the WILL sample, where not already available, to complement the Herschel spectral line information. Approximately half of the sources in the final catalogue were not part of the samples and regions already observed in HCO⁺J = 4−3 by Carney et al. (2016), so such observations were undertaken to confirm the embedded protostellar nature of the sample (see Appendix A.7). The follow-up observations also included maps of ¹²CO J = 3−2 to characterise the entrained molecular outflow and C¹⁸O J = 3−2 to obtain the source velocity and turbulent line-width in the cold envelope.

All but the two WILL sources in Chameleon are observable from the James Clerk Maxwell Telescope (JCMT³) on Mauna Kea, Hawaii. Observations of C¹⁸O J = 3−2 and HCO⁺J = 4−3 were obtained using HARP (Buckle et al. 2009) and the ACSIS autocorrelator at the JCMT as 2′ × 2′ jiggle maps either as part of observing programs M12AN08 and M12BN07 or from the archive where these were already taken as part of other programs. These also included 2′ × 2′ jiggle map observations of all sources in ¹²CO and H¹³CO⁺J = 4−3, while ¹³CO J = 3−2 was obtained simultaneously with C¹⁸O J = 3−2 for those sources that had not been previously observed. In a few cases, the ¹²CO and¹³CO observations were supplemented with cut-outs from the large basket-woven raster maps taken as part of the JCMT Gould Belt survey observations of Perseus, Taurus and Ophiuchus (Curtis et al. 2010; Davis et al. 2010; White et al. 2015).

For the two Chameleon sources, a series of lines were observed with the Atacama Pathfinder EXperiment (APEX⁴) telescope at Llano de Chajnantor, Chile as part of project M0002_90. These consisted of 2′ × 2′ on-the-fly maps of ¹²CO J = 3−2 and 2−1, as well as single pointings of ¹²CO J = 4−3, ¹³CO, C¹⁸O and C¹⁷O J = 3−2 and 2−1, and HCO⁺ and H¹³CO⁺J = 4−3, using the FLASH⁺ (for the 300 GHz and 450 GHz bands) and APEX1 (for the 225 GHz band, Vassilev et al. 2008) receivers.

The initial reduction of the JCMT jiggle maps was performed using the most up-to-date version of the starlink⁵ reduction package orac-dr (Jenness et al. 2015). Similar initial reduction was performed for the APEX data using gildas-class⁶. Following this, all data were (re-)baselined, corrected to the T_mb scale, and re-sampled to a common velocity scale with 0.2 km s^-1 resolution using customised python scripts. A summary of the observed lines, adopted beam efficiencies and typical σ_rms values obtained is presented in Table A.8. ¹²CO emission is detected towards all sources but not all show evidence of outflows (see Sects. 4.1 and A.4 for more details). More details of detections and non-detections in the ¹³CO, C¹⁸O, C¹⁷O, HCO⁺ and H¹³CO⁺ spectra can be found in Sect. A.6.

4. Outflow characteristics and energetics

In this section we present selected characterisation and comparative analysis of the Herschel and ground-based spectral line observations, focusing on the entrained outflow as probed by ¹²CO J = 3−2 and outflow/wind/jet-related shocks traced by PACS [O i] observations and the broader components of the HIFI H₂O and ¹²CO J = 10−9 lines. In this and the following section, the pre-stellar and Class II sources are excluded from all analyses as they do not show strong outflow or envelope signatures (see Sect. A.7 for characterisation of sources).

Details of how the various entrained outflow-related properties (i.e. mass, momentum, energy, force and mass-loss rate, maximum velocity, dynamical time, inclination and radius) were measured are given in Appendix A.4, along with a table of their values for all WILL sources with detected outflows (Table A.6). For consistency, the calculations are performed following the same method as that used by Yıldız et al. (2015) for the WISH sources, thus ensuring consistency between the WISH and WILL measurements.

4.1. Low-J CO emission

One simple, initial question to ask is whether or not the observations are consistent with the common assumption that all embedded protostars have outflows. Overall 34/37 (92%) of the Class 0/I sources in the WILL sample show outflow emission associated with the source in CO J = 3−2 (shown in Fig. A.2). Two of the Class II sources (CHA 01 and TAU 03) also show outflow activity, which is discussed further in Appendix C. Of the three Class 0/I sources without detections, Per 12, a Class 0 source, shows indications in Spitzer images that the outflow is in the plane of the sky (see Fig. 19 in Tobin et al. 2015, and associated discussion). Cha 02, a Class I, is faint or not detected in most tracers, but is detected in [O i] with PACS. W40 01, a Class 0, is also not detected in most PACS lines but shows a faint broad blue-shifted line-wing in the HIFI spectra. Thus the lack of detection for these three sources is likely due to sensitivity, meaning our observations are consistent with the hypothesis that all Class 0/I sources drive a molecular outflow.

An important next step in understanding the mechanism and impact of outflows on star formation is to constrain how the properties of outflows are related to those of the protostar. Figure 5 shows comparisons of the full-width at zero intensity (FWZI), which is the sum of the maximum velocities in the red and blue outflow lobes, mass in the entrained outflow (M_out), the mass entrainment rate in the outflow (Ṁ_out) and the time-averaged momentum or outflow force (F_CO) measured from CO J = 3−2 with L_bol, T_bol and M_env for all WISH and WILL Class 0/I sources. Ṁ_out and F_CO, as calculated quantities, are corrected for the inclination (see Appendix A.4 for details), while we do not correct for inclination for directly measured quantities, such as FWZI and M_out.

The strongest correlation (5.9σ) is between the outflow mass and envelope mass, perhaps unsurprisingly given that the outflow is entrained from the envelope, with M_out/M_env centred around 1%, as shown by the grey dashed line. A significant, though weaker, correlation is also found between outflow mass and luminosity (3.4σ), possibly due to the correlation between L_bol and M_env (see Fig. 1). Since the mass is also an important factor in the calculation of Ṁ_out and F_CO, it is not surprising that both also show significant (i.e. ≥3σ) correlations with M_env and L_bol. In general, all parameters decrease between Class 0 and Class I, as reflected in the significant negative correlations with T_bol seen for FWZI, M_out, Ṁ_out and F_CO (3.9, 3.8, 3.6 and 3.1σ respectively).

Correlations of F_CO with L_bol and M_env have been known for some time, including Cabrit & Bertout (1992) who found the relationship between F_CO and L_bol indicated by the dot-dashed line in Fig. 5 for a sample of Class 0 sources, and Bontemps et al. (1996) who found the relationships indicated by the solid lines for a sample of primarily Class I sources. These have subsequently been confirmed to hold when extended to the high-mass regime for a sample of young protostars in Cygnus-X by Duarte-Cabral et al. (2013) and for a sample of massive young stellar objects (MYSOs) and young H ii regions by Maud et al. (2015). The relationships seen between these variables in the combined WISH and WILL sample are steeper than that found by Cabrit & Bertout (1992) and Bontemps et al. (1996). This may be due to differences in the calculation method (van der Marel et al. 2013), or to the fact that the luminosities were likely overestimated and the F_CO values underestimated in these previous studies due to the larger beam and lower sensitivity of older observations.

At first glance, the lower right panel of Fig. 5 would seem to show a slight offset in the F_CO measurements between the WISH and WILL samples, suggesting that either there is a difference in the measurements or that they come from distinct populations. However, there is no distinct break between the WISH and WILL sources, or between Class 0 and I when considering M_env vs. M_out and Ṁ_out, so the WISH sources are merely the extreme upper end of a continuous distribution. The WISH Class 0 sources were all chosen to be strong outflow sources, and the M_out and Ṁ_out panels of Fig. 5 suggest that they are more prominent due to a larger reservoir of material (i.e. larger M_env), rather than faster outflows as they have similar or even lower FWZI than the WILL sources.

For the Class I sources, there is little difference between the WISH and WILL sources in FWZI or M_out, but the WISH Class I sources tend to have smaller outflows (see Table A.6 for WILL sources and Table 3 in Yıldız et al. 2015, for WISH sources) resulting in larger values for the WISH sources of Ṁ_out and F_CO. This could be because the WISH Class I sources are typically in smaller, more isolated clouds with shorter distances from the protostar to the cloud edge than the Class I sources in WILL.

Let us now consider the physical implications of the main correlations between outflow and source properties. The correlation between F_CO and (current) M_env is often interpreted as the result of an underlying link between envelope mass and the mass accretion rate (Ṁ_acc), which is itself related to the driving of the outflow (Bontemps et al. 1996; Duarte-Cabral et al. 2013). As the central source evolves, M_env and Ṁ_acc decrease, leading naturally to the decrease in F_CO and other outflow-related properties between Class 0 and I sources. The comparatively tight relationship between M_out and M_env further supports this interpretation.

Indeed, the relation between M_out and M_env requires more investigation in its own right. Figure 6 shows a histogram of the fraction of mass in the outflow compared to the envelope (i.e. M_out/M_env), as well as how this varies with T_bol (as a more continuous proxy for source evolution), the mean length of the outflow lobes (R_out), and the strength of the outflow as measured by F_CO. The values of M_out/M_env vary between ~0.1 and 10%, peaking around 1%. The peak is similar between Class 0 and I with no significant trend with T_bol, except that the Class 0 sources extend to larger values. This seems to be related to some Class 0 sources having longer outflows (i.e. larger R_out) and thus have likely entrained additional material from the clump/cloud outside their original envelope. A statistically significant (3.7σ) correlation with F_CO is found, though with more than an order of magnitude spread.

The first impression of the peak value of M_out/M_env being ~10^-2 in the histogram shown in Fig. 6 might be that this is rather low compared to a “typical” star formation efficiency of 30−50% (e.g. Myers 2008; Offner et al. 2014; Frank et al. 2014, and references therein). In order to understand whether or not this value is actually reasonable, let us first assume that the outflow is responsible for removing all of the envelope material that does not end up on the star. In this case, the average mass entrainment rate in the outflow over the Class 0+I lifetime (τ_{Class 0/I}) is given by: $\mathit{Ṁ}\mathrm{out}\mathrm{=}\mathrm{(}\mathrm{1}\mathrm{-}{\mathit{ϵ}}_{\mathrm{sf}}\mathrm{)}\hspace{0.17em}\frac{{\mathit{M}}_{\mathrm{env}}}{{\mathit{\tau }}_{\mathrm{Class}\hspace{0.17em}\mathrm{0}\mathrm{+}\mathrm{I}}}\mathit{,}$(3)where ϵ_sf is the core-to-star formation efficiency, that is, the fraction of the envelope that will end up on the star. The observed mass-loss rate in the outflow, averaged along the flow, is given by: ${\mathit{M}}_{\mathrm{out}}\mathrm{=}\mathit{Ṁ}\mathrm{out}\hspace{0.17em}{\mathit{t}}_{\mathrm{dyn}}\mathit{,}$(4)where t_dyn is the dynamical time of the flow.

It is worth pointing out that t_dyn is not necessarily the age of the source, particularly if outflow activity is time-variable. Indeed, if ejection stops then, after some time, radiative losses and mixing with the ambient cloud material will dissipate all the angular momentum and energy from the flow, meaning that the observed t_dyn is likely a lower limit to the true “age” of total accretion/outflow activity in a given source. However, if we assume that the overall mass outflow rate for a given burst is not significantly different from the average over the lifetime of the main accretion (i.e. Class 0+I) phase, or equally that protostellar outflows have an approximately constant entrainment efficiency per unit length, then we can combine and re-arrange Eqs. (3) and (4) to get: $\frac{{\mathit{M}}_{\mathrm{out}}}{{\mathit{M}}_{\mathrm{env}}}\mathrm{=}\mathrm{(}\mathrm{1}\mathrm{-}{\mathit{ϵ}}_{\mathrm{sf}}\mathrm{)}\frac{{\mathit{t}}_{\mathrm{dyn}}}{{\mathit{\tau }}_{\mathrm{Class}\hspace{0.17em}\mathrm{0}\mathrm{+}\mathrm{I}}}\mathrm{·}$(5)The ratio of t_dyn to τ_{Class 0 + I} effectively expresses the duty cycle of the outflow.

For τ_{Class 0/I} ≈ 0.5 Myr (Dunham et al. 2015; Heiderman & Evans 2015; Carney et al. 2016) and a typical dynamical time for the outflow of approximately 10⁴ yr, the ratio of outflow to envelope mass has a value of ~0.01 if the ϵ_sf = 0.5. Thus, while certainly missing some details, and being affected by variation from source to source and with time, the fact that we find median value for M_out/M_env of approximately 1% is consistent with protostellar outflows having an approximately constant entrainment efficiency per unit length, a core-to-star formation efficiency of approximately 50%, and an outflow duty cycle of order ~5%.

4.2. HIFI water and mid-J CO emission

The water and CO J = 10−9 spectra of many of the WILL sources show broad line wings, indicative of outflow emission, as seen in previous WISH observations (see e.g. Fig. 2 and Kristensen et al. 2012; San José-García et al. 2013, 2016; Mottram et al. 2014). The FWZI of H₂O is strongly correlated with CO J = 10−9 and J = 3−2 (see Fig. 7), with H₂O consistently tracing faster material than these CO transitions, typically by a factor of ~2. The narrower line-widths of CO J = 10−9 with respect to CO J = 3−2 are likely because the CO J = 3−2 FWZI is calculated as the difference of the maximum red and blue velocity offsets anywhere in the 2′ × 2′ maps while for CO J = 10−9 this is measured from a single HIFI spectrum with a 18.4′′ beam centred at the source position, that is, over a smaller region.

Figure 8 shows a comparison of various source and outflow-related properties with the integrated intensity of the H₂O 1₁₀ −1₀₁ (557 GHz) line after scaling to a common distance of 200 pc. A linear scaling is used because the emission is dominated by outflows, which likely fill the beam along the outflow axis but not perpendicular to it (see Mottram et al. 2014,for more details). We are able to confirm the strong correlation found by Kristensen et al. (2012) for the WISH sample alone between the integrated intensity of the water line with its FWZI (at 6.0σ) and M_env (4.9σ). A new correlation is also found with F_CO (4.8σ), firmly showing that water emission is related to, though not tracing the same material as, the entrained molecular outflow. Furthermore, shallower trends of water line intensity with L_bol and inversely with T_bol, hinted at but not significant in the WISH sample (e.g. see Fig. 6 of Kristensen et al. 2012), are now confirmed as statistically significant at 3.7 and 3.3σ respectively.

To further examine the variation of water emission with source evolution, the lower-right panel of Fig. 8 shows the integrated intensity normalised by the source bolometric luminosity (thus minimising the contribution due to source brightness) vs. ${\mathit{M}}_{\mathrm{env}}\mathit{/}{\mathit{L}}_{\mathrm{bol}}^{\mathrm{0.6}}$, which was proposed by Bontemps et al. (1996) as an evolutionary indicator. The clear positive correlation (4.1σ) seen in this panel reinforces the finding that the intensity of water emission decreases as sources evolve, independent of the relationship between integrated intensity and L_bol.

The WILL observations therefore reinforce and confirm the results from WISH: H₂O traces a warmer and faster component of protostellar outflow than the cold entrained molecular outflowing material traced by low-J CO (e.g. Nisini et al. 2010; Kristensen et al. 2012; Karska et al. 2013; Santangelo et al. 2013; Busquet et al. 2014; Mottram et al. 2014; San José-García et al. 2016). In addition, they confirm that the intensity of H₂O is related to envelope mass and the strength of the entrained molecular outflow, and is higher for younger and/or more luminous sources.

4.3. [O I] emission

It has been suggested for some time that emission from [O i] is a good alternative tracer of the mass loss from protostellar systems (e.g. Hollenbach 1985; Giannini et al. 2001). In Class 0/I protostars it is thought to primarily trace the atomic/ionised wind, because most PACS observations are spectrally unresolved and those few sources that do show velocity-resolved emission (e.g. see Nisini et al. 2015) are dominated by the unresolved (≲100 km s^-1) component. While there may be a contribution on-source from the disk, as in more evolved sources (see e.g Howard et al. 2013), [O i] emission in Class 0/I sources is often spatially extended and only fainter off-source by a factor of ~2 compared to the peak position, so the wind likely still dominates.

The first comprehensive surveys of the [O i] 63  μm transition towards samples of YSOs, observed in an 80′′ beam with the Infrared Space Observatory Long Wavelength Spectrometer (ISO-LWS, Swinyard et al. 1998), suggested a link between the mass loss in [O i] and that in CO for Class 0 sources (Giannini et al. 2001). Only a marginal difference was seen in [O i] luminosity between Class 0 and I sources (Nisini et al. 2002), in contrast to the trend in CO. More recent studies by Podio et al. (2012) and Watson et al. (2016) with PACS on Herschel at 9′′ resolution have used [O i] observations to claim trends of decreasing mass loss in the wind between Classes 0, I and II. However, both suffered from low number statistics, and Podio et al. (2012) mixed the same ISO results where no trend was found with early detections from Herschel PACS, which have significantly different beam sizes and observing methods that could induce such changes. For example, the chopping as part of PACS observations can cancel out up to 80% of the large-scale emission that is still detected by ISO (see Appendix E of Karska et al. 2013). The combined WILL, WISH and DIGIT dataset, with consistent observations of a large number of YSOs is ideally placed to help solve this issue.

Figure 9 shows the distribution of [O i] luminosity in the 63  μm line, integrated over the PACS spaxels associated with source outflows, and how this varies with various source parameters for the WILL, WISH and DIGIT samples (see Table A.7 for L[O i] values). Also shown are the measurements from Herschel studies of a number of Class I/II sources in Taurus (Podio et al. 2012), the “Herschel Orion Protostars” survey (HOPS Watson et al. 2016), and the “FU Orionis Objects Surveyed with Herschel” survey (FOOSH Green et al. 2013b) which targeted a number of Flat spectrum and Class II sources that show evidence of FU Ori-type luminosity outbursts. None of the detected sources have a line luminosity below the upper limit for disk sources found by Howard et al. (2013) towards sources in Taurus (4 × 10¹⁷ W m^-2, corresponding to ~2 × 10^-5L_⊙ assuming a distance of 140 pc).

Two primary results stand out from Fig. 9. First, L[O i] is strongly correlated with L_bol but not with M_env, with sources of all evolutionary classification following the overall trend. This is essentially the reverse of the situation found with low-J CO, where the correlation is weak with L_bol and strong with M_env (see e.g. Fig. 5). H₂O shows clear correlations with both M_env and L_bol (see 8), though the relationship is slightly stronger with M_env than L_bol, consistent with it tracing actively shocked outflow material between the entrained outflow, probed by low-J CO, and the wind, probed by [O i].

Second, there is no statistically significant variation in L[O i] with evolutionary stage, either when considering the flat distribution between L[O i] and T_bol or the histogram of L[O i], which shows remarkably similar distributions for Class 0, I or II sources. This is not due to an evolutionary trend being masked by the correlation with L_bol, as shown by the flat distribution in L[O i]/L_bol vs. ${\mathit{M}}_{\mathrm{env}}\mathit{/}{\mathit{L}}_{\mathrm{bol}}^{\mathrm{0.6}}$. There is also no statistically significant correlation with envelope mass or integrated intensity in the H₂O 1₁₀ −1₀₁, which is dominated by the fast, actively shocked component of the molecular outflow.

This apparent contradiction between the evolutionary behaviour of mass-loss indicators, that is, the decrease of CO and H₂O velocity, intensity etc. as sources evolve compared to the invariance of [O i], will be explored and discussed in more detail in the following subsections. It is interesting to note that the FOOSH sources, which are all known to be undergoing luminosity outbursts, are on the upper end of, but consistent with, the distribution of other sources in terms of L[O i] vs. L_bol. Thus, [O i] must react relatively quickly to variations in the mass accretion rate, which has a significant contribution to the observed source luminosity.

4.4. Mass accretion vs. loss

The balance of mass loss vs. accretion is important in revealing the rate at which the central protostar gains mass, as well as what fraction of the initial envelope will become part of the central source, that is, the core to star efficiency of star formation.

Direct measurement of the mass accretion rate is extremely challenging for embedded protostars because the UV, optical and near-IR continuum and lines typically used to do this in more evolved T-Tauri stars (e.g. Ingleby et al. 2013) are too heavily extincted. An approximate estimate can be obtained, however, by rearranging the equation for accretion luminosity, that is, $\mathit{Ṁ}\mathrm{acc}\mathrm{=}\frac{{\mathit{L}}_{\mathrm{acc}}{\mathit{R}}_{\mathrm{\ast }}}{\mathit{G}{\mathit{M}}_{\mathrm{\ast }}}\mathit{,}$(6)with the aid of a number of empirically constrained assumptions.

Firstly, accretion is assumed to generate all the observed bolometric luminosity for Class 0 sources and 50% for Class I sources, in keeping with the range observed in the few cases where this could be measured (L_acc/L_bol = 0.1−0.8: Nisini et al. 2005; Antoniucci et al. 2008; Caratti o Garatti et al. 2012). Next, a typical stellar mass (M_∗) of 0.5 M_⊙ for Class I sources is assumed and 0.2 M_⊙ for Class 0 sources as they are still gaining mass, as in Nisini et al. (2015). The chosen values are for sources that will end up slightly more massive than the peak of the IMF (~0.2 M_⊙, Chabrier 2005). However, as already discussed in Sect. 2.2, our sample is biased towards slightly higher luminosities, and thus presumably stellar masses, than the global distribution, so this assumption is probably not far off. Indeed, these stellar masses are broadly in keeping with several recent mass determinations for similar embedded protostars from disk studies (Tobin et al. 2012; Murillo et al. 2013; Harsono et al. 2014; Codella et al. 2014). Finally, we assume a stellar radius (R_∗) of 4 R_⊙. The calculated values are given in Table A.7 and shown vs. L_bol, T_bol and M_env in the top panels of Fig. 10. The solid line in the upper-left panel shows the relation assumed in the evolutionary models of Duarte-Cabral et al. (2013), that is, Eq. (2) with τ = 3 × 10⁵ yr.

Hollenbach (1985) noted a simple scaling between the [O i] line luminosity at 63  μm and the total mass-flux through the dissociative shock(s) producing it, given by: $\mathit{Ṁ}\mathrm{s}\mathrm{=}{\mathrm{10}}^{-4}\mathit{L}\mathrm{[}\mathrm{O}\hspace{0.17em}{i}\hspace{0.17em}\mathrm{63}\mathit{\mu }\mathrm{m}\mathrm{]}\mathit{.}$(7)For shocks generated by the wind, as is most likely the case for the emission probed by [O i], the mass flux through the shock(s), Ṁ_s, is related to the wind mass-flux, Ṁ_w, by the general formula (see Dougados et al. 2010): $\mathit{Ṁ}\mathrm{s}\mathrm{=}{\mathit{N}}_{\mathrm{s}}\hspace{0.17em}\frac{{\nu }_{\mathrm{s}}}{{\nu }_{\mathrm{w}}\cos \mathrm{\left(}\mathit{\theta }\mathrm{\right)}}\mathit{Ṁ}\mathrm{w}\mathit{,}$(8)where N_s is the number of shocks in the beam, ν_s is the shock speed, and θ is the angle between the normal to the shock front and the wind direction (the 1/cos(θ) term then accounts for the ratio of the shock area to the wind cross section). It may be seen that Ṁ_s = Ṁ_w in the simple case considered by Hollenbach (1985) if we are observing a static terminal shock where the wind is stopped against a much denser ambient medium; in this case, N_s = 1 and ν_s = ν_wcos(θ). This remains valid if the wind is not isotropic but collimated into a jet.

If we are instead observing weaker internal shocks travelling along the jet/wind, then ν_s ≪ ν_wcos(θ) but this will tend to be compensated for by the presence of several shocks in the beam (i.e. N_s> 1), as suggested by the chains of closely spaced internal knots seen in optical jets.

An alternative method for obtaining the average Ṁ_w in this case is to consider that the [O i] emission is approximately uniform along the flow within the aperture, and to divide the emitting gas mass by the aperture crossing time. The derivation of emitting mass requires assumptions on the temperature and electron density, which are somewhat uncertain without also having observations of the [O i] 145  μm transition. However, Nisini et al. (2015) found that the differences are small between this alternative per-unit-length calculation and the Hollenbach (1985) formulation for a terminal static wind shock (i.e. Ṁ_s = Ṁ_w). Hence, although we note that there are some uncertainties involved, we adopt Ṁ_w[O i] = Ṁ_s as given by Eq. (7) to estimate the wind mass-flux from L[O i] for our targets. The calculated values are given in Table A.7.

The ratio of the mass-loss rate in the wind as measured from [O i] using Eq. (7) to the mass accretion rate (i.e. Ṁ_w[O i]/Ṁ_acc) is compared to L_bol, T_bol and M_env in the middle panels of Fig. 10. Ṁ_w[O i]/Ṁ_acc varies from approximately 0.1% to 100% with a median of 13%, in agreement with previous determinations (e.g. Cabrit 2009; Ellerbroek et al. 2013) and in line with theoretical predictions (e.g. Konigl & Pudritz 2000; Ferreira et al. 2006). However, approximately two-thirds of all Class 0 sources lie below 10%.

The lower panels of Fig. 10 show similar comparison using the outflow force as measured from CO J = 3−2. Assuming the entrainment process is momentum conserving: ${\mathit{F}}_{\mathrm{CO}}\mathrm{=}\mathit{Ṁ}\mathrm{w}\hspace{0.17em}{\nu }_{\mathrm{w}}{\mathit{ϵ}}_{\mathrm{ent}}\mathit{,}$(9)where ϵ_ent is the entrainment efficiency. The ratio with the mass accretion rate is then: $\frac{{\mathit{F}}_{\mathrm{CO}}}{\mathit{Ṁ}\mathrm{acc}}\mathrm{=}\frac{\mathit{Ṁ}\mathrm{w}}{\mathit{Ṁ}\mathrm{acc}}\hspace{0.17em}{\nu }_{\mathrm{w}}{\mathit{ϵ}}_{\mathrm{ent}}\mathit{.}$(10)$\frac{\mathit{Ṁ}{}_{\mathrm{w}}}{\mathit{Ṁ}{}_{\mathrm{acc}}}\hspace{0.17em}{\nu }_{\mathrm{w}}$ is expected to be approximately constant due to conservation of angular momentum, with a value close to the Keplerian velocity of the disk at the launching radius (Duarte-Cabral et al. 2013).

We find that F_CO/Ṁ_acc is relatively invariant with L_bol, T_bol and M_env, as shown by values of the Pearson coefficient ρ consistent with 0 (i.e. p< 3σ). Taken together, this suggests that the efficiency of entrainment, ϵ_ent, is not dependent on source properties. The Keplerian velocity for a disk around a 0.2 or 0.5 M_⊙ source is approximately 10−20 km s^-1 at 1 AU, which, for a median value of F_CO/Ṁ_acc of 6.3 km s^-1, suggests values for ϵ_ent of approximately 0.3−0.6. If the wind is launched at larger radii then ϵ_ent could be closer to 1.

L[O i] does not vary with T_bol, M_env or evolutionary stage (see Sect. 4.3 and Fig. 9), so the increase of Ṁ_w[O i]/Ṁ_acc between Class 0 and I with increasing T_bol and with decreasing M_env is caused by the decrease in Ṁ_acc , while Ṁ_w[O i] remains relatively constant. In contrast, the invariance of F_CO/Ṁ_acc is caused by both Ṁ_acc and Ṁ_w decreasing with increasing T_bol and decreasing M_env (see the lower panels of Fig. 5 for variation of F_CO with T_bol and M_env). The reason for the difference in behaviour between these two measures of the ratio of mass loss to mass accretion is discussed further in the following section.

4.5. On the difference between [O I] and CO

The difference in behaviour between the atomic component of the wind (as traced by [O i]) and the entrained molecular outflow (as traced by low-J CO) might seem to be in contradiction with models where the wind is the driving agent of the outflow (see e.g. Arce et al. 2007). Indeed a direct comparison, shown in Fig. 11, suggests that either the wind and outflow are not linked, [O i] is under-estimating the mass loss rate in the wind or F_CO is overestimated. However, there are several factors relating to what component of the system each tracer probes that argue against rushing to such a conclusion.

First, [O i] only traces the atomic component of the wind and/or jet. Jets in Class 0 protostars are known to have a significant molecular component, as identified from high-velocity features (detected in e.g. CO, SiO and/or H₂O, Bachiller et al. 1991; Tafalla et al. 2010; Kristensen et al. 2011) with typical mass-loss rates of approximately 10^-7−10^-5M_⊙ yr^-1 (e.g. Santiago-García et al. 2009; Lee et al. 2010a). These are typically approximately ten times higher than measured from [O i], but the molecular jet component disappears in older sources. This suggests an evolution in composition from molecular to atomic/ionised (see Nisini et al. 2015, for a detailed discussion), most likely due to increasing temperature of the protostar and decreasing density, and thus shielding, in the jet. Such arguments also hold for any wide-angle wind that could be present and contributing to driving the entrained CO outflow. Therefore, while the mass-loss rate due to the wind as a whole will decrease as the source evolves, in line with the decrease in the average mass accretion rate, the mass loss in the atomic component may remain approximately constant due to the shift in the composition of the wind.

Next, the optical depth of the continuum at 63  μm is likely considerable in the inner envelope in Class 0 sources (see e.g. Kristensen et al. 2012), so the observed [O i] flux may be significantly lower than the “true” emission. The continuum optical depth will decrease as the source evolves and M_env decreases, which may also act to counteract the evolution in the mass loss in the wind. However, such an effect should also cause the ratio of the 63  μm to 145  μm [O i] lines to vary with continuum optical depth of the source, and a wavelength-dependent deficit in CO and H₂O transitions. Neither is clearly seen in PACS observations (see e.g. Karska et al. 2013). This is therefore likely a minor effect dominating only for sight-lines directly towards the protostar through the disk.

Finally, there is increasing evidence that episodic or time-variable accretion is important in embedded protostars from the very earliest phases of their evolution (see Dunham et al. 2014; Audard et al. 2014, for recent reviews). Accretion variability provides a consistent explanation for very low luminosity objects (e.g. Dunham et al. 2006), the observed spread and trends in protostellar (e.g. Dunham et al. 2010) and outflow related (e.g. Duarte-Cabral et al. 2013) properties, and luminosity bursts, brighter by at least a factor of ten, have now been observed in at least two embedded sources (Caratti o Garatti et al. 2011; Fischer et al. 2012; Safron et al. 2015). Chains of high-velocity molecular knots or “bullets” observed in Class 0 outflows and jets, with typical spacings of 1000−10 000 AU between minor and major episodes, respectively (e.g. Santiago-García et al. 2009; Lee et al. 2015), past heating of CO₂ ice (e.g. Kim et al. 2012) and the difference between the expected and observed CO snow surface in a number of protostars (Visser et al. 2015; Jørgensen et al. 2015) also provide indirect evidence of outbursts.

The imprint of time-variable accretion will be different for the molecular outflow and atomic wind, leading to differences in their properties. The luminosity will react quickly to any changes in the accretion rate (Johnstone et al. 2013), and thus traces the current or instantaneous activity. Since [O i] is dominated by the wind, it traces material that is closely related to the current accretion state and thus is correlated with luminosity regardless of whether the source is in outburst (e.g. the FOOSH sources) or not (see Fig. 9).

In contrast, the entrained molecular outflow traced by low-J CO, particularly when measured over the full extent of the outflow, is an average of the ejection activity over at least 10³ −10⁵ yr. Indeed, the highest intensity in the entrained outflow as traced by low-J CO is usually offset from the central source. If these spots represent major ejections triggered by accretion bursts, then such episodes should have occurred approximately hundreds to thousands of years ago, leaving enough time for the luminosity and circumstellar material to cool back to pre-burst levels (e.g. Arce & Goodman 2001; Arce et al. 2013). Thus, the mass loss in the molecular outflow is related to the time-averaged mass-accretion rate and may be dominated by any periods of high accretion/ejection during outbursts (see also e.g. Dunham et al. 2006; Lee et al. 2010b).

The decrease of F_CO between Class 0 and I (see Fig. 5) therefore shows that the average mass accretion rate decreases as sources evolve, as originally proposed by Bontemps et al. (1996). The combination of decreasing mass accretion rates and episodic accretion was shown by Duarte-Cabral et al. (2013) to be consistent with the observed relationships between, and spread of, L_bol, M_env and F_CO. In particular, variation of the mass-accretion rate on shorter timescales than the dynamical timescale of the outflow helps to explain why outflow properties are less correlated with L_bol than with M_env (see Fig. 5 and Sect. 4.1). Those sources that show particularly high outflow forces and/or peak emission close to the source position may therefore have recently finished such a burst, or have a higher duty cycle of outburst to quiescent accretion. Thus, the mass-loss rate in the wind from [O i] and in the outflow force measured by low-J CO are not directly correlated because the relationship between the current and time-averaged mass-accretion rate will be different for each source based on a complex combination of the source age, properties and mass-accretion history.

Some combination of the effects discussed above therefore causes the observed lack of correlation between CO and [O i]. As such, [O i] is not necessarily a direct alternative to CO for tracing mass loss and/or entrainment due to the jet/wind/outflow system in protostars, in contradiction to the early findings of Giannini et al. (2001).

5. Envelope

5.1. HCO⁺ vs. H₂O

Instead of outflows, HCO⁺J = 4−3 emission primarily traces cool, high-density envelope material (T ~ 40 K, n_cr = 2 × 10⁷ cm^-3 though the effective density for optically thick emission could be as low as 10⁴ cm^-3, see Shirley 2015), and so is a relatively clean discriminator between young, embedded protostars and pre-stellar or more evolved disk sources (van Kempen et al. 2009; Carney et al. 2016). Spatially compact detections are found in this line towards most of the WILL sample sources, confirming them to be genuine embedded Class 0/I sources, while pre-stellar and Class II sources are either non-detections or show extended emission with no clear peak at the source position (see Appendix A.7 for details).

H₂O emission is a good tracer of warm, relatively dense material in shocks related to protostellar outflows (T ≳ 300 K, n = 10⁵−10⁸ cm^-3, Kristensen et al. 2013; Mottram et al. 2014). Sources with higher luminosities typically have stronger outflows (and thus stronger H₂O emission) and will lead, all other things being equal, to more mass at a given temperature in their envelopes and thus higher intensity in molecular tracers such as HCO⁺. It is therefore not unreasonable to expect that the emission in these two tracers may be related in Class 0/I sources.

Carney et al. (2016) compared the nature of HCO⁺ emission (compact, confused or extended/not detected) and the detection of water in the WILL PACS 179  μm observations. There are 18 sources in common between their sample and WILL: 13 are classified as Class 0/I (i.e. compact HCO⁺J = 4−3), 3 as confused and 2 as Class II (i.e. extended and/or non-detections in HCO⁺). Most sources (14/18) have detected 179  μm water emission, with four Class 0/I sources and one confused source showing extended emission. As both Class II sources were detected, while three Class 0/I sources and one confused source was not, Carney et al. (2016) did not find a clear relationship with evolution between the spatial distributions and detection of HCO⁺ and H₂O.

Considering the H₂O 1₁₀ −1₀₁ HIFI observations (see Fig. 12), for the two Class II sources, TAU 07 is not detected at the 3σ level and TAU 09 shows very weak, narrow emission with an integrated intensity of 0.16 K km s^-1, equivalent to 10^-18 W m^-2, which could include a contribution from the disk (cf. Podio et al. 2012, 2013; Fedele et al. 2013). Thus, while H₂O may be detected in either Class 0/I or II, the origin and intensity of the water emission changes as the source evolves. Strong detection of either HCO⁺ or water is therefore still a good indication of the youth of a protostar.

For the Class 0/I sources, the full WILL dataset enables the relationship between these species to be probed further by comparing the integrated intensity and FWHM of HCO⁺ with the integrated intensity and FWZI of the H₂O 1₁₀ −1₀₁ HIFI observations (see Fig. 12). HCO⁺ intensity is correlated both with L_bol and the intensity of the water line at 4.2 and 3.5σ significance, respectively. There is not a strong, statistically significant relationship between the kinematics of the two lines and no significant line-wings are seen in the HCO⁺ spectra (see Fig. A.3). HCO⁺J = 4−3 therefore seems to primarily trace parts of the envelope that are further from, and thus less disturbed by, the outflow. HCO⁺ can be destroyed through reactions with H₂O (e.g. Jørgensen et al. 2013), so the higher abundance of H₂O in the outflow (X [ H₂O ] or approximately 10^-5 −10^-7: Tafalla et al. 2013; Santangelo et al. 2013; Kristensen et al. 2017) compared to the envelope (10^-8 −10^-11: Mottram et al. 2013; Schmalzl et al. 2014) could be suppressing the HCO⁺ abundance in the outflow. This would explain why it is a poorer tracer of the outflow than might otherwise be expected. The correlation between the intensity of the two lines is therefore likely due to the relation between emission in each line and the source luminosity and structure, assisted by their tracing similar densities.

5.2. Infall signatures

The bulk of the H₂O emission comes from outflow-related shocks that have Gaussian-like profiles in velocity-resolved spectra. Once this contribution is removed, the residual profiles show the remaining water emission and/or absorption associated with the envelope. For the WISH sample, this process revealed seven sources with inverse P-Cygni line profiles indicative of infall and five with regular P-Cygni line profiles indicative of expansion motions in the envelope (Kristensen et al. 2012; Mottram et al. 2013). When the same procedure is performed for the WILL H₂O 1₁₀ −1₀₁ observations, removing the shock emission using the Gaussian decomposition of the profiles from San José-García (2015), six sources (3 Class 0 and 3 Class I) show inverse P-Cygni and two sources (both Class 0) show regular P-Cygni line profiles (see Fig. 13).

The two WILL sources in Serpens South have broad water absorption features, but they are not offset enough to be considered inverse P-Cygni (see Fig. 2). These may trace the large-scale cloud collision identified by Kirk et al. (2013), similar to that in the Serpens Main cloud first identified by Duarte-Cabral et al. (2011) and revealed to be the origin of the strong inverse P-Cygni line profile in water in Serpens-SMM4 by Mottram et al. (2013).

As an optically thick, higher density tracer, HCO⁺J = 4−3 is also sensitive to infall and expansion motions in protostellar envelopes (e.g. Gregersen et al. 1997; Myers et al. 2000). Though single-dish observations do not show absorption below the continuum, they can exhibit asymmetric line profiles, with the peak shifted to either the blue (infall) or red (expansion). The asymmetry between the red and blue peaks can be quantified using the δν parameterisation suggested by Mardones et al. (1997, see Sect. A.6: $\mathit{\delta }\nu \mathrm{=}\frac{{\nu }_{\mathrm{thick}}\mathrm{-}{\nu }_{\mathrm{LSR}}}{{\mathit{FWHM}}_{\mathrm{thin}}}\mathit{,}$(11)where ν_thick is the velocity of the peak emission in an optically thick tracer (in this case HCO⁺) and FWHM_thin is the line width of an optically thin tracer (in this case C¹⁸O). Values above or below 0.25 indicate a shift in the optically thick line of more than a quarter of the optically thin line width and so are considered significant. Thus, values above 0.25 indicate expansion motions while those below −0.25 indicate infall. Values between −0.25 and 0.25 are consistent with the optically thick and thin tracers being in agreement. A histogram of the values calculated for the WILL sample is shown in Fig. 14, with five sources showing blue asymmetry (i.e. δν < −0.25, 3 Class 0 and 2 Class I) and six showing red asymmetry (i.e. δν > 0.25, 3 Class 0 and 3 Class I).

Only one source (PER 08) exhibits non-static line signatures in both water and HCO⁺, however they conflict as the water shows an inverse P-Cygni line profile and the HCO⁺ a red asymmetry, leaving the status of this source uncertain. It is likely that these two tracers probe different radii, and so perhaps infall and expansion dominate in different parts of the envelope.

What is possibly more puzzling is that the vast majority of sources do not show indications of either infall or expansion in either tracer. Two-sample Kolmogorov-Smirnov (K-S) tests were performed comparing the cumulative distributions of source properties (e.g. L_bol, T_bol etc.) for sources that show infall or expansion motions with respect to those that do not, in order to see if any source properties correlate with the detection of infall or expansion signatures. Only the integrated intensity of C¹⁸O 3−2 shows a statistically significant difference (≤1% chance of being drawn from the same distribution) between sources with either an infall or expansion signature and those that do not: sources with higher C¹⁸O line intensity are more likely to show signs of infall.

Sources with clear infall motions are statistically more likely to have higher FWZI in ¹²CO 3−2, that is, broader outflow line-wings, than those showing no or expansion envelope motions (0.9% likelihood of coming from the same distribution). However, there are some sources that have high FWZI but no indication of radial envelope motions. One other result worth noting is that the presence of envelope motions is not more likely for certain values of outflow inclination (see van der Marel et al. 2013, and Sect. A.4 for details of how these were determined), suggesting that the orientation of the protostellar system is not the overriding cause of not detecting infall or expansion in our observations.

Infall must take place in all protostars, at least in the early phases, and at later times it seems unlikely that expansion of the envelope is restricted to a few select sources. Thus the low detection rate of such signatures in both tracers and the lack of a consistent trend with evolutionary Class is puzzling. It may well be that this is an observational effect, caused by small infall motions being lost in the general turbulent field on the large spatial scales that dominate single-dish observations. Mapping the velocity field inside protostellar envelopes with interferometers is likely needed to conclusively understand how material moves radially, and how this varies between different sources and over time (e.g. Yen et al. 2014; Aso et al. 2015; Evans et al. 2015).

6. Evolution of water line profile components from Class 0 to Class I

The intensity and line-width of water emission decreases for WISH sources between Class 0 and I (Kristensen et al. 2012; Mottram et al. 2014), while the rotational temperature of mid-J CO and water excitation conditions do not (Karska et al. 2013; Mottram et al. 2014). Mottram et al. (2014) therefore suggested that the observed evolution in water line intensity and line-width from Class 0 to I was caused by a decrease in the velocity of the wind driving the outflow, due to the increase of the outflow cavity opening angle as proposed in the models of Panoglou et al. (2012), for example, rather than a decrease in density in the H₂O emitting gas. However, the small sample size of the WISH survey when broken down into evolutionary classifications meant that some trends, while hinted at by the data, were not statistically significant.

Inclusion of the WILL sample helps to resolve this issue. For example, Fig. 15 shows a histogram of the components found in HIFI water spectra, updated from that shown by Kristensen et al. (2012) to include the WILL survey sources and now using the more physically motivated nomenclature introduced by Mottram et al. (2014). A broad cavity shock component is observed in almost all Class 0 and I sources. Mottram et al. (2014) argued that this is caused by C-type shocks in the outflow cavity wall, though Yvart et al. (2016) suggest an alternate explanation where this component is formed in a dusty disk wind. Spot shock components, associated with offset J-type shocks either in bullets along the jet or at the base of the outflow (Kristensen et al. 2013; Mottram et al. 2014), are far more likely to be detected in Class 0 than Class I sources. Inverse P-Cygni line profiles associated with infall are more common for younger sources, though the inclusion of the WILL sample means that expansion motions traced by regular P-Cygni line profiles are now approximately equally common in both Class 0 and Class I sources. With the exception of the regular P-Cygni profiles, the evolution of water line profile components found for the combined WISH and WILL samples confirms the conclusion of Kristensen et al. (2012) and Mottram et al. (2014) that the outflows of young Class 0 sources are more energetic and their envelopes more infall-dominated than their more evolved Class I counterparts.

Two-sample K-S tests show that there is less than 2% chance of the Class 0 and I distributions of M_env, M_out, F_CO, Ṁ_out, the integrated intensity of HCO⁺ and H₂O 1₁₀ −1₀₁, and the FWZI of CO J = 3−2 and H₂O 1₁₀ −1₀₁ being drawn from the same distribution, with the values for Class 0 sources being larger on average than those of Class I sources. In particular, this confirms the decrease in F_CO and both the line-width and intensity of water emission with evolution of the central source, reinforcing the direct relation between water emission and outflow/shock activity.

One caveat is that the Class 0 and I sources in WISH and WILL sources are not evenly drawn from the sampled star-forming regions. For example, only Class I or II sources are included in the Taurus star-formation region, while many of the Class 0 sources are in Perseus, which is a much more active star-forming complex. The observed differences and trends between Class 0 and I may therefore be accentuated by environmental differences. However, all the evidence suggests that Class I sources have slower, less powerful outflows and show less sign of strong infall motions in their envelopes than Class 0 sources.

7. Summary and conclusions

This paper has presented a set of Herschel and ground-based follow-up observations, characterisation and initial analysis of a flux-limited sample of Class 0/I YSOs in the Gould Belt. From this comprehensive dataset, combined with observations from the WISH and DIGIT surveys, we are able to conclude that:

• Water line profiles are dominated by emission from the actively shocked regions in outflows, the activity of which decreases in strength (i.e. has lower intensity and FWZI, see Fig. 8) and has fewer J-type shocks (i.e. fewer spot-shock components, see Fig. 15) as sources evolve from Class 0 to I. We also confirm the decrease in the force of the cooler and slower entrained outflowing gas, measured from low-J CO, as young embedded protostars evolve.

• The ratio of mass in the entrained outflow to envelope mass (i.e. M_out/M_env) remains relatively constant between Class 0 and I with a median of approximately 1%, consistent with a core-to-star formation efficiency of approximately 50% and an outflow duty cycle of approximately 5%.

• F_CO/Ṁ_accis relatively constant with L_bol, T_bol and M_env, suggesting that the entrainment efficiency is constant and independent of the power and evolution of the driving source of the flow. The constant value of F_CO/Ṁ_acc implies a median velocity at the wind launching radius of 6.3 km s^-1. This in turn suggests an entrainment efficiency of approximately 30−60% if the wind is launched around 1AU, or close to 100% if it is launched at larger radii.

• L[O i] is strongly correlated with L_bol but not with M_env, in contrast to low-J CO, which is strongly correlated with M_env and more weakly related to L_bol. This suggests that [O i] is more closely related to the current accretion activity while low-J CO traces the average activity over timescales of approximately 10² −10⁴ yr. H₂O is more strongly correlated with L_bol than M_env, but with a smaller difference than for low-J CO, consistent with it tracing actively shocked material between the wind and entrained outflow.

• L[O i] does not vary significantly between Class 0 and I, likely because the molecular to atomic ratio in the wind and jet decreases as the source evolves, as suggested by Nisini et al. (2015). This could be caused by increased temperature and decreased density (and thus shielding) in more evolved sources. [O i] is therefore a poor tracer of the time-averaged mass-loss rate, and thus a poor alternative to CO.

• Infall signatures are predominantly seen in Class 0 sources in both H₂O and HCO⁺ single-dish observations, but with little overlap in detections between the two tracers. However, infall signatures remain elusive in the majority of sources. Thus, while water is a good tracer of infall, it is by no means the universal tracer needed to understand how this proceeds in general.

• The conclusions drawn from the WISH sample hold and become more statistically robust when the combined WISH+WILL+DIGIT sample is analysed.

Further exploitation of these data and this sample can be found for the HIFI data in San José-García (2015) and will be presented for the PACS data in Karska et al. (in prep.). The use of PACS and HIFI in this way is cementing the unique legacy of Herschel on the energetics of star formation and the origin of water in the interstellar medium.

---

Acknowledgments

We thank the anonymous referee for constructive comments that improved the clarity and content of the paper. J.C.M. acknowledges support from grant 614.001.008 from The Netherlands Organisation for Scientific Research (NWO), from the European Union A-ERC grant 291141 CHEMPLAN, and from the European Research Council under the European Community’s Horizon 2020 framework program (2014–2020) via the ERC Consolidator grant “From Cloud to Star Formation (CSF)” (project number 648505). Astrochemistry in Leiden is supported by The Netherlands Research School for Astronomy (NOVA), by a Spinoza grant, by a Royal Netherlands Academy of Arts and Sciences (KNAW) professor prize, and by the CHEMPLAN A-ERC grant. A.K. acknowledges support from the Foundation for Polish Science (FNP) and the Polish National Science Center grant 2013/11/N/ST9/00400. D.F. acknowledges support from the Italian Ministry of Education, Universities and Research project SIR (RBSI14ZRHR). 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; The 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. PACS has been developed by a consortium of institutes led by MPE (Germany) and including UVIE (Austria); KU Leuven, CSL, IMEC (Belgium); CEA, LAM (France); MPIA (Germany); INAF-IFSI/OAA/OAP/OAT, LENS, SISSA (Italy); IAC (Spain). This development has been supported by the funding agencies BMVIT (Austria), ESA-PRODEX (Belgium), CEA/CNES (France), DLR (Germany), ASI/INAF (Italy), and CICYT/MCYT (Spain). This research made use of astropy, a community-developed core python package for Astronomy (Astropy Collaboration et al. 2013), and aplpy, an open-source package for plotting astronomical images with python hosted at http://aplpy.github.com.

References

Appendix A: Property determination

This section presents the various properties of the WILL sources, and details of how they were determined.

Appendix A.1: Spectral energy distributions

In constructing the spectral energy distributions (SEDs) of the WILL sources, the photometric flux densities from the near-IR to 24  μm from 2MASS (Skrutskie et al. 2006) and Spitzer have been collected, where detected, from the latest determinations by Dunham et al. (2015). Flux densities for detections at 450 and/or 850  μm with SCUBA on the JCMT were taken from the catalogue of Di Francesco et al. (2008), and the 1.3 mm detections with MAMBO on the IRAM 30m telescope by Maury et al. (2011) for Aquila, Serpens South and W40 were also included.

All WILL sources lie within the Herschel PACS and SPIRE photometric maps observed in parallel-mode at 70, 160, 250, 350 and 500  μm as part of the Herschel Gould Belt Survey, with many also within the smaller regions observed in PACS-only mode at 100 and 160  μm. Where both parallel and PACS-only mode observations are available at 160  μm we use the data taken in PACS-only mode, as these were observed with a slower scanning speed and thus are of better quality than those data taken in parallel with SPIRE observations. Processed and calibrated mosaics were downloaded from the Herschel Science Archive (HSA⁷, see also the Herschel Gould Belt survey archive for further products⁸). The observation identification numbers for these maps are given in Table B.2.

Aperture photometry was then performed using python routines from the astropy package to extract flux densities for all sources at all available wavelengths. The starting value for the source aperture, inner and outer sky radii for each source at each wavelength was set at 6, 10 and 12 pixels, typically corresponding to 6, 10 and 12 times the beam size. Rings at these radii were then overlaid on images of the data and the radii adjusted to best encompass each source and local background respectively, while excluding nearby sources. Aperture correction factors were calculated and applied for each combination of aperture and sky annuli using the relevant PACS and SPIRE point source function images.

In some cases, additional nearby sources are blended with the primary WILL source at longer wavelengths where the Herschel beam becomes large. For these sources, we specify the longest wavelength where the sources are reliably separated. The flux densities at longer wavelengths are then scaled down by the ratio of the flux densities in the reliable image between the aperture used for that wavelength and the aperture used for the longer wavelength data. For example, if a source is blended at 500  μm, then the flux density at that wavelength is scaled down by the ratio of the flux density at 350  μm to the flux density in the 350  μm image within the region defined by the aperture used for the 500  μm data.

Finally, continuum flux densities from the WILL PACS spectra were also obtained in the same manner as for the WISH sources (see Karska et al. 2013, for details). However, it is sometimes difficult in the PACS spectral maps to separate the emission related to the protostar from surrounding emission. Therefore, in cases where the PACS spectral continuum is significantly higher than obtained from the broad-band photometric maps, the latter is preferred and the PACS spectral continuum flux densities are not included in the SED. The broad-band photometric flux densities for all sources are given in Table A.1, with the PACS spectral continuum flux densities given separately in Table A.2. The SEDs of all sources are shown in Fig. A.1.

For Aquila sources 01−04, the peak PACS continuum is offset from the coordinates in Maury et al. (2011), and so extraction of the SED and ground-based molecular line emission is performed at the peak flux position.

Appendix A.2: Continuum properties

A number of properties are often calculated from YSO SEDs in order to characterise them. While previous determinations have been made for the better-studied members of the WILL sample, variations in method and implementation could lead to biases or systematic effects between sources in different regions. In addition, Herschel data provide significantly higher-resolution images for flux density determination than were previously available between 70 and 250  μm (e.g. Herschel has a FWHM of ~8′′ at 100  μm compared to 3′ × 5′ for IRAS). In order to be consistent with other Herschel surveys, most notably the WISH and DIGIT surveys, we follow the definitions in Dunham et al. (2010) and calculate the bolometric luminosity (L_bol) using: ${\mathit{L}}_{\mathrm{bol}}\mathrm{=}\mathrm{4}\mathit{\pi }{\mathrm{d}}^{\mathrm{2}}{\mathrm{\int }}_{\mathrm{0}}^{\mathrm{\infty }}{\mathit{F}}_{\mathit{\nu }}\mathrm{d}\mathit{\nu ,}$(A.1)the sub-mm luminosity using: ${\mathit{L}}_{\mathrm{smm}}\mathrm{=}\mathrm{4}\mathit{\pi }{\mathrm{d}}^{\mathrm{2}}{\mathrm{\int }}_{\mathrm{0}}^{\mathit{\nu }\mathrm{=}\mathit{c}\mathit{/}\mathrm{350}\mathit{\mu }\mathrm{m}}{\mathit{F}}_{\mathit{\nu }}\mathrm{d}\mathit{\nu ,}$(A.2)and bolometric temperature using: ${\mathit{T}}_{\mathrm{bol}}\mathrm{=}\mathrm{1.25}\mathrm{\times }{\mathrm{10}}^{-11}\frac{{\mathrm{\int }}_{\mathrm{0}}^{\mathrm{\infty }}\mathit{\nu }{\mathit{F}}_{\mathit{\nu }}\mathrm{d}\mathit{\nu }}{{\mathrm{\int }}_{\mathrm{0}}^{\mathrm{\infty }}{\mathit{F}}_{\mathit{\nu }}\mathrm{d}\mathit{\nu }}\mathrm{·}$(A.3)In these equations, F_ν is the flux density at frequency ν, and d is the distance to the source. The integrals were calculated over the available SED flux densities using trapezium integration, which was found by Karska et al. (2013) to provide the most consistent results for the WISH survey sources.

The infrared spectral index (see Lada 1987) is also calculated for those sources where the source is detected in at least three wavelengths between 2 and 25  μm using a least-squares linear fit to the logarithms of wavelength and flux density such that: ${\mathit{\alpha }}_{\mathrm{IR}}\mathrm{=}\frac{\mathrm{d}{\log }_{\mathrm{10}}\mathrm{\left(}\mathit{\lambda }{\mathit{F}}_{\mathit{\lambda }}\mathrm{\right)}}{\mathrm{d}{\log }_{\mathrm{10}}\mathrm{\left(}\mathit{\lambda }\mathrm{\right)}}\mathit{,}$(A.4)where F_λ is the flux density at wavelength λ.

Finally, for the sources in the Aquila Rift we take the envelope mass (M_env) from Maury et al. (2011), calculated from the integrated intensity at 1.2 mm, corrected to the updated distance for that region. In other regions we calculate M_env for those sources with SCUBA 850  μm observations using Eq. (1) from Jørgensen et al. (2009): ${\mathit{M}}_{\mathrm{env}}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\mathrm{0.44}\hspace{0.17em}{\left(\frac{{\mathit{L}}_{\mathrm{bol}}}{\mathrm{1}\hspace{0.17em}{\mathit{L}}_{\mathrm{\odot }}}\right)}^{-0.36}{\left(\frac{{\mathit{S}}_{{\mathrm{15}}^{\mathrm{\prime \prime }}\mathit{,}\hspace{0.17em}\mathrm{850}\mathit{\mu }\mathrm{m}}}{\mathrm{1}\hspace{0.17em}{\mathrm{Jy}\hspace{0.17em}\mathrm{beam}}^{-1}}\right)}^{\mathrm{1.2}}{\left(\frac{\mathrm{d}}{\mathrm{125}\hspace{0.17em}\mathrm{pc}}\right)}^{\mathrm{1.2}}{\mathit{M}}_{\mathrm{\odot }}\mathit{,}$(A.5)where S_{15′′, 850 μm} is the peak SCUBA 850  μm flux density in a 15′′ pixel. This empirical relation was derived from comparison of observed source properties with the results of dust radiative transfer models. The results of these calculations are summarised in Table 1.

Appendix A.3: PACS line properties and detection statistics

Table 3 gives the wavelength ranges covered by the two PACS settings used for the WILL survey. Table A.3 summarises the properties of the principle transitions observed with PACS towards the WILL sample, while Table A.4 indicates which of these are detected towards each WILL source. The overall fraction of sources detected in each transition is also given in the bottom row of Table A.4. The spectral resolution of PACS is not sufficient to separate the CO J = 31−30 transition at 84.41  μm and the OH 84.42  μm transition. Emission at this wavelength is therefore marked as a detection for both lines but could be in only CO or OH. The detection (or not) of neighbouring transitions is likely a reasonable estimate of whether the detected emission is from CO, OH or a blend of the two.

Appendix A.4: Entrained outflow

Molecular outflows, usually detected through observations of low-J¹²CO, are a ubiquitous signpost of ongoing star formation. In order to identify and quantify the properties of the entrained outflowing material associated with WILL sources from ¹²CO J = 3−2 maps, a number of steps were taken, following van der Marel et al. (2013):

• 1.

  The data were resampled to 0.5 km s^-1 to improve the sensitivity to line-wing emission.

• 2.

  Maps of maximum red and blue-shifted velocities were identified in all spectra as the channel where the emission first reaches 1σ_rms.

• 3.

  The outer velocity (ν_out) of each outflow lobe was defined from the maximum velocity maps as the most offset value from the source velocity. The maximum velocity (ν_max) is the absolute difference between ν_out and ν_LSR.

• 4.

  The inner velocity (ν_in) for each lobe was defined by the same approach using a spectrum without any outflow emission, so as to mask out cloud emission.

• 5.

  Integrated maps for the red and blue outflow lobes were created by integrating between the minimum and maximum velocities. Visual inspection and comparison of these maps with the continuum emission was then used to identify the spectra associated with an outflow from the source. Excluding low-velocity emission will lead to an underestimate in the mass and related properties by factors of a few (Downes & Cabrit 2007). However, this is generally preferable to performing an incorrect correction based on poor knowledge of the contribution from the envelope (Cabrit & Bertout 1990).

• 6.

  Each source was assigned an inclination (i) of 10, 30, 50 or 70° between the outflow axis and the line of sight, such that i = 0° is pole-on, through visual inspection of the overlap between the red and blue lobes with each other and the source position.

• 7.

  The radius associated with each outflow lobe (R_CO) was defined as the distance between the source position and furthest pixel containing outflow emission.

• 8.

  The mass in each channel of each pixel (m_ij) was calculated assuming an excitation temperature (T_ex) of 75 K, μ = 2.8 (Kauffmann et al. 2008) and a CO abundance relative to H₂ of 1.2 × 10^-4 (Frerking et al. 1982), consistent with the outflow properties determined by Yıldız et al. (2015) for the WISH sample. Changing T_ex to 100  or 50 K would only raise or lower the mass by a factor of 1.2. No correction for τ is performed, as the optically thick parts of the line near the line centre are excluded and the line wings are typically optically thin in outflows associated with low-mass protostars (see van der Marel et al. 2013; Yıldız et al. 2015).

• 9.

  The physical properties of the outflow were calculated using the separation method (M7 in van der Marel et al. 2013), where the mass, momentum and energy are calculated separately for each lobe on a per-channel basis for each spectrum, then summed over all channels and spectra, while the maximum velocity is used to calculate the dynamical time of the flow, that is: $\begin{array}{ccc}& & {\mathit{M}}_{\mathrm{out}}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathrm{npix}}\sum _{{\nu }_{\mathrm{in}}}^{{\nu }_{\mathrm{out}}}{\mathit{m}}_{\mathit{i,j}}\hspace{0.17em}\mathrm{d}\mathit{v,}\\ & & {\mathit{P}}_{\mathrm{out}}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathrm{npix}}\sum _{{\nu }_{\mathrm{in}}}^{{\nu }_{\mathrm{out}}}{\mathit{m}}_{\mathit{i,j}}\mathrm{|}{\nu }_{\mathit{i}}\mathrm{-}{\nu }_{\mathrm{LSR}}\mathrm{|}\hspace{0.17em}\mathrm{d}\mathit{v,}\\ & & {\mathit{E}}_{\mathrm{out}}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathrm{npix}}\sum _{{\nu }_{\mathrm{in}}}^{{\nu }_{\mathrm{out}}}{\mathit{m}}_{\mathit{i,j}}\mathrm{|}{\nu }_{\mathit{i}}\mathrm{-}{\nu }_{\mathrm{LSR}}{\mathrm{|}}^{\mathrm{2}}\hspace{0.17em}\mathrm{d}\mathit{v,}\\ & & {\mathit{t}}_{\mathrm{dyn}}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\frac{{\mathit{R}}_{\mathrm{CO}}}{{\nu }_{\max }}\mathit{,}\\ & & \mathit{Ṁ}\mathrm{out}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\frac{{\mathit{c}}_{\mathrm{f}}\hspace{0.17em}{\mathit{M}}_{\mathrm{out}}}{{\mathit{t}}_{\mathrm{dyn}}}\mathit{,}\end{array}$and ${\mathit{F}}_{\mathrm{CO}}\hspace{0.17em}\mathrm{=}\hspace{0.17em}\frac{{\mathit{c}}_{\mathrm{f}}\hspace{0.17em}{\mathit{P}}_{\mathrm{out}}}{{\mathit{t}}_{\mathrm{dyn}}}\mathit{,}$(A.11)where c_f is a correction factor to account for inclination, given in Table A.5 and derived from the models of Downes & Cabrit (2007) by van der Marel et al. (2013).

Table A.6 gives the calculated outflow properties for the red and blue lobes separately, including the velocity limits, outflow mass (M_out), momentum (P_out), kinetic energy (E_out), radius, inclination, dynamical time (t_dyn), force (F_CO) and mass outflow rate (Ṁ_out). In some cases, the outflow may extend beyond the 2′ × 2′ coverage of the observations, so these values may be lower limits to the total value. Figure A.2 shows the outflow lobes of the observed outflows associated with the WILL sources overlaid on the Herschel PACS 70 micron maps.

Appendix A.5: Mass accretion and loss

Table A.7 presents the calculated mass accretion rates for all WILL, WISH and DIGIT sources, calculated as discussed in Sect. 4.4 using Eq. (6), as well as the observed luminosity in the [O i] 63  μm transition and the mass-loss rate in the wind derived using the relation from Hollenbach (1985), given in Eq. (7). L [ OI 63μm ] is not given for OPH 01, W40 01 and W40 03−06 as the detections in [O i] towards these sources are almost certainly due to PDR emission.

Appendix A.6: Ground-based line fitting

Table A.8 presents the basic properties of all transitions observed with the JCMT and APEX. Gaussian fitting was performed for the central spectra of all sources for ¹³CO, C¹⁸O, HCO⁺ and H¹³CO⁺ in order to determine line-widths and central velocities. The ν_LSR for each source is defined as the peak position of the fit to the C¹⁸O J = 3−2 central spectrum. For the four sources where this is not detected, three (SCO 01, TAU 03, TAU 07) were also observed and detected in ¹³CO J = 3−2, so the fit to this line is used instead. For the remaining source (TAU 08) we use the velocity derived by Caselli et al. (2002) from N₂H⁺J = 1−0 observations. The integrated intensity of the lines is measured over a window of ±3 FWHM.

For the optically thick HCO⁺ and ¹³CO, we also quantify any blue/red asymmetry by calculating δν as defined by Mardones et al. (1997) using Eq. (11). For most sources, we use the peak of the Gaussian fit to the line, but in the case of five sources (AQU 02, PER 04, SERS 01, SERS 02 and TAU 09), the position of the maximum intensity is used in this calculation for HCO⁺ because strong self-absorption or a broad, fainter second component skews the Gaussian fit away from the maximum intensity.

The results for all sources are presented in Table A.9, and Fig. A.3 shows the central HCO⁺ and C¹⁸O spectra for all sources.

Appendix A.7: Evolutionary classification

Table A.10 presents a summary of the various indicators used to reach the final classification of the evolutionary state of sources in the WILL sample. Firstly, we consider whether or not an entrained molecular outflow is associated with the source in CO J = 3−2, and whether or not there are broad line-wings in the HIFI water and CO J = 10−9 spectra. Sources with all three signatures are likely to be the youngest protostars, with strong, energetic and likely warm outflows. Those without any detected outflow signatures are likely either pre-stellar or more evolved (i.e. Class II) sources.

Next, we follow the method of van Kempen et al. (2009) and Carney et al. (2016) using maps of the molecular emission in HCO⁺J = 4−3 and C¹⁸O J = 3−2 to separate Class 0/I embedded protostars from edge-on Class II disk sources. If both transitions are strong and spatially concentrated then the source is most likely a genuine embedded (i.e. Class 0/I) YSO. If not, the source is either too cold (i.e. pre-stellar) or does not have a significant envelope and so isa more evolved disk source. In the W40 sources, the extended emission in both lines is likely due to the PDR and may mask the presence of compact emission associated with the sources, so they are designated as confused. A few other sources also receive this designation, in line with Carney et al. (2016), if there are multiple sources in the field that cannot be disentangled at the resolution of the observations.

The presence of strong, compact sub-mm continuum emission is also indicative of a young (i.e. pre-stellar or embedded Class 0/I) source, allowing pre-stellar and Class II sources to be distinguished (e.g. André et al. 2000). Finally, the T_bol classification (i.e. T_bol < 70 K corresponds to Class 0, 70 ≤ T_bol < 650 K to Class I) can be used to separate Class 0 and I sources, though on its own this calculation is not always able to correctly separate other highly-reddened sources; for example edge-on Class II disks (Crapsi et al. 2008).

Overall, the WILL sample consists of 23 Class 0, 14 Class I, 8 Class II and 4 potentially pre-stellar sources. Most of the Class II sources are in Taurus, while the pre-stellar sources are all in Aquila/W40. Six sources (one Class 0, one Class II and four pre-stellar) have narrow yet bright ¹²CO J = 10−9 emission (see Fig. 3) suggestive of a PDR along the line of sight to the source. The details of some of the more complex determinations are discussed in Appendix C.

Appendix B: Observation IDs

This section presents the Herschel observation ID numbers for all WILL HIFI and PACS spectral observations (Table B.1), as well as those of the Herschel PACS and SPIRE photometric maps used to extract source photometry (Table B.2).

Appendix C: Discussion of individual cases

This section presents notes on individual sources to explain oddities in the data that bear specific mention.

AQU 01: the additional emission component observed only in the H₂O 1₁₀ −1₀₁ line towards this source (which has the largest beam) is almost certainly due to another source on the edge of the beam. This source is outside the beam for all other HIFI observations.

CHA 01: this source shows a relatively small and weak but detectible outflow in CO J = 3−2, as well as a relatively narrow detection in H₂O with HIFI and detections in [O i], H₂O, CO and OH with PACS. The principle reason this source is not designated as a Class I is the non-detection of HCO⁺. The low velocity of the outflow suggests this may either be a remnant from the Class I phase, or a disk-wind such as recently seen by ALMA in HD 163296 (Klaassen et al. 2013).

CHA 02: the non-detection of water emission in the HIFI observations for this Class I source is unsurprising given that only a few lines are marginally detected in PACS and no outflow is detected in CO J = 3−2. There is a detection in [O i], suggesting that some form of jet or wind is present. This source is therefore probably nearing the Class II phase and simply has too tenuous an envelope for a significant outflow component to still be present and detectible.

OPH 01: this source does not show compact HCO⁺ or C¹⁸O emission, so is not Class 0/I. There is a cold, starless core to the north that increasingly dominates at longer wavelengths (e.g. Sadavoy et al. 2010) and is causing the source to appear more embedded than it really is. The strong, narrow CO J = 10−9 emission profile, along with the PACS detection of CO J = 16−15 are most likely from a PDR or bow-shock, possibly caused by interaction between the infrared source and the cold core if they are spatially associated. Given that the primary infrared source is visible in the mid-IR (e.g. Brown et al. 2013), we classify it as Class II, whilst noting that there is also PDR-like emission.

PER 02: high-resolution BIMA observations of the CO outflows and Spitzer IRAC scattered light observations show that the blue outflow lobe of this source is contaminated by the blue outflow lobe originating from L1448-MM to the south (Kwon et al. 2006; Tobin et al. 2007). We therefore do not report outflow properties as it is impossible to disentangle the two flows. However, the BIMA observations do show activity related

to the targeted source, so it is still classified as having a detected outflow for the purposes of source evolution. The broad absorption in the blue outflow lobe in the ground-state water lines likely takes place against the outflow from L1448-MM, indicating that it is between this source and the observer. In particular, the saturated absorption feature below the continuum in H₂O 1₁₁ −0₀₀ (see San José-García 2015) rules out the possibility that this is instead caused by a combination of emission from the different outflows in the beam. Thus, the source is most likely further away than the outflow from L1448-MM and the two are probably not interacting directly.

PER 04: the non-detection of water for this Class 0 source is slightly surprising given the detection of an outflow in CO J = 3−2. However, the low velocity of the outflow may mean that any shocks are not fast enough to sputter water from the grains or warm enough to lead to efficient gas-phase formation. This may indicate that this source is particularly young or has been in a lower-accretion phase for some time.

PER 06: while the CO J = 3−2 data are consistent with an outflow originating from this source (NGC 1333-IRAS2B), it is impossible to disentangle the contribution of NGC 1333-IRAS2A, something even Plunkett et al. (2013) find impossible at three times higher spatial resolution with CARMA, so we follow those authors in not quoting any outflow property values.

PER 08: at least part of the strange outflow morphology from this source may be from another nearby source or indicate that it is a multiple. However, it is a single star in high-resolution VLA observations (Tobin et al. 2016) and there are no other obvious infrared candidates nearby so for now we attribute all the outflow emission to this protostar.

PER 12: the red outflow lobe of NGC 1333-IRAS4A passes through the H₂O 1₁₀ −1₀₁ beam but not the other transitions that have smaller beam-sizes; hence the detection in only this transition. There is no evidence in the CO J = 3−2 data of an outflow associated with this source. However, Tobin et al. (2015) note that Spitzer images are suggestive of a jet or outflow related to the source (see their Fig. 19), so the non-detection is likely because the outflow is in the plane of the sky rather than because this source does not have an outflow.

PER 22: the morphology of the outflow in CO J = 3−2 is suggestive of two outflows, particularly in the red lobe, with one being approximately north-south and the other east-west. Enoch et al. (2009) find another, more evolved Class I source (Per-emb 55) ~9′′ away, which therefore lies within the Herschel

beam. VLA observations resolve this additional source into a binary (Tobin et al. 2016), suggesting that this is actually a triple if Per-emb 8 and 55 are spatially associated.

SERS 01: this source is significantly offset in both the continuum and molecular lines from the position given by Maury et al. (2011), which is likely why it is relatively weak in the HIFI and PACS observations.

TAU 03: most of the outflow near this source is likely due to a neighbouring, probably younger, source, though there does seem to be a weak flow from the primary target as well.

TAU 06: the outflow observed in CO J = 3−2 seems to have two blue outflow lobes, and the red and blue outflow lobes are not well aligned. Given the low velocity of the outflow, it may be close to the plane of the sky, or there may be multiple outflows in the region. There is no other obvious driving source in the vicinity however, so we assign all the emission to the target.

TAU 07−09: these Class II sources exhibit some line broadening in CO J = 3−2 around the source position. While the line-wings are not at high enough velocity offset (FWZI ≲ 10 km s^-1) to be considered related to a true outflow, this is suggestive of either higher turbulence or a small disk wind (cf. Klaassen et al. 2013) in these sources.

W40 01: the HIFI H₂O spectra show a broad line-wing, in addition to PDR-like absorption close to the ν_LSR, consistent with an outflow, though there is no outflow detected in low-J CO.The source is compact and relatively bright at 70  μm, which led Könyves et al. (2015) to designate this source as protostellar, rather than prestellar, and it is detected in the mid-IR. It is therefore most likely a Class 0 source, though the presence of the W40 PDR complicates the detection and passage of the outflow.

W40 03−06: these sources show little or only very weak continuum emission at wavelengths shorter than 100  μm, while at longer wavelengths, the emission is not particularly compact. They also do not show any signs of outflow and the only line detections are related to emission from the W40 PDR. We therefore designate these sources as potentially prestellar, but note that the PDR makes a reliable classification more difficult.

W40 07: the outflow observed towards this source in CO J = 3−2 is surprisingly strong given that there is no emission detected in H₂O and only a faint, narrow and tentative detection in CO J = 10−9. The low fluxes in the mid-IR would seem to suggest that there is little warm gas around this source, but the sub-mm and mm continuum detections indicate a significant reservoir of cold dust. The bright and compact nature of the 70  μm PACS continuum emission and shape of the SED led Könyves et al. (2015) to designate this source as protostellar. It is therefore most likely a very young protostellar source where the outflow has not become fast and warm enough to be detected in H₂O and high-J CO.

All Tables

All Figures