Galactic cosmic ray propagation: sub-PeV diffuse gamma-ray and neutrino emission

Abstract

The Tibet ASγ experiment just reported their measurement of sub-PeV diffuse gamma-ray emission from the Galactic disk, with the highest energy up to 957 TeV. These diffuse gamma rays are most likely the hadronic origin by cosmic ray (CR) interaction with interstellar gas in the galaxy. This measurement provides direct evidence to the hypothesis that the Galactic Cosmic Rays (GCRs) can be accelerated beyond PeV energies. In this work, we try to explain the sub-PeV diffuse gamma-ray spectrum with different CR propagation models. We find that there is a tension between the sub-PeV diffuse gamma-ray and the local CR spectrum. To describe the sub-PeV diffuse gamma-ray flux, it generally requires larger local CR flux than measurement in the knee region. We further calculate the PeV neutrino flux from the CR propagation model. Even all of these sub-PeV diffuse gamma rays originate from the propagation, the Galactic Neutrinos (GNs) only account for less than ~ 15% of observed flux, most of which are still from extragalactic sources.

Graphical abstract

Keywords

Galactic cosmic ray / diffuse gamma ray / neutrino

1 Introduction

In the all-particle spectrum, the most prominent feature is the so-called “knee” structure at ~ 4 PeV, where the spectrum exhibits a slight steepening, with the change in slope [1–5]. Meanwhile above the knee, a tendency to the heavier nuclei has been noted [6, 7]. It is believed that the cosmic rays (CRs) less than PeV energies are principally of galactic origin. However, which Galactic objects could accelerate PeV CRs, has long been in dispute. The supernova remnants (SNRs) had long been proposed as the sources of Galactic Cosmic Rays (GCRs) [8]. But long gazes into most SNRs associated with the molecular clouds indicate that the cut-off energy in the gamma-ray spectrum is less than 10 TeV, which means the inferred maximum cosmic ray (CR) energy is below ~ 100 TeV.

Massive strides about the PeV sources have been made in recent years. The Galactic center, which embodies a supermassive black hole, Sagittarius A*, with mass up to 4 × 10⁶ solar mass, is regarded as a PeVatron candidate. The observation of anisotropy phase also seems to hint at this [9–11]. The HESS collaboration performed deep observations of the central molecular zone surrounding the Galactic center. They detected a single power-law gamma-ray spectrum up to tens of TeV, without a break or a cutoff. This partly indicates a PeV proton accelerator at the Galactic center [12]. Early in this year, the Tibet ASγ experiment reported their observation of SNR G106.3+2.7 [13]. They successfully detected the 100 TeV gamma-ray photons, which are well correlated with the nearby molecular cloud. Therefore the morphological feature favours a hadronic origin of such very-high energy gamma rays. The young massive star clusters are another possible acceleration sites. Most recently, HAWC published their gamma-ray observation of the Cyngnus Cocoon, a well-known superbubble, with energy ranging from 1 to 100 TeV [14]. The 100 TeV gamma-ray photons may originate from the accelerated CRs in the Cyg OB2, an active star-forming region. These recent discoveries provide hard evidence concerning the galactic origin of PeV CRs.

Apart from searching for point sources, the observations of 100 TeV diffuse gamma-ray emission could also be served as the evidence of PeV CRs origin. At such high energy, the gamma-ray photons are expected to have interactions with the low energy background radiation, so that the extragalactic flux is strongly suppressed. On the other hand, the all-electron spectrum shows a clear steepening above 1 TeV, where the spectral index changes to ~ −3.9 [15, 16]. Coupled with the energy loss during propagation of electrons, the 100 TeV gamma-ray flux from the inverse Compton scattering off electrons could be safely neglected. Such energetic gamma rays are clearly generated from the π⁰-decay process during the hadronic interaction of the CR nuclei with the interstellar medium (ISM). Previously, only CASA-MIA [17] and KASCADE [18] experiments set upper limits on the diffuse gamma rays above 20 TeV. At present, the Tibet ASγ experiment managed to detect the sub-PeV diffuse gamma rays in the Galactic disk, with energy ranging from 100 TeV to 1 PeV [19]. The observed highest energy is up to 957 TeV, very close to 1 PeV. This discovery indicates that the Galactic sources could accelerate CRs to at least ~ 10 PeV. Some subsequent works had been carried out to study the origin of sub-PeV diffuse gamma-ray emission [20–23].

During the p–p interactions of CR with ISM, the generated gamma-ray photons are simultaneously accompanied with the neutrinos, which are the decay products of the π^±s. Thereupon, the high-energy neutrinos are also regarded as a good probe to the hadronic interaction. The astrophysical neutrinos have been detected by the IceCube Neutrino Observatory [24–26]. The overall distribution of the neutrino event samples is consistent with an isotropic distribution, which demonstrates that they are held to be predominantly extragalactic. The directional searches also have found an excess from a starburst galaxy [27] and the neutrino emission associated with a blazar [28, 29]. However, the proportion of Galactic contribution is still uncertain. A separate fit of the Northern and Southern hemisphere signals in the four-year signal shows a preference to a harder spectrum in the Northern hemisphere [30], which could potentially be due to the presence of a softer contribution of the flux from the inner galaxy in the Southern hemisphere.

The sub-PeV diffuse gamma-ray emission effectively traces the spatial distribution of remote CRs, so they could be well applied to testify the available propagation models. Besides, it is helpful to unveil the origin of knee region [31]. In this work, we investigate the propagation origin of sub-PeV diffuse gamma rays. We find, when the sub-PeV diffuse gamma-ray flux is entirely attributed to p–p interactions of CRs, and, based on current under-standings of CR transport and gas distribution, there is a tension between sub-PeV diffuse gamma rays and local CRs measurement. To explain the sub-PeV diffuse gamma rays measured by the Tibet ASγ experiment, the calculated local CR flux inevitably exceeds the observed CR flux above PeV energies. One possibility is that the propagated spectrum close to the Galactic center is harder as for the current propagation models, whereas the propagation around the solar system is still unchanged. Meanwhile these sub-PeV gamma-ray photons from point sources, especially unresolved, may not be entirely subtracted. We further evaluate the Galactic contribution to the diffuse neutrino flux and find at most 15% of observed flux come from GCRs propagation. It should be noted that Lipari and Vernetto [32] have made the prediction for the diffuse gamma rays above hundreds of TeV. But the spatial distribution of CRs is extrapolated according to the local observations. In this work, the whole CR spatial distribution is evaluated by solving the diffusion equation numerically.

2 Model description

2.1 Homogeneous diffusion

In the conventional propagation model, the diffusion process is supposed to be homogeneous and isotropic, so the diffusion coefficient is only a function of rigidity \mathcal{R} = p/Ze, namely,

\begin{array}{}\text{(1)}& D(\mathcal{R})=D_0{\beta }^{\eta }{\left(\frac{\mathcal{R}}{{\mathcal{R}}_0}\right)}^{{\delta }_0}\mathrm{,}\end{array}

where R₀ = 4 GV, β is the particle’s velocity in unit of light speed, η is a phenomenological constant in order to fit the low-energy data, D₀ and δ₀ are constants representing the diffusion coefficient and its high-energy rigidity dependence in the outer halo respectively. Usually, the power index δ₀ is taken from 0.3 to 0.6, as inferred from the fitting of boron-to-carbon ratio [33]. After propagation, the CR spectrum falls off as a single power-law, \varphi \propto {\mathcal{R}}^{-\nu -\delta } , where ν is the power index of CR spectrum at source.

2.2 Spatial-dependent propagation

However the spectral hardening of CR nuclei above ~ 200 GV [34] as well as the anisotropy observation [11] severely challenge the conventional homogeneous propagation model [35, 36]. The spatial-dependent propagation (SDP) was intially introduced to account for the excess of CR nuclei [37]. Afterwards, it is further applied to large-scale anisotropy [38, 39] and diffuse gamma-ray [40] observations. For a comprehensive introduction, one can refer to [41] and [42].

In contrast to the homogeneous diffusion, the entire diffusive halo in the SDP model is split into two zones characterized by diverse diffusion properties, i.e., inner halo (IH) and outer halo (OH). The Galactic disk and its surrounding areas within a few hundred parsecs are called IH, where the diffusion is slower and relevant to the radial distribution of sources. The diffusion in extended regions outside of IH, i.e., OH, is faster and approaches to the conventional propagation. The diffusion coefficient D in the whole region is thus parameterized as

\begin{array}{}\text{(2)}& D_{xx}(r\mathrm{,}z\mathrm{,}\mathcal{R})=D_{0}F(r\mathrm{,}z){\beta }^{\eta }{\left(\frac{\mathcal{R}}{{\mathcal{R}}_0}\right)}^{{\delta }_{0}F(r\mathrm{,}z)}\mathrm{,}\end{array}

where r and z are cylindrical coordinate, the other parameters are same as the conventional propagation model. The spatial dependence function F(r, z) is shown as follows:

\begin{array}{}\text{(3)}& F(r\mathrm{,}z)=\{\begin{array}{ll}g(r\mathrm{,}z)+[1-g(r\mathrm{,}z)]{\left(\frac{z}{\xi L}\right)}^n\mathrm{,}\hfill & |z|\le \xi L\hfill \\ 1\mathrm{,}\hfill & |z|>\xi L\hfill \end{array}\mathrm{,}\end{array}

in which g(r, z) is N_m/[1+f(r, z)], and f(r, z) is the source density distribution. N_m is a constant characterized the particle’s diffusion velocity in the inner halo. L is the halfthickness of the propagation cylinder, and ξL is the halfthickness of the inner halo. The factor {(\frac{z}{\xi L})}^n describes the smoothness of the parameters at the transition between the two halos. Note that the spatial dependence of the diffusion coefficient is phenomenologically assumed. Physically it may be related with the magnetic field distribution, or possibly the turbulence driven by CRs [43]. We adopt the diffusion reacceleration model in this work, with the reacceleration being described by a diffusion in the momentum space. The momentum diffusion coefficient, D_pp, correlates with D_xx via D_{pp}{D}_{xx}=\frac{4p^2{v}_A^2}{3\delta (4-{\delta }^2)(4-\delta )} , where v_A is the Alfvén velocity, p is the momentum, and δ is the rigidity dependence slope of the spatial diffusion coefficient [44]. The spatial distribution of CRs is obtained by numerically solving the diffusion equation with the DRAGON package. The diffuse gamma-ray distribution around the Galactic disk is calculated based on the CRs spatial distribution obtained with DRAGON. The employed gas distribution follows the default setup in the GALPROP.

2.3 Background source distribution

SNRs are considered to be the most plausible candidates for the acceleration of CRs. The spatial distribution of SNRs are approximated as an axisymmetric form parameterized as

\begin{array}{}\text{(4)}& f(r\mathrm{,}z)={\left(\frac{r}{r_{\odot }}\right)}^{\alpha }\exp \left[-\frac{\beta (r-r_{\odot })}{r_{\odot }}\right]\exp \left(-\frac{|z|}{z_s}\right)\mathrm{,}\end{array}

where r_⊙ ≡ 8.5 kpc represents the distance from the Galactic center to the solar system. Parameters α and β are taken to be 1.69 and 3.33 [45]. The density of the SNR distribution decreases exponentially along the vertical height from the Galactic plane, with z_s = 200 pc.

In the calculations below, the injection spectra of all propagation models are assumed to have a power-law plus a high-energy exponential cutoff, i.e.,

\begin{array}{}\text{(5)}& q_i(\mathcal{R})=q_0^i{\left(\frac{\mathcal{R}}{{\mathcal{R}}_{\text{br}}}\right)}^{{\nu }_i}\exp \left(-\frac{\mathcal{R}}{{\mathcal{R}}_{\text{c}}}\right)\mathrm{,}\end{array}

for the i-th composition. Here q₀ is the normalization factor, ν is the spectral incides, {\mathcal{R}}_{\text{br}} is break rigidity fixed at 2.2 GV, {\mathcal{R}}_{\text{c}} is the cutoff rigidity.

2.4 Nearby source

The propagation of particles from the nearby source is calculated using the Green’s function method, assuming a spherical geometry with infinite boundary conditions. Assuming instantaneous injection from a point source, the CR density as a function of space, rigidity, and time can be calculated as

\begin{array}{}\text{(6)}& \varphi (r\mathrm{,}\mathcal{R}\mathrm{,}t)=\frac{q_{\text{inj}}(\mathcal{R})}{{\left(\sqrt{2\pi }\sigma \right)}^3}\exp \left(-\frac{r^2}{2{\sigma }^2}\right)\mathrm{,}\end{array}

where q_{\text{inj}}(\mathcal{R}) is the injection spectrum as a function of rigidity, \sigma (\mathcal{R}\mathrm{,}t)=\sqrt{2D(\mathcal{R})t} is the effective diffusion length within time t. The diffusion coefficient D(\mathcal{R}) takes the solar system value of Eq. (2). The function form of q_{\text{inj}}(\mathcal{R}) is assumed to be power-law with an exponential cutoff, q_{\text{inj}}(\mathcal{R})=q_0{\mathcal{R}}^{-\alpha }\exp (-\mathcal{R}/{\mathcal{R}}_{\text{c}}^{\prime }) . The normalization q₀ is determined through fitting to the CR energy spectra. The distance and age of the local source are set to be d = 330 pc and τ = 3.4×10⁵ years [38], respectively.

3 Results

In this work, to make a detailed study of the propagation origin of sub-PeV diffuse gamma rays, we compare two kinds of common propagation scenarios, i.e., homogeneous diffusion (HD) and spatial-dependent propagation (SDP). We have introduced two SDP models according to the origin of spectral harding. In the model SDP-A, the excess of CR nuclei above 200 GV is regarded as a local effect, which originate from a local SNR. A nearby SNR is also favored in order to describe the evolution of anisotropy amplitude and phase with energy [38, 39, 46]. For model SDP-B, the excesses chiefly originate from the spatial variation of diffusion coefficient. The excess of diffuse gamma-ray emission at Galactic plane above a few GeV has also been well accounted for in model SDP-B [40], due to the spectral hardening of CRs in the whole diffusive halo.

3.1 CR energy spectra

Before calculating the diffuse gamma-ray distribution from the π⁰-decay, the CRs spatial distribution in the Galactic halo have to be evaluated by certain propagation and injection parameters. For this purpose, the local observations of CR energy spectrum have to be fitted first, in order to obtain the propagation and injection parameters. In the HD propagation, the essential propagation parameters includes D₀, δ₀ and L. As for the SDP models, three additional parameters have to be involved, i.e., N_m, ξ and n. Figure 1 shows the fitting to the latest B/C ratio published by AMS-02 experiment [47], and the corresponding propagation parameters are listed in Table 1. Compared with the conventional propagation, the SDP models anticipates a flattening for the B/C ratio above hundreds of GeV.

Given the propagation parameters, we further fit the local CR energy spectra to obtain the source’s injection parameters. The calculated proton and helium spectra are illustrated in Fig. 2, both of which are the major elements of CRs less than 10 PeV. The injection parameters are listed in Table 2. The green lines are the propagated proton and helium spectra in the HD model. For the HD model, we do not taken into account the local source, since the nuclei excess between 200 GV and tens of TeV does not affect the sub-PeV gamma-ray photon productivity. Compared with the observations less than 100 TeV, both proton and helium spectra measured by the KASCADE experiment indicate a visible steepening. In this work, the cutoff rigidity of different elements is assumed to have Z-dependence. To fit KASCADE observations, the cutoff rigidity of proton and helium is set to 7 PV for HD and SDP-A model, but 3.5 PV for SDP-B model.

The SDP-A model is marked as the red lines, in which the dash-dot, dash and solid lines denote the fluxes from background sources, local SNR and sum of them respectively. To reproduce the softening at 20 TeV [48, 49], the cutoff rigidities of local CRs are set to 28 TV. It has been shown that this break is relevant to the peak in the anisotropy amplitude at ~ 10 TeV and the flip of anisotropy phase at ~ 100 TeV. In this SDP model, energetic CRs principally propagate within the IH. Therefore compared with HD, the background spectra (i.e., red dash-dot lines) gradually harden above 10⁴ GeV, so that the calculated fluxes above 10⁴ GeV are higher than HD model. The corresponding sub-PeV diffuse gamma-ray flux is expected to be higher in the SDP model.

The blue lines are the calculated proton and helium fluxes in SDP-B model. Here the parameter N_m is 0.49, which is smaller than 0.51 of SDP-A model, so the propagated spectrum and corresponding B/C ratio are harder than those of SDP-A model. The cutoff rigidity is set to 3.5 PV in order to fit the KASCADE data. But the observed softening structure at around tens of TeV could not be explained, and both proton and helium fluxes inevitably exceed those of HD and SDP-A models above 100 TeV. As we will show below, these extra fluxes are important to account for the observed sub-PeV diffuse gamma rays.

3.2 All-particle spectrum

The diffuse sub-PeV photons principally originate from the CRs at PeV energies, and an amount of ground-based observations of all-particle spectrum have been performed at this energy range with high precision. To calculate all-particle spectrum, we fit other major components, namely C, N, O, Ne, Mg, Si and Fe respectively. All of them conform with the observations at lower energy separately and are extrapolated the model to the knee region. The injection parameters are summarized in Table 2. The cutoff rigidity of nuclei is Z-dependent in order to fit the knee region. It is worth noting that at knee region, which is most relevant to the Tibet sub-PeV gamma-ray observation, the major composition of CR flux is proton and helium, so the fluxes of heavier composition do not take leading role.

In Fig. 3, we compare the calculated all-particle spectra with observations. We could see that the expected all-particle spectrum in HD model is lower than the observations by simply extrapolating the low-energy spectrum assuming a single power-law injection spectrum. This would make the expected diffuse sub-PeV flux well below the observations, as we would see. The all-particle spectrum could be accounted for until ~ 10 PeV in both SDP models. Beyond that energy, the Galactic contribution drops significantly, and these CRs may originate from the extragalactic sources. Furthermore, the CRs above 10 PeV essentially do not contribute significantly to the sub-PeV gamma-ray photons during the interaction with the interstellar medium. Therefore, we do not further study the origin of such high energy CRs. It is worth noting that the all-particle flux in SDP-B model approaches to the upper limits of measurements above PeV energies.

3.3 Diffuse gamma-ray emission

After reproducing the local CR observations, we further calculate the spatial distribution of CRs and the corresponding diffuse gamma-ray spectra. Above ~ 20 TeV, the gamma-ray photons would have interactions with lowenergy Galactic interstellar radiation field (ISRF) and the observed gamma-ray flux is expected to be strongly suppressed [62, 63]. Figure 4 shows the diffuse gamma-ray spectra calculated by the three propagation models. The solid and dash lines are the diffuse gamma-ray spectra with/without gamma-ray attenuation from pair production respectively. As illustrated in the figure, above 100 TeV, the gamma-ray fluxes show observable attenuation.

In HD model, the calculated diffuse gamma-ray fluxes are well below the observations. Compared with HD, the gamma-ray spectra above 1 TeV show noticeable bumps in both SDP models so that the flux has been greatly enhanced correspondingly. This results from the spectral hardening of CRs between 100 TeV and 10 PeV in the entire halo. Meanwhile compared with SDP-B model, the diffuse gamma-ray flux in SDP-A model is lower, which means the CR flux in the whole Galaxy is always lower than SDP-B model above 100 TeV. This is because the excess of CR flux become softening at ~ 20 TeV in SDP-A model, and is attributed to a local effect. Especially close to the direction of anti-Galactic center, i.e., 50° < l < 200°, the gamma-ray flux above 100 TeV is significantly lower than Tibet observation. Furthermore, the flux at 1 TeV is slightly lower than ARGO measurements.

Only the expectations by the SDP-B model could marginally reaches the lower limits of both ARGO-YBJ and Tibet AS+MD observations at 25° < l < 100°. This is due to the enhanced CR flux from 10 TeV to 10 PeV. Meanwhile in SDP-B model, the hardening is regarded as the propagation effect, which enables the CR flux above 10 PeV in the entire diffusive halo to be augmented overall. Close to the Galactic center, the propagated CR energy spectra is harder. But as we have seen that, the CR proton flux exceeds the available observations and the softening at tens of TeV could not be reproduced.

As can be seen, the origin of sub-PeV diffuse gamma-ray flux has a tension with local CR observations more or less. When the local CR observations are reproduced, the calculated gamma-ray flux is inadequate to explain the Tibet observations. But if the sub-PeV gamma-ray flux is accounted for, the required CR flux surely exceeds the local CR energy spectra. This is obvious at 50° < l < 200°, which requires an unexpectedly large CR flux at anti-Galactic direction.

3.4 Diffuse neutrino spectrum

Based on the sub-PeV gamma-ray observation, we could reckon the Galactic contribution below PeV. On average, the p–p collision produces nearly one-third neutral pions and two-thirds charged pions. Each neutral pion decays into a pair of gamma rays, whereas each charged pion decays into two muon neutrinos and one electron neutrino (here we do not distinguish between neutrinos and anti-neutrinos). The initial neutrino flavor ratio is approximately ν_e : ν_µ : ν_τ = 1 : 2 : 0 from charged pion decay. After travelling, the flavor ratio is transformed to ν_e : ν_µ : ν_τ = 1 : 1 : 1 due to the vacuum neutrino oscillation. The typical neutrino energy from charged pion decay is approximately half of the gamma-ray photon from neutral pion decay. The gamma-ray spectrum at source is \text{d}{\text{Φ}}_{\nu }/\text{d}{E}_{\nu }={\varphi }_{\nu }{E}_{\nu }^{-\text{Γ}} , the resulting neutrino spectrum is shifted relative to the gamma-ray spectrum [65], i.e.,

\begin{array}{}\text{(7)}& \frac{\text{d}{\text{Φ}}_{\nu }}{\text{d}{E}_{\nu }}={\varphi }_{\nu }{E}_{\nu }^{-\text{Γ}}={\left(\frac{1}{2}\right)}^{\text{Γ}-1}{\varphi }_{\gamma }{E}_{\nu }^{-\text{Γ}}\mathrm{.}\end{array}

In Fig. 5, we evaluate the diffuse neutrino flux in three propagation models respectively. Compared with the fitting result of IceCube observations, which has a power index of ~ 2.87, the spectrum of Galactic diffuse neutrinos is softer above hundreds of TeV. Because of the cut-off energy of background injection spectrum of SDP-B model is smaller than another two propagation models, it can be seen that the calculated diffuse neutrino flux of SDP-B model is lower than that of the other two models above ~ 500 TeV. While it is larger than the other two models with energy less than ~ 500 TeV. Even in this case, the galactic diffuse neutrinos generated during CRs propagation, only take up to at most 15% of latest observations and are less than the 90% C.L. upper limit for the neutrino flux of the Galactic halo measured by the ANTARES neutrino telescope and IceCube.

4 Conclusion

The argument that the CRs below knee region are accelerated by the Galactic sources, has long been the lack of clear evidence. Most recently, the Tibet ASγ experiment reported their findings of the diffuse gamma ray photons between 100 TeV and 1 PeV, nearby the Galactic disk. These sub-PeV photons are likely the π⁰-decay products from the hadronic interactions between the CR nuclei and interstellar materials. The Galactic origin of CRs with energy up to PeV have a definite conclusion along with observations of Galactic center [12, 69], SNR G106.3+2.7 [13, 70], the massive young cluster Westerlund 1 [69, 71], Cygnus OB2 [69], gamma-ray binary source (LS 5039, PSR B1259-63, LSI +61°303, PSR J2032+4127) [72] as well as Cyngnus Cocoon [14] and so on.

In this work, we study the propagation origin of these sub-PeV diffuse gamma rays by considering three competing propagation models. One is homogeneous diffusion and the other two are SDP models. Compared with the homogeneous diffusion, the SDP-A model could generate more gamma rays above 10 TeV around the disk due to the flattenning propagated CR spectra above tens of TeV. Nonetheless, to meet the observed sub-PeV gamma-ray flux at 25° < l < 100°, the local CR flux in the SDP plus local source model has to slightly exceed the observation of the all-particle spectrum at around 10 PeV. And the calculated gamma rays are still inadequate for the measurements of 1 TeV at 25°< l < 100° and above 100 TeV at 50° < l < 200°. As for the SDP-B model, all above gamma-ray observations could be accounted for. However such a model meets up with severe challenges. The local CR flux would be greatly enhanced, which far exceeds the all-particle measurements. Meanwhile, the softening of CR nuclei at ~ 20 TeV could not be reproduced, much less the observations of large-scale anisotropies.

In the current propagation models, there is a tension between sub-PeV gamma rays observation and local CRs measurements. One possibilities is that the propagated spectrum is harder close to the Galactic center, whereas the propagation around the solar system is still unchanged, as indicated by the analysis of [73] for the radial distribution of diffuse gamma-ray. Therefore the observed sub-PeV gamma rays and local CRs could be simultaneously satisfied by the SDP model. Meanwhile these sub-PeV gamma-ray photons from the Galactic disk may not entirely come from the CR propagation. Most of CR sources and interstellar medium are located around the Galactic disk. There may be some of gamma-ray photons from CR sources, despite that the events within 0.5° from the known TeV sources have been subtracted. This has also been confirmed by the Tibet observation in fact. Above 398 TeV, 4 events at 50° < l < 200° are detected less than 4° from the center of the Cygnus cocoon, which is just proved to be a PeV source by HAWC experiment. We hope more extensive research of diffuse gamma rays and the observation of single CR composition between 100 TeV and 10 PeV with high precision could testify our conclusion.

We evaluate the corresponding neutrino flux from the Galactic comic ray propagation. We find that even if the observed sub-PeV diffuse gamma rays could be accounted for, the Galactic halo could contribute ~ 15% of observed neutrino flux.

References

[1] M. Nagano, T. Hara, Y. Hatano, N. Hayashida, S. Kawaguchi, K. Kamata, T. Kifune, and Y. Mizumoto, Energy spectrum of primary cosmic rays between 10^14.5 and 10¹⁸ eV, J. Phys. G 10(9), 1295 (1984)

[2] M. A. K. Glasmacher, M. A. Catanese, M. C. Chantell, et al., The cosmic ray energy spectrum between 10¹⁴ and 10¹⁶ eV, Astropart. Phys. 10(4), 291 (1999)

[3] M. Aglietta, B. Alessandro, P. Antonioli, F. Arneodo, L. Bergamasco, et al., The cosmic ray primary composition in the “knee” region through the EAS electromagnetic and muon measurements at EAS-TOP, Astropart. Phys. 21(6), 583 (2004)

[4] T. Antoni, W. D. Apel, A. F. Badea, K. Bekk, A. Bercuci, et al., KASCADE measurements of energy spectra for elemental groups of cosmic rays: Results and open problems, Astropart. Phys. 24(1–2), 1 (2005)

[5] M. Amenomori, X. J. Bi, D. Chen, S. W. Cui, Danzengluobu, et al., The all-particle spectrum of primary cosmic rays in the wide energy range from 10¹⁴ to 10¹⁷ eV observed with the Tibet-III air-shower array, Astrophys. J. 678(2), 1165 (2008)

[6] K. H. Kampert and M. Unger, Measurements of the cosmic ray composition with air shower experiments, Astropart. Phys. 35(10), 660 (2012)

[7] R. Aloisio, P. Blasi, I. De Mitri, and S. Petrera, Selected topics in cosmic ray physics, arXiv: 1707.06147 (2017)

[8] W. Baade and F. Zwicky, Cosmic rays from super-novae, Contributions from the Mount Wilson Observatory 3, 79 (1934)

[9] R. Abbasi, Y. Abdou, T. Abu-Zayyad, M. Ackermann, J. Adams, et al., Observation of anisotropy in the Galactic cosmic-ray arrival directions at 400 TeV with IceCube, Astrophys. J. 746(1), 33 (2012)

[10] M. G. Aartsen, K. Abraham, M. Ackermann, J. Adams, J. A. Aguilar, et al., Anisotropy in cosmic-ray arrival directions in the southern hemisphere based on six years of data from the IceCube detector, Astrophys. J. 826(2), 220 (2016)

[11] M. Amenomori, X. J. Bi, D. Chen, T. L. Chen, W. Y. Chen, et al., Northern sky Galactic cosmic ray anisotropy between 10 and 1000 TeV with the Tibet air shower array, Astrophys. J. 836(2), 153 (2017)

[12] HESS Collaboration, Acceleration of petaelectronvolt protons in the Galactic centre, Nature 531(7595), 476 (2016)

[13] The Tibet ASγ Collaboration, M. Amenomori, Y. W. Bao, et al., Potential PeVatron supernova remnant G106.3+2.7 seen in the highest-energy gamma rays, Nat. Astron. 5, 460 (2021)

[14] A. U. Abeysekara, A. Albert, R. Alfaro, et al., HAWC observations of the acceleration of very-high-energy cosmic rays in the Cygnus Cocoon, arXiv: 2103.06820 (2021)

[15] DAMPE Collaboration, Direct detection of a break in the teraelectronvolt cosmic-ray spectrum of electrons and positrons, Nature 552(7683), 63 (2017)

[16] D. Kerszberg for the HESS Collaboration, The cosmic-ray electron spectrum measured with H.E.S.S. (2017)

[17] A. Borione, M. A. Catanese, M. C. Chantell, C. E. Covault, J. W. Cronin, et al., Constraints on gamma-ray emission from the Galactic plane at 300 TeV, Astrophys. J. 493(1), 175 (1998)

[18] W. D. Apel, J. C. Arteaga-Velázquez, K. Bekk, M. Bertaina, J. Blümer, et al., KASCADE-Grande limits on the isotropic diffuse gamma-ray flux between 100 TeV and 1 EeV, Astrophys. J. 848(1), 1 (2017)

[19] M. Amenomori, Y. W. Bao, X. J. Bi, D. Chen, T. L. Chen, et al., First detection of sub-PeV diffuse gamma rays from the Galactic disk: Evidence for ubiquitous Galactic cosmic rays beyond PeV energies, Phys. Rev. Lett. 126(14), 141101 (2021)

[20] R. Y. Liu and X. Y. Wang, Origin of Galactic sub PeV diffuse gamma-ray emission: Constraints from highenergy neutrino observations, Astrophys. J. Lett. 914(1), L7 (2021)

[21] V. Vecchiotti, F. Zuccarini, F. L. Villante, and G. Pagliaroli, Unresolved sources naturally contribute to PeV γ-ray diffuse emission observed by Tibet ASγ, arXiv: 2107.14584 (2021)

[22] S. Koldobskiy, A. Neronov, and D. Semikoz, Pion decay model of the Tibet-ASγ PeV gamma-ray signal, Phys. Rev. D 104(4), 043010 (2021)

[23] P. P. Zhang, B. Q. Qiao, Q. Yuan, S. W. Cui, and Y. Q. Guo, Ultrahigh-energy diffuse gamma ray emission from cosmic-ray interactions with the medium surrounding acceleration sources, Phys. Rev. D 105(2), 023002 (2022)

[24] IceCube Collaboration, Evidence for high-energy extraterrestrial neutrinos at the IceCube detector, Science 342(6161), 1242856 (2013)

[25] M. G. Aartsen, R. Abbasi, Y. Abdou, M. Ackermann, J. Adams, et al., First observation of PeV-energy neutrinos with IceCube, Phys. Rev. Lett. 111(2), 021103 (2013)

[26] M. G. Aartsen, M. Ackermann, J. Adams, J. A. Aguilar, M. Ahlers, et al., Observation of high-energy astrophysical neutrinos in three years of IceCube data, Phys. Rev. Lett. 113(10), 101101 (2014)

[27] M. G. Aartsen, M. Ackermann, J. Adams, J. A. Aguilar, M. Ahlers, et al., Time integrated neutrino source searches with 10 years of IceCube data, Phys. Rev. Lett. 124(5), 051103 (2020)

[28] M. Aartsen, M. Ackermann, J. Adams, J. A. Aguilar, M. Ahlers, et al., Multimessenger observations of a flaring blazar coincident with high-energy neutrino IceCube-170922A, Science 361(6398), eaat1378 (2018)

[29] M. Aartsen, M. Ackermann, J. Adams, J. A. Aguilar, M. Ahlers, et al., Neutrino emission from the direction of the blazar TXS 0506+056 prior to the IceCube-170922A alert, Science 361(6398), 147 (2018)

[30] M. G. Aartsen, K. Abraham, M. Ackermann, J. Adams, J. A. Aguilar, et al., A combined maximum-likelihood analysis of the high energy astrophysical neutrino flux measured with Ice Cube, Astrophys. J. 809(1), 98 (2015)

[31] Y. Q. Guo, H. B. Hu, Q. Yuan, Z. Tian, and X. J. Gao, Pinpointing the knee of cosmic rays with diffuse PeV γ-rays and neutrinos, Astrophys. J. 795(1), 100 (2014)

[32] P. Lipari and S. Vernetto, Diffuse Galactic gamma-ray flux at very high energy, Phys. Rev. D 98(4), 043003 (2018)

[33] Q. Yuan, S. J. Lin, K. Fang, and X. J. Bi, Propagation of cosmic rays in the AMS-02 era, Phys. Rev. D 95(8), 083007 (2017)

[34] O. Adriani, G. C. Barbarino, G. A. Bazilevskaya, R. Bellotti, M. Boezio, et al., PAMELA measurements of cosmicray proton and helium spectra, Science 332(6025), 69 (2011)

[35] P. Blasi and E. Amato, Diffusive propagation of cos mic rays from supernova remnants in the Galaxy (II): Anisotropy, J. Cosmol. Astropart. Phys. 2012(1), 11 (2012)

[36] W. Liu, X. J. Bi, S. J. Lin, B. B. Wang, and P. F. Yin, Excesses of cosmic ray spectra from a single nearby source, Phys. Rev. D 96(2), 023006 (2017)

[37] N. Tomassetti, Origin of the cosmic-ray spectral hardening, Astrophys. J. Lett. 752(1), L13 (2012)

[38] W. Liu, Y.-Q. Guo, and Q. Yuan, Indication of nearby source signatures of cosmic rays from energy spectra and anisotropies, J. Cosmol. Astropart. Phys. 2019(10), 010 (2019)

[39] B.-Q. Qiao, W. Liu, Y.-Q. Guo, and Q. Yuan, Anisotropies of different mass compositions of cosmic rays, J. Cosmol. Astropart. Phys. 2019(12), 007 (2019)

[40] Y. Q. Guo and Q. Yuan, Understanding the spectral hardenings and radial distribution of Galactic cosmic rays and Fermi diffuse γ rays with spatially-dependent propagation, Phys. Rev. D 97(6), 063008 (2018)

[41] Y. Q. Guo, Z. Tian, and C. Jin, Spatial-dependent propagation of cosmic rays results in the spectrum of proton, ratios of P/P, and B/C, and anisotropy of nuclei, Astrophys. J. 819, 54 (2016)

[42] W. Liu, Y. H. Yao, and Y. Q. Guo, Revisiting the spatially dependent propagation model with the latest observations of cosmic-ray nuclei, Astrophys. J. 869(2), 176 (2018)

[43] P. Blasi, E. Amato, and P. D. Serpico, Spectral breaks as a signature of cosmic ray induced turbulence in the galaxy, Phys. Rev. Lett. 109(6), 061101 (2012)

[44] E. S. Seo and V. S. Ptuskin, Stochastic reacceleration of cosmic rays in the interstellar medium, Astrophys. J. 431, 705 (1994)

[45] G. Case and D. Bhattacharya, Revisiting the Galactic supernova remnant distribution, Astron. Astrophys. Suppl. 120, 437 (1996)

[46] M. Ahlers, Deciphering the dipole anisotropy of Galactic cosmic rays, Phys. Rev. Lett. 117(15), 151103 (2016)

[47] M. Aguilar, L. A. Cavasonza, G. Ambrosi, et al., Precision measurement of the boron to carbon flux ratio in cosmic rays from 1.9 GV to 2.6 TV with the alpha magnetic spectrometer on the International Space Station, Phys. Rev. Lett. 117(23), 231102 (2016)

[48] Y. S. Yoon, T. Anderson, A. Barrau, N. B. Conklin, S. Coutu, et al., Proton and helium spectra from the CREAM-III flight, Astrophys. J. 839(1), 5 (2017)

[49] Q. An, R. Asfandiyarov, P. Azzarello, P. Bernardini, X. J. Bi, et al., Measurement of the cosmic ray proton spectrum from 40 GeV to 100 TeV with the DAMPE satellite, Sci. Adv. 5(9), eaax3793 (2019)

[50] M. Aguilar, D. Aisa, B. Alpat, A. Alvino, G. Ambrosi, et al., Precision measurement of the proton flux in primary cosmic rays from rigidity 1 GV to 1.8 TV with the alpha magnetic spectrometer on the International Space Station, Phys. Rev. Lett. 114(17), 171103 (2015)

[51] M. Aguilar, L. A. Cavasonza, B. Alpat, et al., Observation of the identical rigidity dependence of He, C, and O cosmic rays at high rigidities by the alpha magnetic spectrometer on the International Space Station, Phys. Rev. Lett. 119(25), 251101 (2017)

[52] E. Atkin, V. Bulatov, V. Dorokhov, N. Gorbunov, S. Filippov, et al., First results of the cosmic ray NUCLEON experiment, J. Cosmol. Astropart. Phys. 07, 020(2017)

[53] W. D. Apel, J. C. Arteaga-Velázquez, K. Bekk, M. Bertaina, J. Blümer, et al., KASCADE-Grande measurements of energy spectra for elemental groups of cosmic rays, Astropart. Phys. 47, 54 (2013)

[54] M. G. Aartsen, M. Ackermann, J. Adams, J. A. Aguilar, M. Ahlers, et al., Cosmic ray spectrum and composition from PeV to EeV using 3 years of data from IceTop and IceCube, Phys. Rev. D 100(8), 082002 (2019)

[55] J. C. Arteaga-Velázquez, HAWC measurements of the energy spectra of cosmic ray protons, helium and heavy nuclei in the TeV range, arXiv: 2108.03208 (2021)

[56] J. R. Hörandel, On the knee in the energy spectrum of cosmic rays, Astropart. Phys. 19(2), 193 (2003)

[57] M. G. Aartsen, R. Abbasi, M. Ackermann, J. Adams, J. A. Aguilar, et al., Cosmic ray spectrum from 250 TeV to 10 PeV using IceTop, Phys. Rev. D 102(12), 122001 (2020)

[58] R. Alfaro, C. Alvarez, J. D. Álvarez, R. Arceo, J. C. Arteaga-Velázquez, et al., All particle cosmic ray energy spectrum measured by the HAWC experiment from 10 to 500 TeV, Phys. Rev. D 96(12), 122001 (2017)

[59] G. Di Sciascio, Measurement of the Cosmic Ray Energy Spectrum with ARGO-YBJ, arXiv: 1408.6739 (2014)

[60] A. D. Panov, J. H. Jr Adams, H. S. Ahn, G. L. Bashinzhagyan, et al., Energy spectra of abundant nuclei of primary cosmic rays from the data of ATIC-2 experiment: Final results, Bull. Russ. Acad. Sci, Physics 73(5), 564 (2009)

[61] E. V. Atkin, V. L. Bulatov, O. A. Vasiliev, A. G. Voronin, N. V. Gorbunov, et al., Energy Spectra of Cosmic-Ray Protons and Nuclei Measured in the NUCLEON Experiment Using a New Method, Astron. Rep. 63(1), 66 (2019)

[62] J.-L. Zhang, X.-J. Bi, and H.-B. Hu, Very high energy γ ray absorption by the Galactic interstellar radiation field, Astron. & Astrophys. 449, 641 (2006)

[63] I. V. Moskalenko, T. A. Porter, and A. W. Strong, Attenuation of very high energy gamma rays by the milky way interstellar radiation field, Astrophys. J. 640(2), L155 (2006)

[64] B. Bartoli, P. Bernardini, X. J. Bi, P. Branchini, A. Budano, et al., Study of the diffuse gamma-ray emission from the Galactic plane with ARGO-YBJ, Astrophys. J. 806(1), 20 (2015)

[65] M. D. Kistler and J. F. Beacom, Guaranteed and prospective Galactic TeV neutrino sources, Phys. Rev. D 74(6), 063007 (2006)

[66] R. Abbasi, M. Ackermann, J. Adams, et al., The IceCube high-energy starting event sample: Description and flux characterization with 7.5 years of data, arXiv: 2011.03545 (2020)

[67] S. Adrián-Martínez, A. Albert, M. André, M. Anghinolfi, G. Anton, et al., Constraints on the neutrino emission from the Galactic Ridge with the ANTARES telescope, Phys. Lett. B 760, 143 (2016)

[68] M. G. Aartsen, M. Ackermann, J. Adams, J. A. Aguilar, M. Ahlers, et al., Constraints on Galactic neutrino emission with seven years of IceCube data, Astrophys. J. 849(1), 67 (2017)

[69] F. Aharonian, R. Yang, and E. de Oña Wilhelmi, Massive stars as major factories of Galactic cosmic rays, Nat. Astron. 3(6), 561 (2019)

[70] P. Cristofari, The hunt for pevatrons: The case of supernova remnants, Universe 7(9), 324 (2021)

[71] A. M. Bykov, D. C. Ellison, P. E. Gladilin, and S. M. Osipov, Ultrahard spectra of PeV neutrinos from supernovae in compact star clusters, Mon. Not. R. Astron. Soc. 453(1), 113 (2015)

[72] A. M. Bykov, A. E. Petrov, M. E. Kalyashova, and S. V. Troitsky, PeV photon and neutrino flares from Galactic gamma-ray binaries, Astrophys. J. Lett. 921(1), L10 (2021)

[73] R. Yang, F. Aharonian, and C. Evoli, Radial distribution of the diffuse γ-ray emissivity in the Galactic disk, Phys. Rev. D 93(12), 123007 (2016)

[74] A. W. Strong and I. V. Moskalenko, Propagation of cosmic-ray nucleons in the galaxy, Astrophys. J. 509(1), 212 (1998)

[75] A. W. Strong, I. V. Moskalenko, and O. Reimer, Diffuse continuum gamma rays from the galaxy, Astrophys. J. 537(2), 763 (2000)

[76] C. Evoli, D. Gaggero, D. Grasso, and L. Maccione., Cosmic ray nuclei antiprotons and gamma rays in the galaxy: A new diffusion model, J. Cosmol. Astropart. Phys. 10, 018 (2008)

[77] C. Evoli, D. Gaggero, A. Vittino, G. Di Bernardo, M. Di Mauro, et al., Cosmic-ray propagation with DRAGON2 (I): Numerical solver and astrophysical ingredients, J. Cosmol. Astropart. Phys. 02, 015 (2017)

Rights & permissions

Higher Education Press 2022