Detection of extended gamma-ray emission in the vicinity of Cl Danks 1 and 2

ABSTRACT

We report the detection of high-energy gamma-ray emission towards the G305 star-forming region. Using almost 15 yr of observation data from Fermi Large Area Telescope, we detected an extended gamma-ray source in this region with a significance of |$\sim 13 \sigma$|⁠. The gamma-ray radiation reveals a clear pion-bump feature and can be fitted with the power-law parent proton spectrum with an index of |$-2.5$|⁠. The total cosmic ray (CR) proton energy in the gamma-ray production region is estimated to be of the order of |$10^{49}\ \rm erg$|⁠. We further derived the CR radial distribution from both the gamma-ray emission and gas distribution and found it roughly obeys the |$1/r$| type profile, though a constant profile is not ruled out. This is consistent with other similar systems and expected from the continuous injection of CRs by the central powerful young massive star cluster Danks 1 or Danks 2 in this region. Together with former detections of similar gamma-ray structures, such as Cygnus cocoon, Westerlund 1, Westerlund 2, NGC 3603, and W40, the detection supports the hypothesis that young massive star clusters are CR accelerators.

1 INTRODUCTION

The origin of cosmic rays (CRs) remains an open question in astrophysics. The recent advances in gamma-ray astronomy indicate that the young massive star clusters (YMCs) can be an alternative type CR sources (Aharonian, Yang & de Oña Wilhelmi 2019) in addition to supernova remnants (SNRs). Furthermore, several independent observations also hint at such a hypothesis. For example, recent findings by Binns et al. (2016) suggest that a significant portion of CRs may be accelerated in young OB star clusters and associated superbubbles, as indicated by measurements of |${ }^{60} \mathrm{Fe}$| abundance in CRs. Additionally, studies on the Galactic diffuse gamma-ray emission reveal that the radial distribution of CRs aligns more closely with OB stars rather than SNRs, as reported by Acero et al. (2016) and Yang, Aharonian & Evoli (2016). Superbubbles, fueled by supernova explosions and collective stellar winds, are capable of providing the necessary kinetic energy to generate the observed flux of CRs locally, as proposed by Parizot et al. (2004). These structures are expected to emit gamma-rays due to the interaction of freshly accelerated CRs with the surrounding gas. In recent days, YMCs have become popular candidates for CR accelerators in several studies. The Cygnus cocoon is a typical example of such sources, despite the Cygnus OB system not being the most potent young star cluster in our Galaxy. Ackermann et al. (2011) have studied the Fermi-LAT data of the Cygnus cocoon and provided an illustration for investigating the early stages of CR youth within a superbubble setting before they integrate into the older Galactic population. Other YMC examples include Westerlund 2 (Yang, de Oña Wilhelmi & Aharonian 2018), W40 (Sun et al. 2020b), W43 (Yang & Wang 2020), Carina (Ge et al. 2022), G25 (Sun, Yang & Wang 2020a), NGC 3603 (Yang & Aharonian 2017), NGC 6618 (Liu, Yang & Chen 2022), and NGC 2244 (Liu, Liu & Yang 2023). Characteristics of CRs in such regions are well studied and described, making YMCs important sites for studying CRs.

The G305 complex (G305.4+0.1) is identified as one of the most massive star-forming structures in the Galaxy. It is located in the Scutum-Crux arm (⁠|$l=305.4^\circ$|⁠, |$b=+0.1^\circ$|⁠; Clark & Porter 2004) at an estimated distance of |$\sim 4~\rm kpc$|⁠, and was delineated by both mid- and far-infrared, submillimetre, and radio emission (Davies et al. 2012). The presence of methanol masers and ultracompact H ii regions (G305 H ii complex) indicates ongoing massive star formation in the G305 complex. It has the form of a large trilobed cavity with a maximum extent of |$\sim 34~\rm pc$|⁠, delineated by both mid- and far-infrared, submillimetre, and radio emission. The existence of embedded massive young stellar objects and a substantial reservoir of cold molecular gas available for future activity is observed in the G305 complex. Similar properties are also observed in the star-forming complex associated with NGC 3603 and 30 Dor.

The young clusters Danks 1 and 2 are situated at the core of the star-forming complex. The age of Danks 1 is |$1.5_{-0.5}^{+1.5}~\mathrm{Myr}$|⁠, and the age of Danks 2 is |$3_{-1}^{+3}~\mathrm{Myr}$| (Davies et al. 2012). The complex’s overall morphology strongly suggests several epochs of sequential star formation, which have been initiated and sustained by the influence of these two central clusters (Davies et al. 2012). The G305 complex’s morphology strongly suggests multiple epochs of sequential star formation, a process initiated and sustained by the influence of the two central clusters. Photometric investigations of these clusters were undertaken by Bica et al. (2004) and Baume, Carraro & Momany (2009), albeit somewhat impeded by significant visual extinction and densely populated fields. Regarding the distance to Danks 1 and Danks 2, we assume they are located at a distance of |$4~\rm kpc$|⁠, the same distance as the G305 complex (Davies et al. 2012).

In this paper, we analysed the 15 yr of Fermi-LAT data in this region and found hard gamma-ray emission associated with radius from the centre YMCs Danks 1 and Danks 2. We also derived the radial distribution of CR energy density, which is consistent with the former results. This paper is organized as follows. In Section 2, we studied the distribution of three phases of gas content around Danks 1 and Danks 2. In Section 3, we analysed Fermi-LAT data and studied the spatial and spectral distribution of gamma-ray emissions. In Section 4, we derived the CR spectrum utilizing analysis results of Section 3 and Section 2, and investigated the distribution of CR energy density over radius. Finally, in Section 5, we discussed our results along with former works and drew our conclusions.

2 GAS CONTENT AROUND DANKS 1 AND DANKS 2

We investigated three different gas phases, i.e. the H|$_{2}$|⁠, the neutral atomic hydrogen (H i), and the H ii, in the vicinity of Danks 1 and Danks 2. We followed the method in section 3 of Ge et al. (2022) to derive the gas column density map of these three gas phases scaled to H column density.

First, we used the CO composite survey data (Dame, Hartmann & Thaddeus 2001) to trace the |$\rm H_{2}$|⁠. We took |$N({\rm H_{2}}) = X_{\rm CO} \times W_{\rm CO}$| (Lebrun et al. 1983), where |$X_{\rm CO}$| is the |$\rm H_{2}/CO$| conversion factor that is chosen to be |$\rm 2.0 \times 10^{20}\ cm^{-2}\ K^{-1}\ km^{-1}\ s$| as suggested by Dame et al. (2001) and Bolatto, Wolfire & Leroy (2013).

Then for the H i density, we used the data cube of the H i|$\rm {4\pi }$| survey (HI4PI), which is a 21-cm all-sky data base of Galactic H i (HI4PI Collaboration 2016). Under an optically thin assumption, the H i column density can be estimated as:

where |$T_{\rm B}$| is the brightness temperature of the H i emission. However, the optical thin approximation can no longer be used for Danks 1 and 2, as their positions are close to the Galactic disc. So we estimated the H i column density using the equation,

where |$T_{\rm bg} \approx 2.66\ \rm K$| is the brightness temperature of the cosmic microwave background radiation at 21 cm. In addition, in the case when |$T_{\rm B} \gt T_{\rm s}-5\ \rm K$|⁠, we truncate |$T_{\rm B}$| to |$T_{\rm s}-5\ \rm K$|⁠, where |$T_{\rm s}$| is chosen to be 150 K.

The radial velocity range of G305 star-forming complex was estimated to be |$[-50.6, -24.9]~\rm km~s^{-1}$| by Hindson et al. (2013). The total gas mass and mean column density of the H atom in these two phases are shown in Table 1. The integrated maps of the H atom column density in these two phases are shown in the left and middle panels of Fig. 1.

G305 region is one of the largest ionized gas complexes in our Galaxy, thus the ionized gas can also make a significant contribution to the total gas mass. For the H ii column density, we used the Planck free–free map (Planck Collaboration XIII 2016). First, we transformed the emission measure into free–free intensity using the conversion factor in table 1 of Finkbeiner (2003). Then, we calculate the H ii column density from the intensity (⁠|$I_{\nu }$|⁠) of free–free emission by using equation (5) of Sodroski et al. (1997),

in which the frequency |$\nu = \rm 353\ GHz$| and the electron temperature |$T_{e} =\rm 8000\ K$|⁠. The H ii column density is inversely proportional to the effective density of electrons |$n_{\rm e}$|⁠. Thus, we assumed an effective density |$10\ \rm cm^{-3}$|⁠, which is the value suggested in Sodroski et al. (1997) for the region inside the solar circle. This could introduce considerable uncertainty to the H ii column density. The H ii accounts for approximately 30 per cent of the total H column density in our work. Consequently, the uncertainty in the H ii column density can lead to a 10 per cent to 20 per cent uncertainty in the final estimation of the total H column density. The total gas mass and mean column density of the H atom in H ii phase are shown in Table 1. The map of H atom column density in H ii phase is shown in the right panel of Fig. 1.

3 FERMI-LAT DATA ANALYSIS

To study the gamma-ray emission in the vicinity of Danks 1 and 2, we used the latest Fermi-LAT Pass 8 data from 2008 August 4 (MET 239557417) until 2023 May 22 (MET 706481557), and used the Fermitools from Conda distribution¹ for the data analysis. Considering the extended feature of the gamma-ray emission in the vicinity of Danks 1 and 2, we chose a 20|$^\circ$||$\times$| 20|$^\circ$| square region centred at the position of Danks 1 and 2 (R.A.  = 198.222|${^\circ }$|⁠, Dec.  = −62.681|${^\circ }$|⁠) as the region of interest (ROI). This choice also accounts for the location of Danks 1 and 2 near the Galactic disc, where accurate modelling of the diffuse emission is crucial. We utilized gtselect to choose ‘source’ class events (evtype = 3 and evclass  = 128) with zenith angles below 90|$^{\circ }$| to eliminate background contamination from the Earth’s limb. The good time intervals were determined by applying the expression |$\rm (DATA\_QUAL \gt 0) \&\& (LAT\_CONFIG == 1)$|⁠. The analysis of the SOURCE events was conducted using the instrument response functions P8R3_SOURCE_V3.

A standard binned analysis was carried out following the official tutorial of binned likelihood analysis.² The source models were generated using make4FGLxml.py,³ which included sources from the LAT 12-Year Source Catalogue (4FGL-DR4; Ballet et al. 2023) within the ROI enlarged by 10|$^{\circ }$|⁠, the Galactic diffuse background emission (gll_iem_v07.fits), and the isotropic emission background (iso_P8R3_SOURCE_V3_v1.txt). The spectral parameters of sources within 10|${^\circ }$| from the centre and the normalization factor of the diffuse and isotropic background were set free.

In our ROI, there are several sources in the 4FGL catalogue. Their positions are shown in Fig. 2. During the likelihood analysis, we tested the results when the sources within 1 degree that were classified as ‘c’ sources were removed. Such classification indicates that they require caution in interpretation or analysis.

3.1 Spatial analysis

We first used the events from 1 to |$500~\rm GeV$| to study the spatial distribution of the gamma-ray emission in the ROI. The gamma-ray counts map in the |$3{^\circ }\times 3{^\circ }$| region of the vicinity of Danks 1 and 2 is shown in the left panel of Fig. 2. The clusters mentioned in Section 1 are marked in the counts map. There are three sources in 4FGL-DR4 in the ROI with corresponding identified sources (4FGL J1305.5-6241 and Kes 17, 4FGL J1315.9-6349c and PMN J1315-6354, HESS J1303-631), which are shown in blue marks. Other sources in 4FGL-DR4 are shown in green crosses.

With the sources in the 4FGL-DR4 catalogue, we conducted a binned likelihood analysis. We then added a point source (pt) at the position of Danks 1 and 2 (R.A.  = 198.222|${^\circ }$|⁠, Dec.  = −62.681|${^\circ }$|⁠) to account for the emission near Danks 1 and 2. We also tested a radial Gaussian template (disc) centred at the same position. We calculated the log-likelihood −|$\log ({\cal L})$| assuming different sigma of the Gaussian disc and found that the best fit is achieved when sigma equals 0.5|$^{\circ }$| (−|$\log ({\cal L})$| = 21447). The significance of this added source, corresponding to 2 degrees of freedom for a |$\chi ^2$| distribution, was calculated to be |$\sim 12\sigma$|⁠. We also tested another spacial template generated using the total gas column density map derived in Section 2. To build this template, we normalized the three different gas phases, i.e. the H|$_{2}$|⁠, the neutral atomic hydrogen (H i), and the H ii, to the column density of H atom (N_H), and summed them together. The likelihood fitting result using this gas template (−|$\log ({\cal L})$|21446) is similar to the Gaussian disc. Furthermore, we utilized the Akaike Information Criterion (AIC) to determine which models best suited the data. The AIC is computed by subtracting twice the logarithm of the likelihood from twice the number of free parameters (k) in the model, expressed as |$-2 \log ({\cal L})+2 \mathrm{k}$|⁠. The comparison of fitting results of different models can be seen in Table 2. So in the later analysis, we used the radial Gaussian disc as the spacial template for the gamma-ray excess around Danks 1 and 2. Then, we generated the model maps with this source model using gtmodel, and obtained the residual counts map by subtracting the contribution of the diffuse background and sources other than the Gaussian disc, which is shown in the middle panel of Fig. 2. We also calculated the test statistic (TS) map of this model using gttsmap. The TS value is defined as |$2\left(\log \mathcal {L}/\mathcal {L}_0\right)$|⁠, where |$\mathcal {L}_0$| is the likelihood of null hypothesis and |$\mathcal {L}$| is the likelihood with a test source included at a given position. The TS map is shown in the right panel of Fig. 2.

3.2 Spectral analysis

To study the gamma-ray spectrum from the vicinity of Danks 1 and 2, we divided the data from |$200\rm ~MeV$| to |$500\rm ~GeV$| into eight energy bins evenly distributed in logarithmic space. We tested several methods for selecting energy bins by dividing the data into 8, 9, and 10 energy bins. We obtained similar gamma-ray spectra in most cases. However, for certain energy bin configurations, the likelihood fitting for the lowest energy bin failed. Ultimately, we selected a method that ensures successful fitting in most energy bins. For the data in each energy bin, we conducted likelihood analysis with the source model in Section 3.1 and derived the total flux of the region near Danks 1 and 2. Systematic uncertainties of the flux for each energy bin were also estimated using the method described on the Fermi official website.⁴ We applied photon weights during the likelihood analysis for the energy bins below and around 1 GeV to avoid systematics induced by diffuse emission. The derived gamma-ray spectrum is shown in Fig. 3. The total gamma-ray luminosity was estimated to be |$\sim 1.0\times 10^{35}~\rm erg~s^{-1}$|⁠. Besides, we calculated the 95 per cent flux upper limit for energy bins with the source test statistic (TS) value less than 4, those are shown in inverted triangles. We note that the energy flux of the first energy bin shows a clear hint of a low energy break, which is consistent with the pion-bump feature in the hadronic scenario of gamma-ray productions.

4 CR CONTENT IN THE VICINITY OF DANKS 1 AND 2

Due to the clear pion-bump feature, we assumed that the gamma-ray emission from Danks 1 and 2 is mainly generated from the proton–proton interaction between CR protons and the gases. Under such assumptions, we studied the CR content around Danks 1 and 2 utilizing gamma-ray data from Section 3 and gas distributions derived from Section 2.

4.1 CR spectral analysis

We assumed the CR protons follow the momentum distribution function: |$F=f_ip^{\Gamma }$|⁠, where F is the flux density of protons in the unit of |$\mathrm{cm}^{-2}\, \mathrm{s}^{-1}\, \mathrm{GeV}^{-1}$|⁠, p is the proton momentum in the unit of GeV c⁻¹, and |$\Gamma$| stands for the spectral index. Note that during the likelihood analysis of Fermi-LAT data, the emission of the Galactic diffuse background emission has been taken into account (see Section 3). So the following results of the CR density in the vicinity of Danks 1 and 2 is actually an excess from additional CR sources.

Under this assumption, we used the gamma-ray production cross-section from Kafexhiu et al. (2014) to perform a likelihood fitting to the gamma-ray data. We also calculated the gas mass in the gas template we used in Section 3.1 so that we can derive the absolute CR fluxes. We found that the best-fitting CR proton spectral index |$\Gamma$| is 2.36|$\pm$|0.09. The corresponding gamma-ray spectrum is shown in the dotted line in Fig. 3.

4.2 The radial distribution of CR protons

Using both the gas and gamma-ray spatial distributions, we can derive the spatial distribution of CRs, which contain unique information on the injection and propagation of CRs. If the CRs are injected continuously from the central source as in some other YMC systems (Aharonian et al. 2019), then in regions near the source, the energy distribution of CR protons/electrons |$F(R) \sim \frac{1}{R D}$| (Yang & Liu 2022), in which R is the radial distance and D is the diffusion coefficient of CRs. Consequently, CR energy density |$\Omega _\mathrm{CR}(R)$| would also be proportional to 1/R, i.e.

where |$R_0$| is the size of the emission region taken to be |$56\rm ~pc(\sim$|0.8|$^\circ$|⁠), and |$\omega _0$| is the CR energy density at distance |$R_0$|⁠. Due to the projection effect, the 2D radial distribution function of CR energy density can be calculated as follows:

in which r is the projected radial distance. We can then derive:

To test whether the CR energy density in the vicinity of Danks 1 and 2 follows the distribution described above, we divided the emissive region into four regions: a central disc with a radius of |$14\rm ~pc$| (approximately 0.2|$^\circ$|⁠), and three concentric rings extending from 14 to |$28\rm ~pc$| (approximately 0.2|$^\circ$| to 0.4|$^\circ$|⁠), 28 to |$42\rm ~pc$| (approximately 0.4|$^\circ$| to 0.6|$^\circ$|⁠), and 42 to |$56\rm ~pc$| (approximately 0.6|$^\circ$| to 0.8|$^\circ$|⁠), respectively. Considering the angular resolution of Fermi-LAT, to obtain sufficient statistical power, choosing 0.2|$^\circ$| as the ring width is preferable. If wider rings are selected, there would be only three or fewer rings, making it difficult to analyse the radial distribution of cosmic ray energy density. We also tested adding one more ring and found that the ring extending 56 to |$70\rm ~pc$| (approximately 0.8|$^\circ$| to 1.0|$^\circ$|⁠) provides insufficient flux to be modelled. We replaced the radial Gaussian disc with these four flat spacial templates and applied the same spatial analysis with that of Section 3.1 to derive the gamma-ray flux of each region. Next, assuming the spectral shape of CRs is uniform in the whole region, we calculated the CR energy density of each region. Applying maximum-likelihood estimation, we then derived that |$\omega _0=0.1~\rm eV~cm^{-3}$|⁠. We calculated |$\chi ^2$| of each distribution and obtained |$\chi ^2/\rm d.o.f=4.07/3$| for the uniform profile and |$\chi ^2/\rm d.o.f=1.30/3$| for the |$1/r$| profile. From the statistics, the |$1/R$| distribution is favoured marginally. The comparison of the fitting results of two assumptions is shown in Fig. 4. We also calculated the local CR energy density by integrating the local CR proton flux measured by AMS-02 (Aguilar et al. 2015), which is also shown in Fig. 4 as a comparison. In addition, we integrated equation (4) and derived the total CR energy |$E_\mathrm{CR}$|⁠:

where |$R_{\rm dif}$| is the diffusion length of the CRs injected by the central accelerator. It should be noted that |$R_{\rm dif}$| can be significantly larger than |$R_{0}$|⁠, which represents the radius of the gamma-ray emission region. Since the size of the measured gamma-ray emission is limited by both the gas distribution and the gamma-ray instrument sensitivity, the CRs can indeed occupy a much larger volume. Assuming the age of Danks 1 of |$T\sim 1.5~\rm Myr$| and a diffusion coefficient of |$D \sim 10^{28}~\rm cm^2\,s^{-1}$|⁠, |$R_{\rm dif}$| can be estimated as |$R_{\rm dif} \sim 2 \sqrt{DT} \sim 400~\rm pc$|⁠, which is much larger than the gamma-ray emission region detected in GeV band.

4.3 Leptonic scenario

From the gamma-ray spectrum, we cannot formally rule out a leptonic scenario. In this section, we considered the leptonic scenario, in which the observed GeV gamma-ray emission is generated through two potential mechanisms: inverse-Compton (IC) scattering of relativistic electrons off low-energy seed photons, or non-thermal bremsstrahlung radiation produced by the relativistic electrons interacting with the surrounding matter in the regions around Danks 1 and 2. In the calculations of the photon field for IC processes, we accounted for the cosmic microwave background (CMB) radiation using naima package (Zabalza 2015) based on the model by Khangulyan, Aharonian & Kelner (2014). In this context, Danks 1 and 2 are situated within an H ii region, where ionizing massive stars significantly enhance the optical and UV fields, leading to increased IC emissions (Liu & Yang 2022). The target particle density is assumed to be |$113~\rm cm^{-3}$| for relativistic bremsstrahlung. This is estimated from the gas map in Section 2. It is important to note that variations in ambient matter density will merely scale the contributions to the bremsstrahlung spectrum. For the relativistic bremsstrahlung spectrum, we used the parametrization in Baring et al. (1999). To fit the lower energy break in the gamma-ray spectrum, it is necessary to establish a corresponding break in the spectrum of parent electrons. Consequently, we adopted a broken power-law distribution for the relativistic electrons,

where A, |$\alpha _1$|⁠, |$\alpha _2$|⁠, |$E_\mathrm{b}$| are treated as free parameters for the fitting. As is shown in Fig. 5, both the IC and bremsstrahlung model can fit the observable data.

For IC model, the derived parameters for the electrons are |$\alpha _1=1.3_{-1.4}^{+1.2}$|⁠, |$\alpha _2=4.0_{-0.4}^{+0.5}$|⁠, |$E_\mathrm{b}=359_{-102}^{+142}~\rm GeV$|⁠, and the total energy of the electrons above |$2\rm ~GeV$| is estimated to be around |$1.3\times 10^{48}\rm ~erg$|⁠. For bremsstrahlung, |$\alpha _1=1.1_{-1.0}^{+0.8}$|⁠, |$\alpha _2=2.6_{-0.2}^{+0.3}$|⁠, |$E_\mathrm{b}=1.9_{-0.9}^{+2.2}~\rm GeV$|⁠, and the total energy of the electrons above |$2\rm ~GeV$| is estimated to be around |$1.3\times 10^{46}\rm ~erg$|⁠. We cannot completely exclude a leptonic origin for this source. Bremsstrahlung dominates the radiation mechanism in such a dense region. The good spatial coincidence of the gamma-ray emission with the H ii gas region also favours a hadronic or bremsstrahlung origin. However, in the bremsstrahlung scenario, the derived parent proton index undergoes a significant change from 1.1 to 2.6 around an energy of approximately |$2~\rm GeV$|⁠. This abrupt transition is highly unusual and cannot be reconciled with any established mechanisms.

5 DISCUSSION AND CONCLUSION

In this paper, by analysing almost 15 yr of Fermi-LAT data, we found extended gamma-ray emission in the vicinity of Danks 1 and 2, which can be modelled with a radial Gaussian disc. From the derived spectral energy distribution (SED) of the extended gamma-ray emission, we found a bump feature, which is a strong hint that the gamma-ray emissions are produced via the pion-decay process in the interaction of CR protons with ambient gas. A power-law distribution of CR protons with an index of |$\sim 2.5$| can typically reproduce the gamma-ray SED.

We found that the Gaussian model is as good as the gas distribution model in modelling the gamma-ray emission, which may indicate that the distribution of gamma-ray emission is influenced by both the radial distribution of CRs and the spatial distribution of gas. We obtained the gamma-ray fluxes at various radii by applying the spatial templates of one disc and three rings, then calculated the variation of CR energy density across different radii. We found that the CR energy density around Danks 1 and 2 decreases with radius and its radial distribution is generally consistent with |$1/r$| profile as other similar structures in Aharonian et al. (2019) considering the projection effect, though a constant profile is not completely ruled out. This also implies a continuous injection of CRs as in other YMC systems.

Given the extended emission and hard spectrum, the gamma-ray emission in this region is more likely produced by the CR injected by the YMC interacting with the ambient gas, which is similar to other gamma-ray emission YMCs (Aharonian et al. 2019). Danks 1 and Danks 2 can both be the CR accelerators. Indeed, both clusters are powerful enough to accelerate the CRs needed to explain the observed gamma-ray emissions. Danks 1 is associated with six WR stars, while Danks 2 is associated with one WR star (Celli et al. 2024). The winds of these WR stars alone can provide the kinetic power to accelerate the CRs, and dozens of OB stars in these two clusters can provide similar power. The total kinematic power can be estimated roughly as |$(3-5) \times 10^{38} \rm erg\,s^{-1}$|⁠. Thus, taken into account the age of about 1.5 Myr, the total kinetic energy injected is more than |$10^{52} ~\rm erg$|⁠. On the other hand, the total CR energy derived in Section 4.2 is about |$10^{50}~\rm erg$| even if we assume the diffusion of CRs near this region is as fast as in the Galactic plane, which may imply an acceleration efficiency of at most 1 per cent. Indeed, the diffusion coefficient near the CR sources is expected to be much smaller (Aharonian et al. 2019), and in this case, the acceleration efficiency can be orders of magnitude smaller.

The wind power in Danks 1 and Danks 2 is as strong as in other powerful YMCs, such as Cygnus OB2 (Ackermann et al. 2011; Aharonian et al. 2019) and NGC 3603 (Yang & Aharonian 2017), but the gamma-ray luminosity of about |$10^{35} \rm erg\,s^{-1}$| is nearly one order of magnitude smaller, which also imply a lower CR acceleration efficiency in this region. Furthermore, this system is similar to W40 (Sun et al. 2020b), regarding the very young age (about 2 Myr) and high gas density (⁠|$n\sim 100~\rm cm^{-3}$|⁠). The gamma-ray emission in W40 also reveals a spectral index of about |$-2.4$| (Sun et al. 2020b), which is similar to the results we derived here, but significantly softer than other YMCs (⁠|$2.1-2.2$|⁠) (Aharonian et al. 2019). For W40 the derived total CR energy budget is of |$10^{47}~\rm erg$| which also indicates a quite low acceleration efficiency than in other systems. So W40 and Danks 1/Danks 2 may represent a new sub-type of YMC systems, in which CRs are accelerated in a lower efficiency and softer spectrum.

Nevertheless, we cannot completely exclude a leptonic origin for this source. The gamma-ray spectra can also be explained by IC or bremsstrahlung emission. And the required total electron energy above |$2\rm ~GeV$| is |$1.3\times 10^{48}\rm ~erg$| for IC emission or |$1.3\times 10^{46}\rm ~erg$| for bremsstrahlung, respectively. Our model on the radial distribution of CRs is also limited due to the detector’s performance. Hence future gamma-ray observations with higher angular resolution and multiwavelength study are crucial to understanding the particle acceleration and confinement in these systems. In this regard the Einstein Probe (Yuan et al. 2022) which is already launched and the forthcoming CTA (Cherenkov Telescope Array Consortium 2019) can provide unique information to pinning down the radiation mechanism and shed light on the mystery of the origin of CRs.

ACKNOWLEDGEMENTS

Rui-zhi Yang is supported by the NSFC under grant 12 041 305 and 12393854. Ruizhi Yang gratefully acknowledge the support of Cyrus Chun Ying Tang Foundations. Bing Liu acknowledges the support from the NSFC under grant 12103049.

DATA AVAILABILITY

The Fermi-LAT data used in this work is publicly available, which is provided online by the NASA-GSFC Fermi Science Support Center.⁵ We make use of the CO data⁶ to derive the H|$_{2}$|⁠. The data from Planck legacy archive⁷ are used to derive the H ii. The H i data are from the HI4PI⁸

Footnotes

REFERENCES