© ESO, 2015

1. Introduction

It is now well established that the formation of protostars powers fast jets surrounded by wide-angle winds that impact the high-density parent cloud. These jets interact through shock fronts that compress and heat the ambient gas while driving the formation of low-velocity bipolar outflows. The low-excitation lines of CO have been the privileged tool of the study of protostellar outflows for decades (see e.g. Bachiller 1996). Many studies on CO outflow emission conducted from the ground have reported the presence of warm gas with a kinetic temperature in the range 50–200 K (see e.g. Hatchell et al. 1999a; Stark et al. 2004; van Kempen et al. 2006, 2009a,b; Gomez-Ruiz et al. 2013).

However, the atmospheric transmission at submillimeter and far-infrared wavelengths severely limits the CO transitions that can be observed from the ground above J = 7–6. As a consequence, it is unclear whether any higher temperature component is present in the outflowing gas. In addition, the origin of the warm gas component is not well constrained and several processes have been proposed, such as outflow shocks, and irradiation processes associated with shocks or from the protostar (see e.g. van Kempen et al. 2009b). As suggested by Visser et al. (2012), the gas heating mechanism may depend on the evolutionary stage of the protostar itself.

This situation has changed with the advent of the Herschel Space Observatory¹. For the first time, it has been become possible to observe in a systematic way the CO ladder up to very high energy levels (J_up = 16) at high spectral resolution. This allows a much more accurate determination of the physical conditions and the origin of the high-temperature gas in protostellar outflows. Recently, Yildiz et al. (2013) have carried out a systematic study of the high-J lines of ¹²CO and its isotopologues up to J_up = 10 (E_up up to 300 K) towards a sample of 26 deeply embedded low-mass Class 0 and Class I young stellar objects as part of the Water In Star-forming regions with Herschel (WISH) key program (van Dishoeck et al. 2011). The aim was to constrain the physical characteristics (excitation, kinematics, column density) of the warm gas in low-mass protostellar envelopes and compare them with the colder gas traced by lower excitation lines. The authors have reported CO excitation temperatures in the range 60–130 K. The outflow of the star-forming region L1157 is one of the very few outflows whose emission was investigated at high spectral resolution in the high-excitation lines of CO, up to J_up = 16 (E_up = 751.8 K), as part of the Chemical HErschel Spectral Surveys (CHESS) key program (Ceccarelli et al. 2010). Lefloch et al. (2012) and Busquet et al. (2014) showed evidence for temperature stratification inside the outflow, with the presence of temperature components above 210 K, which are best studied in the high-excitation transitions of CO J_up ≥ 10.

The heterodyne receiver GREAT² on board SOFIA is the only heterodyne instrument that permits observation of the highly excited CO transitions J_up = 11–22, at an angular resolution of 10″–20″, hence comparable to the ground-based submillimeter telescopes.

The intermediate-mass Class 0 source Cep E-mm (L = 100 L_⊙; Lefloch et al. 1996, Chini et al. 2001) in the Cepheus OB3 association at a distance of 730pc (Sargent 1977) drives an exceptionally powerful and luminous molecular outflow, which has long been studied in the rovibrational and pure rotational lines of H (Eisloeffel et al. 1996; Noriega-Crespo et al. 1998, 2004, 2014). Both the protostar and the outflow have been studied in the millimeter to far-infrared wavelength range since their discovery (Eislöffel et al. 1996; Lefloch et al. 1996; Noriega-Crespo et al. 1998, 2004b; Hatchell et al. 1999b; Moro-Martin et al. 2001; Froebrich et al. 2003; Smith et al. 2003). Lefloch et al. (1996) reported the presence of high-velocity bullets of molecular gas in the outflow. The southern lobe is terminated by the Herbig-Haro object HH377, detected at optical wavelengths (Ayala et al. 2000).

Giannini et al. (2001) and Moro-Martin et al. (2001) analyzed the (unresolved) far-infrared emission lines of CO and H₂O in Cep E, as observed with ISO. They found that the high-J CO emission (J_up = 14–25) can be accounted for by an extended gas component with density and temperature in the range 4 × 10⁴–4 × 10⁶ cm^-3) and 200–1000 K, respectively. Lefloch et al. (2011) reported the presence of water-rich shocked clumps of gas in the protostellar jet of Cep E, in which they estimated typical densities n(H₂) ≃ 10⁶ cm^-3 and temperatures T_kin ~ 200 K from the detection of the masing water line 3_1,3–2_2,0 emission at 183.3 GHz. Recent observations of HH377 at high spectral resolution with GREAT/SOFIA by Gómez-Ruiz et al. (2012) showed that the emission of the CO J = 12–11 and 13–12 lines in those bullets was consistent with gas at a temperatures higher than 100 K and density higher than 10⁴ cm^-3. Interestingly, they discovered an additional high-velocity bullet in the high-J line emission. However, the compact size of the outflow (~30″ in each lobe, see Figs. 1–2) with respect to the angular resolution of the observations (20″) made it difficult to pinpoint the origin of the CO emission.

Fig. 1 Cep E outflow emission observed in the CO J = 2–1 line with the PdBI at 1″ resolution (contours) and with IRAC/Spitzer H24.5 μm line (color). Each panel shows the emission associated with each of the three outflow components identified: (left) outflow cavity walls emission, integrated between −15.2 and −8.7 km s-1 in the southern lobe (blue contours) and between −2.2 and −8.7 km s-1 in the northern lobe (white contours); (center) jet emission integrated between −135 and −110 km s-1 and between +40 and +80 km s-1 towards the southern (blue contours) and northern lobes (white contours), respectively; (right) HH377 bowshock emission integrated between −77 and −64 km s-1 (blue contours). First contour and contour interval are 20% and 10% of the peak intensity, respectively.

Fig. 2 Left: CO J = 2–1 velocity-integrated emission of the protostellar jet observed with the PdBI. The emission is integrated between −135 and −110 km s-1 and between +40 and +80 km s-1 towards the southern and northern lobes, respectively. Right: the three kinematic components in the southern lobe of the Cepheus E outflow as seen by the PdBI at 1″ resolution in the CO J = 2–1 line: a) outflow cavity emission between −15 and −12 km s-1 (color); b) HH377 bowshock emission between −77 and −64 km s-1 (yellow contours); c) jet emission between −135 and −110 km s-1 (blue contours). First contour and contour interval are 20% and 10% of the peak intensity, respectively. The white stars mark the location of the positions observed with SOFIA (−8.5″, −23″) and Herschel (HSO), IRAM, and JCMT (−12″, −18″). The beam width at half-power of SOFIA in the CO J = 16–15 and of HIFI in the J = 5–4 are drawn by a solid and a dotted circle, respectively.

In this article, we report on a multi-transition study of the CO rotational ladder from millimeter to far-infrared wavelengths of the southern lobe of the Cep E-mm outflow. The data were collected with the IRAM 30 m telescope (J = 1–0) and Plateau de Bure interferometer (PdBI; J = 2–1), the James Clerk Maxwell Telescope (JCMT; J = 3–2, J = 4–3), Herschel (J = 5–4, 9–8, 14–13), and SOFIA (J = 12–11, 13–12, 16–15). This is the most complete coverage of the CO ladder ever carried out at high spectral resolution towards a protostellar outflow. We used a map of the J = 2–1 line obtained with the PdBI at 1″ resolution to characterize the origin of the features detected with single-dish observations. The observations are described in Sect. 2. We describe in Sect. 3 the various components of the outflow and their specific spectral signatures, detected both with the PdBI and the single-dish telescopes: the protostellar jet, the bowshock HH377, and the outflow cavity. We present in Sect. 4 an analysis of the physical conditions in each outflow component. The properties of the jet and the origin of the high-excitation line emission are discussed in Sect. 5. Our conclusions are summarized in Sect. 6.

2. Observations

We describe below the CO line observations carried out with each of the following observatories: SOFIA, IRAM 30 m, IRAM PdBI, JCMT, and Herschel. Line intensities are expressed in units of antenna temperature corrected for atmospheric attenuation (for ground-based observations) . The observation and telescope parameters are summarized in Table 1. The calibration uncertainty is within 10% in the millimeter bands (J_up = 1, 2) and up to 20% in the submillimeter bands (J_up = 3, 4, 5, 9, 12, 13, 14, 16).

2.1. SOFIA

The observations in the Cep E star-forming region were conducted with the GREAT spectrometer (Heyminck et al. 2012). The position of the HH377 object, coinciding with the terminal bowshock of the outflow was targeted, with coordinates , (see Fig. 2). The receiver was tuned to the frequency of the CO lines J = 13–12 at 1496.923 GHz and J = 16–15 at 1841.346 GHz, respectively, during the Cycle 1 SOFIA flight on 2 November 2013. The emission of the CO line J = 12–11 at 1381.995 GHz was observed towards the same position in July 2011 during Basic Science flights (Gomez-Ruiz et al. 2012).

The receiver was connected to a digital FFT spectrometer (Klein et al. 2012) providing a bandwidth of 1.5 GHz, and resulting in the spectral resolutions given in Table 1. The observations were performed in double beam-switching mode, with an amplitude of 80″ (or a throw of 160″) at the position angle of 135°, and a phase time of 0.5 s. The nominal focus position was updated regularly against temperature drifts of the telescope structure. The pointing was established with the optical guide cameras to an accuracy of ~5″. The data were calibrated with the KOSMA/GREAT calibrator (Guan et al. 2012), removing residual telluric lines, and further processed with the GILDAS/CLASS software³. We caution that the illumination of the telescope changed between Cycle 0 and Cycle 1. The telescope efficiencies have been recently revised. We adopted the values displayed on the SOFIA website⁴ (see also Table 1).

2.2. IRAM 30 m

The emission of the Cep E outflow in the CO transitions J = 1–0 at 115.271 GHz, J = 2–1 at 230.538 GHz and ¹³CO J = 2–1 at 220.399 GHz was observed with the IRAM 30 m telescope in January 2011. Deep integrations were carried out using the EMIR receivers towards HH377, targeting the position α(2000) = δ(2000) = +61°42′08″. The observations were carried out in Wobbler Switching Mode, with a throw of 3′, and using the broadband EMIR receivers. We used the FTS spectrometers in its 200 kHz resolution mode, corresponding to velocity resolutions of 0.6 km s^-1, for the J = 1–0 observations; the WILMA spectrometer with 2MHz resolution was used for the observations of the ¹²CO and ¹³CO J = 2–1, providing a velocity resolution of about 2.7 km s^-1. The data were processed using the GILDAS/CLASS software.

We obtained a fully sampled map of 2′ × 2′ centered on the protostellar source CepE-mm with the receiver A230 in July 1996, using the On-The-Fly mode. The map of the CO J = 2–1 emission integrated in channels of 3.25 km s^-1 is displayed in Figs. A.1, A.2. This permitted us to obtain the spectrum of the CO J = 2–1 line at 230.538 GHz towards the same position as the J = 1–0 observations in a beam of 20″ size (HPBW), hence comparable to the resolution of the Herschel CO J = 9–8 and the IRAM CO J = 1–0 observations. We present in Fig. 3 the spectrum of the CO J = 2–1 emission in a beam of 20″ half-power beamwidth (HPBW). We adopted the forward and beam efficiency values F_eff and B_eff, as monitored by IRAM⁵ (see Table 1).

2.3. PdBI

The J = 2–1 line emission in the Cep E protostellar outflow was observed with the five-antenna array of the Plateau de Bure Interferometer in March 1997 in configurations C and D, and again in January and February 1998 in configurations A and B. In order to cover both lobes of the outflow, a mosaic of seven fields overlapping at half power was observed, covering an area of approximately 24″ × 90″. The size of the synthesized beam is 1.07″ × 0.87″ and makes a PA of 73deg. We used a spectral correlator configured to give 128 channels of 2.5 MHz (3.25 km s^-1 at the frequency of the CO J = 2–1 line), allowing coverage of the full velocity range of outflow emission. Calibration was carried out following standard procedures, using the GILDAS/CLIC software⁶.

Fig. 3 Montage of all the CO transitions observed towards HH377. Line profiles are drawn in thick black. Spectra enhanced with a magnifying factor in order to enhance the high-velocity emission are drawn in thin red. Fluxes are expressed in units of antenna temperature , corrected for atmospheric absorption (for ground-based observations). The ambient cloud velocity (vlsr = −10.9 km s-1), the peak velocity of the protostellar jet (v = −125 km s-1), and the bowshock HH377 (v = −68 km s-1) are marked by dashed blue lines. On the right hand side of each panel we have indicated the observatory and the HPBW in arcsec.

Our IRAM 30 m map was used to provide the short spacings for our PdBI map. We checked that the flux of the jet (both blue- and redshifted components) and of HH377 were unchanged to better than 97% and 87%, respectively, when including the short spacings, i.e. they were fully recovered by the interferometer. By comparison, 52% of the flux from the southern outflow cavity was recovered by the interferometer.

We present here the interferometric data uncorrected for the missing flux, which therefore affects only the outflow cavity, with the motivation of highlighting the structure of the outflow itself, as can be seen in Figs. 1, 2. We note that the global properties of the outflow cavity (mass, density, temperature, column density, momentum) were all derived from the single-dish observations. Hence, the flux filtered out by the interferometer does not affect the qualitative analysis presented here.

2.4. JCMT

Observations of the CO lines J = 3–2 345.796 GHz and J = 4–3 461.041 GHz were carried out at the JCMT in June 1997. The outflow emission was mapped at 6″ sampling in the J = 3–2 line over an area 42″ × 60″, centered on the position of the protostar Cep E-mm. The map of the CO J = 3–2 emission integrated in channels of 3.25 km s^-1 is shown in Figs. A.3–A.4. At the frequency of the J = 3–2 line, the HPBW of the telescope is ≈ 15″, and the main-beam efficiency of the telescope is 0.67. This map permitted us to obtain the spectrum of the CO J = 3–2 line in a beam of 20″ (HPBW). A deep integration was performed towards the position α(2000) = 23^h03^m11.1^sδ(2000) = +61°42′08.0″, at offset position (−12″−18″) with respect to the protostar Cep E-mm. The observations were carried out in Wobbler Switching Mode with a throw of 3′ using receiver RxB3 and the Digital Autocorrelating Spectrometer (DAS) in its 625 kHz resolution mode corresponding to a velocity resolution of ≈ 0.54 km s^-1. The data were converted from FITS to GILDAS format, and subsequently reduced with CLASS.

2.5. Herschel

The emission in the CO transitions J = 5–4 576.268 GHz, J = 9–8 1036.912 GHz, J = 14–13 1611.794 GHz, and ¹³CO J = 5–4 550.296 GHz was observed with the HIFI spectrometer on board the Herschel Space Observatory on 30 January and 2 February 2013, in the course of the program SPECSO, dedicated to spectroscopy of outflow shocks (PIs: Lefloch & Benedettini). The position targeted was α(2000) = 23^h03^m11.1^sδ(2000) = +61°42′08.0″, which is offset from the SOFIA position by −3.5″ in RA and −5″ in declination (see Fig. 2). The Obs_Id of the J = 5–4, 9–8, 14–13 observations are 1342262584, 1342262801, 1342262555, respectively.

Observations were carried out in double beam switching mode with a throw of 3 arcmin. The receiver was tuned in double sideband and the Wide Band Spectrometer was used, providing a spectral resolution of 1.1 MHz, which was subsequently degraded to reach a final velocity resolution of 0.5 km s^-1. The telescope parameters (F_eff, B_eff, HPBW) were adopted from Roelfsema et al. (2012; see Table 1). The data were processed with the ESA-supported package Herschel interactive processing environment (HIPE⁷). FITS files from level 2 data were then created and transformed into GILDAS format for baseline subtraction and subsequent data analysis.

3. Results

Our large-scale maps of the CO J = 2–1 and 3–2 emission reveal a double outflow structure centered on the protostellar core (Figs. A.1–A.4): the compact (size ~1′), collimated high-velocity outflow, oriented in the southwest-northeast direction, and an extended (size ~3′), less collimated, lower velocity outflow, propagating in the southeast-northwest direction. The relative postion angles of both outflows differ by 52 deg. We do not find evidence for source multiplicity down to 1″ (730 AU) in our PdBI map. This double outflow structure was reported for the first time by Ladd & Hodapp (1997). These authors noticed that the high-velocity outflow is associated with strong line emission from vibrationally excited molecular hydrogen, while the large-scale outflow structure shows no associated molecular hydrogen emission. In this work, we present a detailed study of the collimated, high-velocity outflow, hereafter “the outflow”.

3.1. Multiple outflow components

Our CO J = 2–1 interferometric observations reveal three spatial components in the protostellar outflow: the molecular jet, the outflow cavity, and the terminal shock HH377. In Fig. 1, we superpose the emission of each outflow component on an IRAC/Spitzer (channel 2) image at 2″ resolution of the H emission (Noriega-Crespo et al. 2004). IRAC channel 2 at 4.5 μm is an excellent tracer of the H pure rotational emission, since three of the brightest lines, 0–0 S(11) 4.18 μm, 0–0 S(10) 4.40 μm, and 0–0 S(9) 4.18 μm fall within its passband. Figure 2 shows a superposition of the emission of the three components, contributing to the emission of the southern outflow lobe. A montage of all the CO transitions observed towards Cep E with IRAM 30 m, JCMT, Herschel, and SOFIA is presented in Fig. 3. Each of the components identified in the PdBI map appears to have a specific spectral signature.

The molecular jet powered by the protostar is detected between v_lsr = −140 and −100 km s^-1 (Fig. 3). The emission is detected as a line typically 15 km s^-1 wide and well separated from the bipolar outflow wing emission. Our PdBI map shows that the jet is approximately 21″ long (1.5 × 10⁴ AU, beam-deconvolved), with a transverse size of , or 1200 AU (beam-deconvolved). The jet emission terminates at the offset position (−9″; −23″), at upstream of the nominal position of HH377 (Fig. 1). There is no evidence of highly excited H emission, as traced by IRAC, along the jet.

The blueshifted wing of the bipolar outflow is detected from the ambient cloud velocity (v_lsr = −10.9 km s^-1; Lefloch et al. 1996) up to −100 km s^-1 in the low-J transitions of CO (Figs. 1–3). In the PdBI map, the emission arises from the outflow cavity walls, which appear to tightly surround the protostellar jet. The emission is distributed into two imbricated shells, which provides evidence for multiple ejection. The bulk of the emission is detected up to 10 km s^-1 from the ambient cloud velocity (Fig. 3). Comparison with the IRAC/Spitzer image shows an excellent match between the CO J = 2–1 and H line emission from the outflow cavity walls, as traced by IRAC2 and the PdBI. The outflow cavity can be modeled by an ellipse with a major (minor) axis of 22″ (10″) (beam-deconvolved) with a PA of 70deg; see Fig. 1. In what follows, we refer to the outflow cavity as the gas that is contained within the cavity walls, which is not part of the molecular jet, such as detected in our PdBI map.

The terminal shock HH377 is seen downstream of the molecular jet, close to the offset position (−9″, ); it is detected in the velocity range −90; −50 km s^-1, i.e. at intermediate velocities between the outflow cavity and the jet. At the PdBI, in CO J = 2–1, it is detected as a bright component downstream of the molecular jet, with a typical size of ≈ 4″ (beam-deconvolved). The CO J = 2–1 emission region coincides with the rear of HH377, as observed with IRAC/Spitzer (Fig. 1). This suggests that the CO J = 2–1 line traces the post-shocked gas layer, whereas the pure H lines trace the high-excitation gas in the shock.

In the same velocity range, two weaker components are detected in addition to HH377. Both display a compact size of a few arcsec. The first coincides with the clump of shocked gas reported by Lefloch et al. (2011) in the inner envelope of the protostar; the second coincides with the infrared knot BI reported by Gomez-Ruiz et al. (2012).

The emission from HH377 appears as a bright “bump” in the highest excitation lines of CO (J_up = 16), whose intensity is similar to that of the entrained gas at v = −20 km s^-1. The intensity of the HH377 emission decreases in the lower-excitation lines of CO, and appears to be fully dominated by the entrained gas in the transitions J_up ≤ 4 (Fig. 2). This is consistent with highly-excited emission from a small region.

3.2. Dynamical age

Noriega-Crespo et al. (2014) have measured the proper motions in the jet of Cep E using IRAC images containing mostly a mix of emission of several H₂ lines obtained over a period of 6 years. They derived tangential velocities of the northern and southern outflow lobes of 62 ± 29 and 94 ± 26 km s^-1, respectively.

Adopting the value measured in our CO observations (see above) as jet radial velocity, we derive both the jet velocity V_j and the jet inclination with respect to the line of sight. We find a jet velocity of +107 km s^-1 and −148 km s^-1) in the northern and southern lobes, respectively. We find in a consistent way that the jet makes the same angle Φ = 40deg with respect to the line of sight in both outflow lobes. The difference in the jet velocity could arise from a density gradient in the cloud, with a higher density towards the north.

We have estimated the dynamical age of the terminal bowshock HH377, located at 26″ (0.092 pc) from the protostar, as τ_dyn = d/V_jsinΦ where d is the distance of HH377 to the location of the protostar; it comes τ_dyn = 960 yr. The small size of both outflow lobes and the young dynamical age suggest that the jet is probably at an early evolutionary stage.

4. Outflow physical conditions

4.1. Methodology

The determination of the physical conditions in the outflow is made difficult by the presence of the second large-scale outflow (see above). As discussed in Sect. 3, the CO emission of both outflows overlap between +2 and −17 km s^-1 in the central 30″ around the protostar (Figs. A.1 to A.4). We note that the CO J = 2–1 emission from the large-scale outflow is detected up to −30 km s^-1. Therefore, it is difficult to disentangle the respective contributions of both outflows to the entrained gas (the outflow cavity) in the central 30″ around the protostar. The absence of H emission associated with the large-scale outflow suggests, however, low excitation in the latter, so that the contamination most likely affects only the low-J transitions. In order to accurately derive the excitation conditions in the CO cavity gas, we have used the line fluxes integrated between −50 and −17 km s^-1, outside the velocity interval in which the outflows overlap.

The CO gas excitation temperature and column density were first derived from a simple rotational diagram analysis of the line fluxes, after correcting the intensities for the main-beam filling factor (see Goldsmith & Langer 1999). The CO fluxes from the molecular jet and the outflow cavity emission were obtained from direct integration in the corresponding velocity range (see Table 2). The emission from the shock was obtained after subtracting the contribution of the entrained gas. As the shock emission is dominated by the entrained gas emission in the transitions J_up = 1–4, only the emission in the transitions (J_up ≥ 5) is reported here and taken into account in the analysis.

For each component of the outflow (jet, HH377, cavity), we estimated the solid angle of the emission encompassed by the telescope beam at the frequency of each CO transition, from the intersection between the half-power contour of the J = 2–1 emission mapped at the PdBI and the telescope main-beam modeled by a disk of diameter equal to the HPBW, centered at the targeted position. We then derived the corresponding beam filling factor from the ratio of the source solid angle to the telescope main beam in order to estimate the intrinsic brightness temperature of the three outflow components in each CO transition. The results of this procedure were successfully compared to the value obtained from a χ² minimization of the CO line emission modeled with a radiative transfer code in the large velocity gradient approximation, as explained in Sect. 4.3; we feel confident that the procedure is valid. The integrated intensities are indicated in Tables 1 and 2, while the results are given in Table 3.

4.2. Line opacity

We used the ¹³CO lines J = 2–1 and 5–4 in order to constrain the opacity of the transitions of the main isotopologue in the outflow. The spectra of the ¹³CO transitions are presented in the bottom panels of Fig. 4. The line profiles peak at ambient velocity v_lsr = −10.9 km s^-1. The linewidth is Δv ≃ 2.7 km s^-1. Emission is detected at the 3σ level in the entrained gas up to −50 km s^-1. We note that the absorption features detected between −5 and 0 km s^-1 in the J = 2–1 spectrum are caused by contamination from the reference position adopted in the observations (see Sect. 2.2).

Fig. 4 Bottom: 13CO J = 2–1 and 5–4 emission detected towards HH377. The dashed lines (in blue) mark the peak velocity of the cloud emission (v = −10.9 km s-1); we have indicated the observatory and the telescope HPBW in arcsec on the top right hand side. Top: intensity line ratio of 13CO/12CO J = 2–1 and 5–4 as a function of velocity, between −70 and +5 km s-1.

The ratio of the ¹²CO to ¹³CO line intensities is shown in the top panel of Fig. 4. A similar behavior is observed for both transitions. The ratio is strongly peaked at the ambient cloud velocity with values ≈ 0.4 and 0.26 for the J = 2–1 and J = 5–4, respectively. These values refer to the ambient cloud molecular gas, however, and not to the entrained gas of the outflow. As the J = 5–4 line ratio is sampled at 1 km s^-1 resolution, it appears much more strongly peaked than the J = 2–1 data, which are sampled at 3 km s^-1 only. At v = −14 km s^-1 and −17 km s^-1, we measure a ratio of 0.10 and 0.032 for the J = 2–1, respectively, and 0.047 and 0.032 for the J = 5–4. Adopting a standard abundance ratio [¹²CO]/[¹³CO] = 75 (Wilson & Rood 1994), these ratios correspond to line opacities of and 2.1 ( and 2.2), respectively. Both line ratios decrease sharply with increasing velocity, and vary little for velocities beyond −20 km s^-1, with typical values of 0.015 and 0.22 for the J = 2–1 and 5–4 transitions, respectively. Overall, both transitions can be regarded as optically thin with and over this velocity range. The ratio appears to decrease again for velocities greater than −40 km s^-1, corresponding to lower optical depths. The signal-to-noise ratio of the ¹³CO data is not sufficient to accurately estimate the line opacity.

Therefore, the line opacity decreases strongly with increasing velocity, a result similar to that obtained by van Kempen et al. (2009a). We conclude that the J = 2–1 and 5–4 line emissions from the outflow cavity are optically thin, except for velocities close to ambient, between −15 and −10.9 km s^-1.

Because of the spatial overlap between the two outflows from Cep E, with rather different excitation conditions, it must be kept in mind that these opacities may not reflect the optical depth of the emission from the high-velocity outflow. This point is addressed in Sect. 4.4.3.

Fig. 5 Rotational diagram of the CO emission in the three components of the outflow: the jet (top), the terminal shock (middle), and the outflow cavity (bottom).

4.3. Rotational diagram analysis

In a first step, we estimated the CO excitation conditions in each outflow component from a simple rotational diagram analysis. The CO line fluxes corrected for the main-beam filling factor are given in Table 2. The rotational diagrams are displayed in Fig. 5. The CO column densities derived in this analysis are therefore values averaged over the source emitting region.

In the jet (top panel of Fig. 5), the flux distribution of the CO lines can be approximately modeled by two components. A component with a rotational temperature T_rot ≃ 70 K and column density N(CO) = 9 × 10¹⁶ cm^-2 accounts for the emission in the range J_up = 1–5. It fails, however, to account for the emission of the higher-excitation lines J_up ≥ 9. The emission of this line is well fit by a second gas component of higher excitation with T_rot = 240 K and of lower column density N(CO) ≃ 2 × 10¹⁶ cm^-2.

In the outflow cavity, a similar pattern with two excitation components is observed (bottom panel of Fig. 5). The flux distribution of the low-J CO lines is well described by a gas component with T_rot = 55 K and N(CO) = 2.5 × 10¹⁷ cm^-2, whereas the emission of the higher-J lines arises from with gas with T_rot = 200−320 K and N(CO) ≃ (5.6−7.6) × 10¹⁶ cm^-2. The fit to the high-J lines suffers large uncertainties and we cannot exclude the presence of higher excitation gas components. We speculate that this emission most likely arises from the shocked gas in the cavity walls, as suggested by the good match between the spatial distributions of the H pure rotational lines and the CO J = 2–1 line emission, obtained with Spitzer and the PdBI, respectively.

The two regimes of excitation detected in the jet and the outflow cavity intersect at E_up ≃ 150–200 K; hence, the millimeter (J_up ≤ 3) and submillimeter (J_up ≥ 9) transitions of CO probe physically different regions both in the jet and in the entrained gas. Therefore, great care must be taken when modeling the emission of the different CO transitions as it can introduce severe bias in the analysis (Gomez-Ruiz et al. 2012).

Last, in HH377, only one gas component of column density N(CO) ≃ 1.3 × 10¹⁷ cm^-2 and T_rot = 320 K is detected.

4.4. Large velocity gradient analysis

In a second step, we have better constrained the CO column density, the H density, and the kinetic temperature in the outflowing gas, using the radiative transfer code MADEX in the large velocity gradient (LVG) approximation. We used the CO-H collisional coefficients of Yang et al. (2010) and built a grid of models with density between 10⁴ cm^-3 and 10⁷ cm^-3 and temperature between 30 K and 1000 K to determine the region of minimum χ² as a function of density and temperature. We adopted a typical linewidth Δv of 13 and 15 km s^-1 for the jet and the shock, respectively, as derived from a simple Gaussian fit to the CO lines. Since the CO line profiles from the outflow cavity are far from a Gaussian line, we analyzed the emission by velocity intervals between −20 and −17 km s^-1, and from −50 to −20 by interval of 5 km s^-1. The physical conditions of the three outflow components are summarized in Table 2.

4.4.1. The jet

We first attempted, but failed, to find a reasonable fit of all the CO line fluxes with one single set of parameters (n(H₂), T_kin, N(CO)). In agreement with our rotational diagram analysis (see top panel of Fig. 5) which suggest two regimes of excitation conditions, we modeled the emission of the low- (J = 1–0 to 5–4) and high- (J = 9–8, 12–11, 13–12, 14–13, 16–15) excitation transitions separately. The results are shown in the top panel of Fig. 6.

In the low-excitation gas, the best-fit solution (χ² = 0.54) is obtained for T_kin = 100 K and n(H₂) = 5 × 10⁴ cm^-3, and N(CO) = 10¹⁷ cm^-2, close to the local thermodynamical equilibrium (LTE) determination (see Sect. 4.2). We find acceptable solutions for gas kinetic temperatures in the range 80–100 K and H jet density ≈ (0.5−1.0) × 10⁵ cm^-3. Solutions with a much higher temperature ≥ 200 K and lower density ≤ 2 × 10⁴ cm^-3, although possible, are less likely, as indicated by the higher values of χ², on the order of 0.8. modeling predicts that all the transitions are optically thin; the highest opacities are found for the J = 5–4 and 6–5 transitions with τ ≃ 0.4. This is consistent with our results on the opacity of the J = 2–1 transition in Sect. 4.1.

Fig. 6 χ2 distribution of LVG slab models for the three components of the outflow. Top: the jet: contour levels at 0.6, 0.7, 0.8 and at 0.55, 0.65, 0.75 are drawn for the low-excitation component () and the high-excitation component (), respectively. Middle: HH377: contour levels are 1.25, 1.5, and 2.5 (). Bottom: the cavity: only emission between −20 and −17 km s-1 is considered. Contour levels at 0.08, 0.12, 0.15 and 0.6, 1.2, 2.4 are drawn for the high- and low-excitation components, respectively ( and 0.6, respectively).

A second component at higher temperature is needed in order to account for the emission of the high excitation lines (see Fig. 6). The best-fit solution (χ² = 0.44) is obtained for a gas column density N(CO) = 1.5 × 10¹⁶ cm^-2 at T_kin = 650 K and n(H₂) = 3 × 10⁴ cm^-3. The physical conditions are less well constrained as solutions with T_kin in the range 400–750 K and n(H₂) in the range (0.5−1.0) × 10⁶ cm^-3 are also possible. All the lines are extremely optically thin (τ ≪ 1). For solutions in the high-temperature regime, the gas density is comparable to the density of the lower-excitation jet component, ≈ 10⁵ cm^-3. As discussed below in Sect. 5, we propose that the high-excitation CO emission arises from internal shocks along the jet. Indeed, complementary observations from specific shock tracers in the jet favor solutions in the high-density–low-temperature regime (Lefloch et al., in prep.).

4.4.2. HH377

The CO fluxes of the transitions between J = 5–4 and J = 16–15 are well reproduced by one gas component of uniform density and temperature. The best-fit solution (χ² = 1.1) is found for a gas column density N(CO) ≃ 10¹⁷ cm^-2 at T_kin = 460 K and density n(H₂) ≃ 10⁶ cm^-3. We note that solutions with T_kin temperature as high as ≈ 1000 K, and a density as low as ≈ 2.5 × 10⁵ cm^-3, with χ² values ≤ 2.5, can also account for the CO line fluxes. Good solutions (χ² ≤ 1.3) are found for T_kin in the range 400–500 K and n(H₂) in the range (1−2) × 10⁶ cm^-3. All the CO transitions are predicted to be optically thin with τ ≪ 1.

Based on observations of the [OI] lines λ6300, 63μm and [SII] line λλ6717 / 31 lines, Moro-Martin et al. (2001) concluded that HH377 is associated with a dissociative J-type shock with an estimated pre-shock density of ~10⁵ cm^-3 (Ayala et al. 2000). Dissociative shock models (Hollenbach & McKee 1989; Gusdorf et al. 2015) predict the CO emission to arise from the molecular reformation zone in the post-shocked gas, with temperatures in the range 400–500 K. Hence, both previous work and shock models favor solutions with a low temperature (≃400–500 K) and high density ((1−2) × 10⁶ cm^-3).

4.4.3. The cavity

We first considered the flux of the transitions J = 1–0, 2–1, 3–2, and 9–8 between −20 and −17 km s^-1. These transitions were observed towards the same position and at a common angular resolution of ≈ 20″. We ran our LVG code with a χ² minimization procedure varying the main-beam filling factor, CO column density, temperature, and density; the minimum χ² is obtained for a filling factor 0.17−0.20. This value is very close to the determination obtained independently by means of the interferometric map of the Plateau de Bure, as explained in Sect. 4.1 (Table 2); we are confident in the method adopted.

Because the intensity distribution between −50 and −14 km s^-1 is far from a Gaussian distribution, we decided to analyze the CO line intensities integrated over narrower velocity intervals: between −14 and −17, −17 and −20, −20 and −25, −25 and −30, −30 and −35, −35 and −40, and −50 and −40 km s^-1. Again, we found that it is possible to account for the emission of all the transitions from J = 1–0 up to J = 9–8 with one single component. The best fitting solution is rather good in the low velocity range (χ² ≃ 0.45). The best fitting solution cannot account for the intensity of the high-excitation lines J = 12–11 to 16–15. We note that the J = 9–8 flux lies at the intersection of the two regimes of excitation found in the rotational diagram (see Fig. 5). Therefore, we have searched for the best fitting solutions which allow each excitation regime to be reproduced separately, J = 1–0 to 5–4, and J = 12–1 to 16–15, respectively, with the additional constraints that a) the resulting flux for the J = 9–8 was consistent with the observational determination, and b) the opacity of the CO J = 5–4 agrees with the observations. The physical parameters for each velocity interval are summarized in Table 3 and the best fitting solutions for each velocity interval are displayed in Figs. B.1–B.2. We show in Fig. B3 the best fitting solutions for both components and the total fit to the full spectral line flux distribution for each velocity interval.

The emission of the low-excitation gas is well described by a gas layer with a kinetic temperature in the range 55 K to 85 K and density of (1−8) × 10⁵ cm^-3 (Fig. B.1). The CO column density decreases with velocity from 10¹⁷ cm^-2 between −20 and −17 km s^-1 down to 3 × 10¹⁶ cm^-2 between −50 and −40 km s^-1. The H density is reasonably well constrained in the low-velocity range up to −30 km s^-1, with n(H₂) ≃ (1−4) × 10⁵ cm^-3. If we extrapolate the excitation conditions derived between −14 and −17km to the lower velocity gas between −14 and −11 km s^-1, we estimate a typical column density 2 × 10¹⁷ cm^-2. The total gas column density is then N(CO) ≃ 7 × 10¹⁷ cm^-2.

We give in Table 3 the opacities of the J = 2–1 and 5–4 transitions, as predicted by our simple modeling. The agreement with the observational determination is found satisfying for the J = 5–4 transition, with on the order of 1 between −20 and −30 km s^-1. On the contrary, our modeling predicts the J = 2–1 transition to be optically thin even in the low-velocity range, whereas observations indicate between −14 and −20 km s^-1. At higher velocities, the agreement is reasonable. As mentioned in Sect. 4.2, the emission of the low-excitation lines is contaminated by the contribution of the large-scale, lower-excitation outflow. The physical conditions of this outflow were derived by Ladd & Hodapp (1997) from CO J = 1–0 and 2–1 observations, who reported a typical temperature T_kin ≃ 20 K and column density N(CO) ≃ 4 × 10¹⁵ cm^-2. Assuming a typical outflow density of a few 10⁴ cm^-3, LVG calculations of the CO emission shows that the gas emits mainly in the lines J_up = 1 to 4. As a result, in the velocity range −11 to −20 km s^-1, the emission of both the CO and ¹³CO lines could be biased by the lower-excitation emission from the large-scale outflow, whereas the higher-excitation lines, J = 5–4 and higher, would not be affected as they arise mainly from the compact, high-excitation outflow. In all cases, the modeling presented here is very crude; an approach taking into account the variations of temperature and density as a function of location inside the outflow cavity would be far more realistic. Such modeling would be at the cost of many additional parameters, which cannot be constrained by the observations discussed here owing to their limited angular resolution.

The emission of high-excitation lines is accounted for by gas at temperatures ranging from ≃ 500 K up to 1500 K, and density in the range (1−5) × 10⁶ cm^-3 (see Fig. 6; also Fig. B.2). The CO column density decreases with velocity from 6.5 × 10¹⁵ cm^-2 between −20 and −17 km s^-1 down to 8.5 × 10¹⁵ cm^-2 between −50 and −40 km s^-1. Proceeding as above and extrapolating the excitation conditions derived between −14 and −17km to the lower velocity gas between −14 and −11 km s^-1, we obtain a total column density 6 × 10¹⁶ cm^-2. The emission of this component is predicted to be highly optically thin in all the CO transitions: τ ≪ 1. Our χ² minimization analysis suggests a gradient of excitation in the outflowing gas as higher-temperature solutions are favored when the velocity range increases (see Table 4 and Fig. B.2). Complementary observations from high-density gas tracers will be useful to better constrain the density and the temperature distributions in the entrained gas. The best fitting temperature and CO column density are in good agreement with the results of Moro-Martin et al. (2001), based on ISO/LWS observations of the CO transitions J = 14–13 up to J = 33–32 at 80″ resolution towards the southern lobe of Cep E. These authors found T_kin = 1200 K and N(CO) = 10¹⁷ cm^-2 from a large velocity gradient analysis. These results are also consistent with the analysis of the pure H rotational lines by the same authors.

Fig. 7 CO velocity-integrated intensity diagrams predicted by non-stationary CJ-type shock models for various pre-shock densities: nH= 104 (upper left panel), 5 × 104 cm-3 (lower left panel), 105 (upper right panel), 5 × 105 cm-3 (lower right panel), and G. The observational data are drawn with filled (empty) squares for a filling factor of 1 (of 0.25). The shock velocity and age are given in each panel, as is the color code.

4.5. Mass and momentum budget

We have estimated the mass contained in the different outflow components as follows.

For the jet, we integrated the CO J = 2–1 emission in the southern jet lobe, as mapped with the PdBI, using the following relation under the hypothesis of LTE: (see Bachiller et al. 1990). This relation is valid only in the optically thin regime and can be applied safely to the jet, as its emission is optically thin (see Sect. 4.1). Adopting T_ex = 70 K (see Sect. 4.3), we obtain 0.02M_⊙ for the jet mass. Given the dynamical age estimated in Sect. 3.2, the mass-loss rate associated with the jet is then Ṁ = 2 × 10^-5 M_⊙yr^-1.

For the outflow cavity, the interval of integration was decomposed into two velocity intervals:

• Between −100 and −50 km s^-1, we applied the above relation to the CO J = 2–1 emission which we corrected for the line opacity τ¹². To do so, we modeled the variations of τ¹² as a function of velocity by a simple synthetic profile and introduced the correcting factor τ¹²/ (1−exp(−τ¹²)) in the previous relation. We find a mass of 0.02M_⊙.

• Between −50 and −11 km s^-1, we adopted the total CO column density obtained in our LVG analysis (Sect. 4.4) and adopted the size of the outflow cavity, as estimated with the PdBI. The mass of the outflowing gas between −50 and −20 is ≈ 0.30M_⊙.

Hence, we estimate a total mass of ≃ 0.32M_⊙ for the outflow cavity gas. Our estimate of the outflow cavity mass is somewhat larger by a factor 6–10 to the previous determinations by Moro-Martin et al. (2001) and Ladd & Hodapp (1997), who found 0.05M_⊙ and 0.03M_⊙, respectively. Such discrepancies can easily be reconciled. These authors adopted a CO excitation temperature of 20 K, much lower than the 55 K estimated in the present work, and neglected any correction for line opacity effects. The discrepancy with the estimate of Moro-Martin et al. (2001) is then ≈ 2, which is satisfying taking into account that their estimate is based on ¹³CO interferometric observations of the outflow, uncorrected for flux filtering and primary beam attenuation.

We proceeded in the same way in order to estimate the radial momentum associated with the jet and the outflow cavity. In the jet, we found 1.7M_⊙ km s^-1. This is very similar to a crude estimate based on the mass of the jet and the mean jet velocity with respect to the ambient cloud. In the outflow cavity, the total radial momentum carried away amounts to ≈ 2.6M_⊙ km s^-1, a value very similar to that in the jet. This very good agreement suggests that the jet momentum is transferred to the entrained gas in the low-velocity outflow with a high efficiency. It supports the idea that the Cep E outflow is driven by the molecular jet through the bowshock HH377, in agreement with the model of Raga & Cabrit (1993).

For the terminal bowshock HH377, the mass and momentum estimates are more uncertain as the J = 2–1 line appears to be dominated by the entrained gas. Hence, we obtained a rough estimate of the mass from the source-averaged CO column density derived from our LVG analysis and the typical size of 5″ estimated from our interferometric data. We found 0.006M_⊙. An estimate of the momentum of 0.4M_⊙ km s^-1 was obtained by adopting a mean value of 60 km s^-1 with respect to the ambient cloud (see also Fig. 3).

Fig. 8 Comparison of the best fitting CJ-type shock models and their predictions for the CO velocity-integrated flux distribution of the protostellar jet. The observational data are drawn with filled (empty) squares for a filling factor of 1 (of 0.25). Upper panel: a relative filling factor of 1 is assumed. Solutions with a pre-shock density nH = 5 × 104 (filled diamonds) and 105 cm-3 (filled circles), νs = 20−30 km s-1, and age in the range ~220–450 yr, with the color code indicated in the figure are found to provide an equally good fit. Lower panel: a relative filling factor of 0.25 is assumed. Solutions with a pre-shock density nH = 105 cm-3 (empty circles), νs = 25–30 km s-1, and age in the range ~595–740 yr, with the color code indicated in the figure, provide a good fit to the data.

5. The jet structure

In the jet, the high-J CO transitions (J_up ≥ 9) trace mainly the warm gas component (Fig. 5). Our PdBI map of the CO 2–1 emission reveals the presence of bright emission knots along the jet (Figs. 1 and 2). The main-beam solid angle of the high-J CO lines (J ≥ 9) encompasses two such knots, peaking at δ = −16″ and δ = −20″, respectively, which present a typical size of 2″–4″ (see Figs. 1, 2). Variations of the centroid of emission and of the line profiles are also observed, on the order of 10 km s^-1, showing evidence of the presence of internal motions and a complex small-scale structure inside the jet. Interestingly, the presence of warm (T_kin ≃ 200 K), dense (≃ 10⁶ cm^-3) clumps of gas of arcsec size, associated with H₂O and SiO emission was reported by Lefloch et al. (2011) along the protostellar jet. These previous results suggest that the emission of the high-J CO lines could arise from shocked clumps of gas along the jet.

We tested this hypothesis by comparing the CO flux spectral distribution of the hot jet component with the predictions of a small grid of shock models, following a method extensively presented in e.g. Gusdorf et al. (2011, 2012). In this simplified approach, clumps of shocked gas propagating along the jet are compared to one-dimensional shock models. This approach is justified by the lack of information on the internal structure of the jet. The fluxes of the low-J transitions were estimated from our LTE modeling (Sect. 4.2). The grid consists of steady-state J-type and C-type models from Gusdorf et al. (2015) and also from a set of non-stationary models (Lesaffre et al. 2004) with input parameters adapted to the present region: pre-shock nuclei density n_H = 10⁴, 5 × 10⁴, 10⁵, 5 × 10⁵ cm^-3; magnetic field (transverse to the shock direction) given by G; and shock velocity ν_s = 20, 25, 30 km s^-1, well in the range of internal motions detected in the PdBI data. The non-stationary shock models present combined characteristics of both stationary C- and J-type models, and are consequently called “CJ-type” shock models. For these models, the age is an additional input parameter; we chose to explore an age less than 1000 years, the dynamical age of HH 377.

The observational CO spectral flux distribution reaches its maximum at J_up = 9. There are no J-type shock models that can account for the CO flux distribution: the corresponding shock layers are too thin to generate sufficient CO emission to match the data. We searched for a CJ-type solution to fit the CO line fluxes given in Table 2, i.e. we assumed that the high-excitation component has the same filling factor as the low-excitation gas, such as that mapped with the PdBI map (relative filling factor of 1). The model grid is displayed in Fig. 7 and the observational fluxes are drawn with filled squares. We find that CJ-type shocks propagating into a pre-shock gas of density n_H = (0.5−1.0) × 10⁵ cm^-3, corresponding to n(H₂) = (2.5−5) × 10⁴ cm^-3 and B = (0.2−0.3) mG, with a velocity v_s ≃ 20–30 km s^-1 match rather well the observed CO fluxes (see bottom left and upper right panel in Fig. 7, and top panel in Fig. 8). Such a pre-shock density is consistent with our determination of the density in the bulk material of the jet (Table 2). The pre-shock density is rather well constrained: for pre-shock gas densities lower than n_H = 10⁴ cm^-3 (upper left panel in Fig. 7), CJ-type shocks fail to produce bright enough emission in the high-J CO lines. On the contrary, higher pre-shock gas densities n_H = 5 × 10⁵ cm^-3 (bottom right panel in Fig. 7) tend to produce CO lines much brighter than those observed, by a factor of about 5. For a filling factor of 1, the best agreement is obtained for CJ-type shocks at a relatively early age of 220–450 yr, as is shown in the upper panel in Fig. 8. Our models predict that the density of the CO emitting layer lies in the range from a few 10⁵ to 10⁶ cm^-3, which is consistent with the physical conditions derived from our LVG analysis (Sect. 4.3).

In what precedes, we assumed that the shocked gas has the same filling factor as the bulk material of the CO jet. The size of the shock is not well constrained and it could be smaller than the size of the 70 K component, which dominates the J = 2–1 jet emission. In this case, the CO line emission from the shock would be brighter than estimated above (Table 2). In Fig. 7 we have superposed the CO line fluxes to our model predictions under the assumption that the relative filling factor (FF) is 0.25 (empty square symbol) and have compared these data points to our models. A filling factor of 0.25 corresponds to a typical size 4″–5″ for one component, or ≃ 2″ for two components of similar excitation. Such sizes are compatible with the size of the CO knots in the jet, such as those mapped with the PdBI. We have also chosen to display the best fitting results for this value of the filling factor in the lower panel of Fig. 8. This figure shows that slightly older shock models are then favored (age 595–740 yr), with a pre-shock density more tightly constrained to n_H = 10⁵ cm^-3.

Complementary or additional CO line observations (J_up = 6, 7, 11, 12 are possible from the ground or from SOFIA) would help to more tightly constrain our shock models. Interferometric observations at subarcsec scale of shock tracers such as SiO (see e.g. Lefloch et al. 1996, 2011) would help to constrain the size of these shocks and to provide a precise structure of the protostellar jet.

6. Conclusions

We have carried out an unprecedented multi-transition study at high spectral resolution of the CO emission in the southern outflow lobe of the intermediate-mass protostar Cep E-mm, based on observations of transitions in the range = 1–16. PdBI observations at arcsec scale of the J = 2–1 line show that the “high-velocity molecular bullets” reported in a previous work are actually tracing a protostellar jet and the terminal bowshock HH377 where the jet impacts the molecular cloud. The dynamical age of the jet is short, ≈ 960 yr, suggesting that the jet is still at an early evolutionary stage.

We have constrained the physical conditions in three outflow components from an LTE and LVG analysis:

• In the jet, of size 21″ × 1.7″, the low-J CO emission is dominated by a gas layer of column density N(CO) = 9.0 × 10¹⁶ cm^-2, density n(H₂) = (0.5−1) × 10⁵ cm^-3 and temperature T_kin = 80–100 K; the high-J CO emission is well accounted for by a gas layer of lower column density N(CO) = 1.5 × 10¹⁶ cm^-2, higher density (n(H₂) = (0.5−1.0) × 10⁶ cm^-3) and higher temperature (T_kin = 400–750 K).

• In the terminal bowshock HH377, of size 4.5″, the CO emission is well described by one single gas component of dense gas (n(H₂) = (1−2) × 10⁶ cm^-3) with temperature in the range 400–500 K and column density N(CO) = 10¹⁷ cm^-2.

• In the outflow cavity, of size 22″ × 10″, the low-J CO emission is dominated by a gas layer of column density N(CO) ≃ 7 × 10¹⁷ cm^-2, density n(H₂) ≃ (1−8) × 10⁵ cm^-3, and temperature T_kin ≃ 50−85 K; the high-J CO line emission is accounted for by a gas layer of lower column density N(CO) = 6 × 10¹⁶ cm^-2, higher density n(H₂) = (1−5) × 10⁶ cm^-3, and higher temperature T_kin = 500–1500 K. Evidence of a gradient of excitation in the high-excitation gas suggests a complex thermal structure in the outflow cavity.

Overall, the jet and the outflow cavity appear to carry away the same amount of angular momentum, which supports Cep E as a convincing case of jet bowshock driven outflow.

The presence of hot CO gas inside the jet can be consistently explained by the presence of internal CJ-shocks propagating at a velocity of ~20–30 km s^-1 along the jet. This implies that the jet is magnetized with typical magnetic field intensity (transverse to the shock direction) on the order of 200–300 μG. Interferometric observations at subarcsec scale of specific shock tracers with NOEMA would help constrain the structure of the protostellar jet and the bowshock HH377.

Acknowledgments

We thank the referee Dr. Tim van Kempen and the Editor Malcolm Walmsley for detailed comments which have helped to improve the manuscript. We also thank the SOFIA operations and the GREAT instrument teams, whose support has been essential for the GREAT accomplishments, and the DSI telescope engineering team. Based on observations made with the NASA/DLR Stratospheric Observatory for Infrared Astronomy, with the IRAM 30 m telescope and the IRAM Plateau de Bure Interferometer. SOFIA Science Mission Operations are conducted jointly by the Universities Space Research Association, Inc., under NASA contract NAS2-97001, and the Deutsches SOFIA Institut, under DLR contract 50 OK 0901. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). The James Clerk Maxwell Telescope is operated by the Joint Astronomy Centre on behalf of the Science and Technology Facilities Council of the United Kingdom, the Netherlands Organisation for Scientific Research, and the National Research Council of Canada. This work was partly supported by the CNRS program “Physique et Chimie du Milieu Interstellaire”, by the PRIN INAF 2012 – JEDI and by the Italian Ministero dell’Istruzione, Università e Ricerca through the grant Progetti Premiali 2012 – iALMA A.I. Gómez-Ruiz is supported by Consejo Nacional de Ciencia y Tecnología, through the program Cátedras CONACYT para Jóvenes Investigadores.

References

Online material

Appendix A: Complementary single-dish CO observations

Fig. A.1 Map of the CO J = 2–1 emission integrated in channels of 3.25km from −13.5 to +113.1 km s-1, as observed with the IRAM 30 m telescope. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (10.7″) is indicated by the white circle.

Fig. A.2 Map of the CO J = 2–1 emission integrated in channels of 3.25km from −137.0 to −10.3 km s-1, as observed with the IRAM 30 m telescope. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (10.7″) is indicated by the white circle.

Fig. A.3 Map of the CO J = 3–2 emission integrated in channels of 3.25km from −13.5 to +113.1 km s-1, as observed with the JCMT. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (15″) is indicated by the white circle.

Fig. A.4 Map of the CO J = 3–2 emission integrated in channels of 3.25km from −137 to −10.3 km s-1, as observed with the JCMT. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (15″) is indicated by the white circle.

Appendix B: Excitation conditions in the outflow cavity

Fig. B.1 χ2 distribution of LVG slab models for different velocity intervals in the low-excitation gas of the outflow cavity. The considered velocity intervals are indicated in each panel. The minimum value and the contour levels are indicated in each panel.

Fig. B.2 χ2 distribution of LVG slab models for different velocity intervals in the high-excitation gas of the outflow cavity. The considered velocity intervals are indicated in each panel. The minimum value and the contour levels are indicated in each panel.

Fig. B.3 Best fitting solutions of the CO spectral line flux distribution as a function of the rotational number of the transition from our LVG analysis with the MADEX radiative transfer code. The fits to the low- and high-excitation components are shown by black solid and dashed lines, respectively. The fit to the total flux is shown in red.

All Tables

All Figures

Fig. 1 Cep E outflow emission observed in the CO J = 2–1 line with the PdBI at 1″ resolution (contours) and with IRAC/Spitzer H24.5 μm line (color). Each panel shows the emission associated with each of the three outflow components identified: (left) outflow cavity walls emission, integrated between −15.2 and −8.7 km s-1 in the southern lobe (blue contours) and between −2.2 and −8.7 km s-1 in the northern lobe (white contours); (center) jet emission integrated between −135 and −110 km s-1 and between +40 and +80 km s-1 towards the southern (blue contours) and northern lobes (white contours), respectively; (right) HH377 bowshock emission integrated between −77 and −64 km s-1 (blue contours). First contour and contour interval are 20% and 10% of the peak intensity, respectively.

In the text

Fig. 2 Left: CO J = 2–1 velocity-integrated emission of the protostellar jet observed with the PdBI. The emission is integrated between −135 and −110 km s-1 and between +40 and +80 km s-1 towards the southern and northern lobes, respectively. Right: the three kinematic components in the southern lobe of the Cepheus E outflow as seen by the PdBI at 1″ resolution in the CO J = 2–1 line: a) outflow cavity emission between −15 and −12 km s-1 (color); b) HH377 bowshock emission between −77 and −64 km s-1 (yellow contours); c) jet emission between −135 and −110 km s-1 (blue contours). First contour and contour interval are 20% and 10% of the peak intensity, respectively. The white stars mark the location of the positions observed with SOFIA (−8.5″, −23″) and Herschel (HSO), IRAM, and JCMT (−12″, −18″). The beam width at half-power of SOFIA in the CO J = 16–15 and of HIFI in the J = 5–4 are drawn by a solid and a dotted circle, respectively.

In the text

Fig. 3 Montage of all the CO transitions observed towards HH377. Line profiles are drawn in thick black. Spectra enhanced with a magnifying factor in order to enhance the high-velocity emission are drawn in thin red. Fluxes are expressed in units of antenna temperature , corrected for atmospheric absorption (for ground-based observations). The ambient cloud velocity (vlsr = −10.9 km s-1), the peak velocity of the protostellar jet (v = −125 km s-1), and the bowshock HH377 (v = −68 km s-1) are marked by dashed blue lines. On the right hand side of each panel we have indicated the observatory and the HPBW in arcsec.

In the text

	Fig. 4 Bottom: 13CO J = 2–1 and 5–4 emission detected towards HH377. The dashed lines (in blue) mark the peak velocity of the cloud emission (v = −10.9 km s-1); we have indicated the observatory and the telescope HPBW in arcsec on the top right hand side. Top: intensity line ratio of 13CO/12CO J = 2–1 and 5–4 as a function of velocity, between −70 and +5 km s-1.
In the text

	Fig. 5 Rotational diagram of the CO emission in the three components of the outflow: the jet (top), the terminal shock (middle), and the outflow cavity (bottom).
In the text

Fig. 6 χ2 distribution of LVG slab models for the three components of the outflow. Top: the jet: contour levels at 0.6, 0.7, 0.8 and at 0.55, 0.65, 0.75 are drawn for the low-excitation component () and the high-excitation component (), respectively. Middle: HH377: contour levels are 1.25, 1.5, and 2.5 (). Bottom: the cavity: only emission between −20 and −17 km s-1 is considered. Contour levels at 0.08, 0.12, 0.15 and 0.6, 1.2, 2.4 are drawn for the high- and low-excitation components, respectively ( and 0.6, respectively).

In the text

Fig. 7 CO velocity-integrated intensity diagrams predicted by non-stationary CJ-type shock models for various pre-shock densities: nH= 104 (upper left panel), 5 × 104 cm-3 (lower left panel), 105 (upper right panel), 5 × 105 cm-3 (lower right panel), and G. The observational data are drawn with filled (empty) squares for a filling factor of 1 (of 0.25). The shock velocity and age are given in each panel, as is the color code.

In the text

Fig. 8 Comparison of the best fitting CJ-type shock models and their predictions for the CO velocity-integrated flux distribution of the protostellar jet. The observational data are drawn with filled (empty) squares for a filling factor of 1 (of 0.25). Upper panel: a relative filling factor of 1 is assumed. Solutions with a pre-shock density nH = 5 × 104 (filled diamonds) and 105 cm-3 (filled circles), νs = 20−30 km s-1, and age in the range ~220–450 yr, with the color code indicated in the figure are found to provide an equally good fit. Lower panel: a relative filling factor of 0.25 is assumed. Solutions with a pre-shock density nH = 105 cm-3 (empty circles), νs = 25–30 km s-1, and age in the range ~595–740 yr, with the color code indicated in the figure, provide a good fit to the data.

In the text

	Fig. A.1 Map of the CO J = 2–1 emission integrated in channels of 3.25km from −13.5 to +113.1 km s-1, as observed with the IRAM 30 m telescope. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (10.7″) is indicated by the white circle.
In the text

	Fig. A.2 Map of the CO J = 2–1 emission integrated in channels of 3.25km from −137.0 to −10.3 km s-1, as observed with the IRAM 30 m telescope. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (10.7″) is indicated by the white circle.
In the text

	Fig. A.3 Map of the CO J = 3–2 emission integrated in channels of 3.25km from −13.5 to +113.1 km s-1, as observed with the JCMT. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (15″) is indicated by the white circle.
In the text

	Fig. A.4 Map of the CO J = 3–2 emission integrated in channels of 3.25km from −137 to −10.3 km s-1, as observed with the JCMT. The position of the protostar is marked by a red star. The positions targeted by GREAT/SOFIA and HIFI/Herschel are marked by blue stars. The size of the telescope HPBW (15″) is indicated by the white circle.
In the text

	Fig. B.1 χ2 distribution of LVG slab models for different velocity intervals in the low-excitation gas of the outflow cavity. The considered velocity intervals are indicated in each panel. The minimum value and the contour levels are indicated in each panel.
In the text

	Fig. B.2 χ2 distribution of LVG slab models for different velocity intervals in the high-excitation gas of the outflow cavity. The considered velocity intervals are indicated in each panel. The minimum value and the contour levels are indicated in each panel.
In the text

	Fig. B.3 Best fitting solutions of the CO spectral line flux distribution as a function of the rotational number of the transition from our LVG analysis with the MADEX radiative transfer code. The fits to the low- and high-excitation components are shown by black solid and dashed lines, respectively. The fit to the total flux is shown in red.
In the text