Detection of Gamma-Rays from the Protostellar Jet in the HH 80–81 System

Radio jets/outflows have been observed in protostars, which are possibly driven by the accretion in the formation of the stars (e.g., Rodríguez et al. 1990; Martí et al. 1993; Anglada et al. 1996; Rodriguez 1996; Garay et al. 2003; Rodríguez et al. 2005). In contrast to the relativistic jets in active galactic nuclei (AGNs), the protostellar jets move at much smaller velocities which are typically from 100 to 1000 km s⁻¹ (see Anglada et al. 2018, for a review). The radio emission is usually dominated by the free–free emission from the thermal motions of electrons, which is characterized by a positive spectral index and no linear polarization (e.g., Anglada et al. 2018).

Among the protostellar jets, the HH 80–81 jet (located at a distance of 1.7 kpc) is an intriguing object (e.g., Carrasco-González et al. 2010). The central source in the HH 80–81 system is IRAS 18162-2048, which is identified as a massive B-type protostar (e.g., Carrasco-González et al. 2012). The entire radio jet system extends up to ∼5 pc (e.g., Martí et al. 1993; Carrasco-González et al. 2010; Rodríguez-Kamenetzky et al. 2017). In the central region of the jet, the spectral index of the radio emission is positive, which suggests a dominating thermal emission in the region (e.g., Rodriguez 1996; Rodríguez-Kamenetzky et al. 2017). However, in some knots as well as in the lobes, the radio emissions show negative spectral indices, which suggests an additional non-thermal component in these regions (e.g., Martí et al. 1993; Carrasco-González et al. 2010; Rodríguez-Kamenetzky et al. 2017). In particular, the radio emission from the knots located ∼0.5 pc from the central source is found to be linearly polarized (Carrasco-González et al. 2010). This clearly confirms the non-thermal origin of the radio emission (Carrasco-González et al. 2010). The non-thermal radio emission is believed to be the synchrotron radiation of relativistic electrons in magnetic field.

Particle acceleration and γ-ray production in the protostellar jet have been studied (e.g., Araudo et al. 2007; Bosch-Ramon et al. 2010; Munar-Adrover et al. 2013; Rodríguez-Kamenetzky et al. 2016, 2017, 2019). The relativistic electrons in the jet would inverse-Compton (IC) scatter optical-ultraviolet photons to GeV γ-ray energies. In a dense material environment, γ-rays could be produced through the relativistic Bremsstrahlung process. If the protons in the jet are also accelerated accompany with the acceleration of the electrons, the inelastic proton-proton (pp) interaction also produce γ-rays. However, the detection of γ-rays from protostellar jet is still lacking. Motivated by the above arguments, we analyze the Fermi-LAT data in the direction of IRAS 18162-2048 to search for γ-rays in the HH 80–81 system.

In this work, we select 10 yr data observed by Fermi-LAT (from 2008 August 4 to 2018 August 4), covering the energies from 100 MeV to 300 GeV. The Pass 8 SOURCE class events within 14° × 14° region of interest centered at the position of IRAS 18162-2048 are used (Atwood et al. 2013). We use instrument response function P8R3_SOURCE_V2. The Fermi-LAT fourth source catalog (4FGL; The Fermi-LAT collaboration 2019) based on the 8 yr data is used to construct the background model. Galactic and extragalactic diffuse components are modeled by gll_iem_v07.fits and iso_P8R3_SOURCE_V2_v1.txt, respectively. We employ the Fermi-tools to perform the analysis. To reduce the contamination from the Earth limb, we exclude the events within zenith angles of >90°. Standard binned likelihood analysis is performed to fit the free parameters in the model, and the obtained best-fit model is used in next procedure.

2.1. TS maps and SEDs

We use gttsmap to create 5° × 5° Test-Statistic (TS) map. It is obtained by moving a putative point source through a grid of locations on the sky and maximizing the likelihood function (−logL) at each grid point.

In Figure 1, one can see a significant γ-ray excess in the direction of IRAS 18162-2048. We add a point-like source at the position of IRAS 18162-2048 to describe this excess. The spectrum of the source is assumed to be a power-law form,

$$\begin{eqnarray}&&\displaystyle \frac{{dN}}{{dE}}{\left(E\right)={F}_{0}\left(\displaystyle \frac{E}{0.1\ \,\mathrm{GeV}}\right)}^{-{{\rm{\Gamma }}}_{\gamma }},\end{eqnarray} \tag{ 1 }$$

where the normalization F₀ and photon index Γ_γ are free parameters. Under such assumption, the tool gtlike gives TS value of 101 for the emission above 100 MeV (panel a). The TS value is reduced to 25 for the emission above 200 MeV (panel b), while above 300 MeV it is only 5 (panel c). From the residual map (panel d), one can find that the γ-ray emission can be well described by a point-like source with a power-law spectrum. The γ-ray excess has Γ_γ = 3.53 ± 0.11 and the flux above 100 MeV F_γ = (5.2 ± 0.4) × 10⁻⁸ photons cm⁻² s⁻¹.

Zoom In Zoom Out Reset image size

We then construct spectral energy distribution (SED). Energies are divided into logarithmic-equivalent bins, and in each bin we perform likelihood analysis to calculate the flux. In this step, only the prefactors of all sources are left free (i.e., fix the spectral shapes to those derived before). The SED of the γ-ray emission is shown in Figure 2 (black data). Again, one can see that the majority of emission is below 1 GeV.

Zoom In Zoom Out Reset image size

To investigate the effect of the fluctuation of the diffuse Galactic background (the so-called systematic uncertainty) on the γ-ray signal, we change the normalization of the Galactic background component to do the analysis. The TS value above 100 MeV changes to 76 when the diffuse Galactic background is enhanced by 3%, and it is 144 when the background is reduced by 3%. The derived systematic uncertainty on F_γ is ${\sigma }_{\mathrm{syst}}{=}_{-0.77}^{+0.45}\times {10}^{-8}\ \mathrm{photons}\,{\mathrm{cm}}^{-2}\ {{\rm{s}}}^{-1}$. The corresponding SEDs are shown in Figure 2 with red and blue circles, respectively. One can find that the effect of the systematic uncertainty on the SED is negligible.

To get accurate position of this γ-ray emission, we run gtfindsrc to optimize the source location. This is performed using the photons above 200 MeV, with the consideration of a narrower point-spread function and lower computation complexity. The best-fit position is (R.A. = 2748, Decl. = −208; J2000) with the 1σ error circle radius of 028. IRAS 18162-2048 is well located in this region.

2.2. Extension

2.3. Variability

We obtain the flux variation during the observation. Light curve is one-year binned (see the top panel in Figure 3). Since the γ-ray emission is weak, in each time bin we fix the spectral shapes of background sources (i.e., only free their prefactors). The χ² under constant hypothesis is 14.02 with 8° of freedom (d.o.f.). Its corresponding p-value is 0.08 (∼1.4σ significance), suggesting no significant variability. The variation of Γ_γ is fitted well by a constant of 3.3 ± 0.2 (see the bottom panel in Figure 3).

Zoom In Zoom Out Reset image size

For the γ-ray excess region, we exclude known possible γ-ray sources within the error circle of 028 (this circle is obtained from previous gtfindsrc performing), by cross-checking known database.

Totally 315 objects are retrieved in this circle from SIMBAD database. ⁵ We carefully check these objects and find that most of them are alternative names of the HH 80–81 system. Among the rest of the objects (most are stars), no potential γ-ray emitter is found. The only γ-ray candidate 3FGL J1819.5-2045c ($4.5^{\prime}$ away) has been removed in the Fermi-LAT fourth source catalog. In 4FGL, the nearest known γ-ray source, 4FGL J1818.1-2000, is about 08 away, which has little impact on the γ-ray excess (see Figure 1).

The 5th Roma blazar catalog (Massaro et al. 2009, 2015) contains 3561 sources, and the nearest one from IRAS 18162-2048, 5BZQJ1833-2103, has an offset of 34. By checking the Australia Telescope National Facility Pulsar Catalog (Manchester et al. 2005), no pulsar is found in this circle. The Molonglo Reference Catalog which is one of the largest homogeneous catalogs of radio sources (Large et al. 1981, 1991) does not cover regions within 3° of the Galactic equator.

Relativistic particles in a protostellar jet can produce γ-rays through a variety of processes, including IC scattering, relativistic Bremsstrahlung, and pp interaction (e.g., Araudo et al. 2007; Bosch-Ramon et al. 2010). The present observations cannot determine the radiative mechanism for the γ-rays. Nevertheless, the γ-ray spectrum that we obtained can put constraints on the characteristics of the relativistic particles in the jet.

The characteristic cooling time of pp interaction in the hydrogen medium with number density n₀ is written as (e.g., Aharonian 2004)

$$\begin{eqnarray}&&{t}_{{pp}}={({n}_{0}{\sigma }_{{pp}}{fc})}^{-1}\approx 5\times {10}^{7}{({n}_{0}/1\,{\mathrm{cm}}^{-3})}^{-1}\ \mathrm{yr}.\end{eqnarray} \tag{ 2 }$$

Here, an average cross-section at high energies of about 40 mb is used, and the coefficient of inelasticity f = 0.5 is adopted (assuming that on average the proton loses about half of its energy per interaction). For a typical n₀ ∼ 1000 cm⁻³ in HH 80–81 (e.g., Anglada et al. 2018), we derive t_pp ∼ 50,000 yr. The ratio of the mean energy of the produced γ-ray to the energy of the incident proton is ∼0.05 (Kelner et al. 2006). Therefore, the energy of the proton that produces the observed 1 GeV photons through pp interaction is ∼20 GeV. The spectral index of the protons distribution s_p is same as the γ-ray photon index (Kelner et al. 2006), i.e., s_p ≈ 3.5.

The cooling time of electrons due to the Bremsstrahlung losses is (e.g., Aharonian 2004)

$$\begin{eqnarray}&&{t}_{\mathrm{br}}\approx 4\times {10}^{7}{({n}_{0}/1\,{\mathrm{cm}}^{-3})}^{-1}\ \mathrm{yr}.\end{eqnarray} \tag{ 3 }$$

With n₀ ∼ 1000 cm⁻³, we get t_br ∼ 40,000 yr. The energy of the electron that produces the observed 1 GeV photons through Bremsstrahlung is 1 GeV (Blumenthal & Gould 1970). The spectral index of the electron distribution is s_e ∼ Γ_γ ≈ 3.5. In this case, the maximum energy of synchrotron photon is ${E}_{\mathrm{syn},\ \max }{(\mathrm{eV})\approx 2\,\times \,{10}^{-8}(1\ \mathrm{GeV}/{m}_{e}{c}^{2})}^{2}B$. Considering a magnetic field strength of 0.1 mG (Carrasco-González et al. 2010), we have ${E}_{\mathrm{syn},\ \max }(\mathrm{eV})\sim {10}^{-5}$. This energy range is covered by the current radio observations, and the photon index of the observed radio spectrum is 1.3 (Martí et al. 1993). From the relativistic electrons distribution, the photon index of the synchrotron emission is expected as (s_e + 1)/2 ≈ 2.3 (Blumenthal & Gould 1970). It is noted that this predicted synchrotron spectrum at ∼0.01 meV is inconsistent with the radio observations.

The cooling time of relativistic electrons due to the IC scattering is (e.g., Aharonian 2004)

$$\begin{eqnarray}&&{t}_{\mathrm{IC}}=\displaystyle \frac{3{m}_{e}{c}^{2}}{4{\sigma }_{{\rm{T}}}c\gamma {U}_{{\rm{s}}}}\approx (1/\gamma )({U}_{{\rm{s}}}/\ \mathrm{erg}\,{\mathrm{cm}}^{-3})\ \mathrm{yr},\end{eqnarray} \tag{ 4 }$$

where γ is the Lorentz factor of the relativistic electrons and U_s is the energy density of seed photon field. The energy of the scattered photons E₁ is written as E₁ ≈ γ² E_s (Blumenthal & Gould 1970), where γ is the Lorentz factor of the relativistic electrons and E_s is the energy of the seed photons. Protostars are detected as strong infrared sources. Assuming E_s ∼ 0.01 eV (corresponding to the temperature T ∼ 100 K) (Bosch-Ramon et al. 2010), we can derive γ ≈ 3 × 10⁵ with E₁ = 1 GeV. Considering U_s ∼ 10⁻¹² erg cm⁻³ (Bosch-Ramon et al. 2010), we have t_IC ∼ 2 × 10⁶ yr for the electrons that produce 1 GeV photons. The energy of synchrotron photon is calculated by E_syn(eV) ≈ 2 × 10⁻⁸ γ² B, where B is strength of magnetic field. With B = 0.1 mG, the maximum energy of the synchrotron photons is ≈0.2 eV. There is no observed data at this energy range. In this scenario, s_e is predicted as s_e = 2Γ_γ − 1 ≈ 6. This very steep spectrum indicates that there is a cut-off at γ_c (γ_c < 3 × 10⁵) in the electron distribution.

As mentioned above, the observed γ-rays can be produced through IC scattering of relativistic electrons or inelastic pp interaction. It should be noted that the jet phase timescale in massive protostars is roughly 40,000 yr (Guzmán et al. 2012). This lifetime is comparable to t_pp , but significantly shorter than the characteristic t_IC for the electrons that produce 1 GeV photons.

The kinetic luminosity of the jet can be calculated by

$$\begin{eqnarray}&&{L}_{j}=\displaystyle \frac{1}{2}{\dot{M}}_{j}{v}_{j}^{2},\end{eqnarray} \tag{ 5 }$$

where ${\dot{M}}_{j}$ is mass-loss rate in the jet and ${v}_{j}^{2}$ is jet velocity. Taking ${\dot{M}}_{j}={10}^{-6}\ {M}_{\odot }\ {\mathrm{yr}}^{-1}$ and v_j = 1000 km s⁻¹ (e.g., Anglada et al. 2018), we have L_j ≈ 3 × 10³⁵ erg s⁻¹. The γ-ray luminosity (>100 MeV) is L_γ ≈ 5 × 10³³ erg s⁻¹ which is 0.017% of L_j . The energy budget is reasonable.

The detection of the γ-rays from HH 80–81 not only confirms again that there is a population of relativistic particles in the protostellar jet, but also makes the protostellar jet an additional high-energy astrophysics laboratory. Evidence for strong shocks has been found in the jet of HH 80–81 (e.g., Rodríguez-Kamenetzky et al. 2019), which are produced by the interaction of the jet with the environment. Diffusive Shock Acceleration (DSA) is therefore considered as the mechanism that accelerates particles to high energies (e.g., Padovani et al. 2015; Araudo et al. 2021). However, details on the acceleration mechanism are still unclear. The jet properties such as jet ionization fraction and jet temperature which can be probed through detection of molecular emission lines are crucial for investigating the acceleration mechanism (e.g., Araudo et al. 2021). Moreover, the detection of non-thermal emissions from radio to X-ray energies is helpful to determine the magnetic field in the jet, which is crucial for investigating particle acceleration mechanism and radiative processes. Follow-up intensive multi-band observations are necessary to reveal the nature of the HH 80–81 jet.

The clear detection of linearly polarized radio emission from the HH 80–81 jet (Carrasco-González et al. 2010) confirms the existence of relativistic particles in the protostellar jet. γ-rays from the HH 80–81 jet are therefore expected. We search for γ-rays from the HH 80–81 system using 10 yr Fermi-LAT data. A significant γ-ray excess is found toward the direction of HH 80–81 system. After excluding other possible γ-ray candidates, we suggest that the stable and point-like γ-ray emission is contributed by the HH 80–81 system, which is likely the non-thermal γ-rays from the HH 80–81 jet. The spectrum of the γ-ray emission extends to 1 GeV with a photon index of 3.5. Moreover, from the aspects of emission mechanism and energy budget, the observed γ-ray emission can be produced by the HH 80–81 jet. Our result suggests the protostellar jet as an efficient particle accelerator, which is proposed in theoretical studies (e.g., Padovani et al. 2015, 2016).

We thank the anonymous referee for his/her comments that helped us to improve the quality of our paper. We thank Dr. Ruizhi Yang for helpful discussions. We acknowledge financial supports from the National Natural Science Foundation of China (NSFC-11803081, NSFC-U1931114, NSFC-U2031205 and NSFC-12163006) and the joint foundation of Department of Science and Technology of Yunnan Province and Yunnan University (2018FY001(-003)). The work of D.H. Yan is also supported by the CAS Youth Innovation Promotion Association and Basic research Program of Yunnan Province (202001AW070013).