On the relevance of prompt neutrinos for the interpretation of the IceCube signals

TL;DR

This paper emphasizes the importance of prompt neutrinos in interpreting IceCube signals, proposing a model to understand flavor-specific spectra and highlighting the need for consistent data analysis across different datasets.

Contribution

It introduces a phenomenological model for high-energy neutrino production via cosmic ray collisions and discusses its implications for IceCube data interpretation.

Findings

Prompt neutrinos can explain spectral differences between datasets

A specific dataset is crucial for detecting the prompt component

Consistent modeling is essential for accurate interpretation

Abstract

The IceCube collaboration has discovered a new, cosmic component of high-energy neutrinos. Although neutrino oscillations suggest that the cosmic neutrino spectrum is almost the same for every neutrino flavor, the attempts to reconstruct it, based on different analyses, lead to different energy spectra below 100 TeV. In this work, we propose a phenomenological model that, assuming collisions between cosmic rays and hadrons as the production mechanism of high-energy neutrinos, yields quantitative expectations for each neutrino flavor. We discuss the detectability of the prompt component of the atmospheric neutrino spectrum, pointing out the most relevant dataset, which has to be muon-neutrino depleted and to cover the energy region 10-100 TeV. We argue that the prompt component can cause the spectral difference between the High-Energy-Starting-Event (HESE) and the through-going muon datasets. Finally, we point out the need for adopting a consistent model for the interpretation of the data, stressing that a separate treatment of the different datasets is, by converse, a suboptimal procedure.

Tables

Table 1. Table 2 : 68 % percent 68 68\% CL intervals relative to the yearly rate of cascade-like events in the energy range 1   TeV times 1 TeV 1\text{\,}\mathrm{T}\mathrm{e}\mathrm{V} and 10   PeV times 10 PeV 10\text{\,}\mathrm{P}\mathrm{e}\mathrm{V} in IceCube, as computed with eq. ( 3.1 ), due to the different components of the three neutrino fluxes.

Component		${\Gamma }_{{\nu }_e}$	${\Gamma }_{{\nu }_{\mu }}$	${\Gamma }_{{\nu }_{\tau }}$	${\Gamma }_{\text{tot}}$
Conventional		160 – 210	420 – 570	0	580 – 780
Prompt		20 – 30	3 – 5	2 – 3	25 – 40
Cosmic		10 – 40	2 – 6	5 – 20	15 – 65

Table 2. Table 3 : The flux normalizations, in units of ​ 10 − 18   GeV − 1   cm − 2   s − 1   sr − 1 times E-18 times gigaelectronvolt 1 centimeter 2 second 1 steradian 1 {10}^{-18}\text{\,}{\mathrm{GeV}}^{-1}\text{\,}{\mathrm{cm}}^{-2}\text{\,}{\mathrm{s}}^{-1}\text{\,}{\mathrm{sr}}^{-1} , and slopes deriving from the astrophysical best fits to the cascade (C) and starting tracks (ST) samples from [ 31 ] , to the 6-years HESE sample from [ 35 ] and to the through-going muons (TM) sample from [ 2 ] . The numerical values of the flux normalizations for the cascade and starting tracks best fits have been obtained from figure 10 of [ 31 ] as they are not explicitly reported in the paper.

dataset		${\Phi }_{\text{astro}}$		${\gamma }_{\text{astro}}$
C		${2.2}_{-0.5}^{+0.6}$		$2.62\pm 0.08$
ST		${1.6}_{-1.0}^{+1.6}$		${2.43}_{-0.30}^{+0.28}$
HESE		$2.46\pm 0.8$		${2.92}_{-0.29}^{+0.33}$
TM		${0.90}_{-0.27}^{+0.30}$		$2.13\pm 0.13$

Table 3. Table 4 : The natural parameters obtained from Monte Carlo sampling the oscillation parameters from [ 30 ] according to their distribution.

ordering	$P_0$	$P_1$	$P_2$
NO	$0.113\pm 0.006$	${0.0345}_{-0.012}^{+0.010}$	${0.0075}_{-0.0038}^{+0.0045}$
IO		${0.0285}_{-0.057}^{+0.010}$	${0.008}_{-0.006}^{+0.005}$

Table 4. Table 5 : The survival/oscillation probabilities obtained from Monte Carlo sampling the oscillation parameters from [ 30 ] according to their distribution.

ordering	$P_{ee}$	$P_{e\mu }$	$P_{e\tau }$	$P_{\mu \mu }$	$P_{\mu \tau }$	$P_{\tau \tau }$
NO	$0.56\pm 0.01$	${0.25}_{-0.01}^{+0.02}$	${0.19}_{-0.02}^{+0.01}$	${0.37}_{-0.02}^{+0.01}$	$0.381\pm 0.005$	${0.43}_{-0.01}^{+0.02}$
IO		${0.23}_{-0.03}^{+0.04}$	${0.21}_{-0.04}^{+0.04}$	${0.39}_{-0.04}^{+0.04}$	$0.381\pm 0.006$	${0.40}_{-0.03}^{+0.03}$

Equations

$$\frac{d\Phi_{p,\text{He}}^{\text{exp2-cut}}}{dE}=N_{p,\text{He}}\left(\frac{E}{$10\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$}\right)^{-\gamma_{p,\text{He}}}\exp\left[-\left(\frac{E}{Z_{p,\text{He}}R_{\text{knee}}}\right)^{2}\right]$$

$$\frac{d\Phi_{p,\text{He}}^{\text{delta-slope}}}{dE}=N_{p,\text{He}}\left(\frac{E}{$10\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$}\right)^{-\gamma_{p,\text{He}}}\times\begin{dcases}1&E\leq Z_{p,\text{He}}R_{\text{knee}}\\ \left(\frac{E}{$10\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$}\right)^{-\alpha}&E>Z_{p,\text{He}}R_{\text{knee}}\end{dcases}$$

$$\frac{d\Phi_{\text{eg}}}{dE}=N_{\text{eg}}\left(\frac{E}{$100\text{\,}\mathrm{P}\mathrm{e}\mathrm{V}$}\right)^{-2.7}$$

$$\frac{d{\Phi }_{\text{tot}}^k}{dE}=\sum _{i=p,\text{He}}\frac{d{\Phi }_i^k}{dE}+\frac{d{\Phi }_{\text{eg}}}{dE}$$

$$\frac{d\Phi}{dE}=\Phi_{\text{astro}}\left(\frac{E}{$100\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$}\right)^{-\gamma_{\text{astro}}}$$

$$\Phi_{\text{astro}}^{\text{LW}}=2\times 10^{\pm 0.5}\times${10}^{-18}\text{\,}{\mathrm{GeV}}^{-1}\text{\,}{\mathrm{cm}}^{-2}\text{\,}{\mathrm{s}}^{-1}\text{\,}{\mathrm{sr}}^{-1}$\qquad\gamma_{\text{astro}}^{\text{LW}}=2.15\pm 0.10$$

$$\Phi_{\text{astro}}^{\text{IC},\mu}=0.90^{+0.30}_{-0.27}\times${10}^{-18}\text{\,}{\mathrm{GeV}}^{-1}\text{\,}{\mathrm{cm}}^{-2}\text{\,}{\mathrm{s}}^{-1}\text{\,}{\mathrm{sr}}^{-1}$\qquad\gamma_{\text{astro}}^{\text{IC},\mu}=2.13\pm 0.13$$

$$\frac{d\Phi_{\nu_{\mu}}}{dE}=0.90^{+0.30}_{-0.27}\times 10^{-18}\left(\frac{E}{$100\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$}\right)^{-2.14\pm 0.08}\;${\mathrm{GeV}}^{-1}\text{\,}{\mathrm{cm}}^{-2}\text{\,}{\mathrm{s}}^{-1}\text{\,}{\mathrm{sr}}^{-1}$$

$

$$\mathcal{L}(\mathbf{v})=\frac{1}{2\pi \sqrt{\det {\Sigma }^2}}\exp \left[-\frac{1}{2}(\mathbf{v}-{\mathbf{v}}^{\text{best}}{)}^T{\Sigma }^{-2}(\mathbf{v}-{\mathbf{v}}^{\text{best}})\right]$$

$$\mathbf{v}=\left(\begin{array}{c}{\Phi }_{\text{astro}}\\ {\gamma }_{\text{astro}}\end{array}\right){\Sigma }^2=\left(\begin{array}{cc}{\sigma }_{\Phi }^2& \rho {\sigma }_{\Phi }{\sigma }_{\gamma }\\ \rho {\sigma }_{\Phi }{\sigma }_{\gamma }& {\sigma }_{\gamma }^2\end{array}\right)$$

$$⟨\frac{d\Phi }{dE}⟩={\int }_{{\Phi }_{\text{astro}}^{\text{best}}-5{\sigma }_{\Phi }}^{{\Phi }_{\text{astro}}^{\text{best}}+5{\sigma }_{\Phi }}d{\Phi }_{\text{astro}}{\int }_{{\gamma }_{\text{astro}}^{\text{best}}-5{\sigma }_{\gamma }}^{{\gamma }_{\text{astro}}^{\text{best}}+5{\sigma }_{\gamma }}d{\gamma }_{\text{astro}}\mathcal{L}({\Phi }_{\text{astro}},{\gamma }_{\text{astro}})\frac{d\Phi }{dE}$$

$$\delta \left(\frac{d\Phi }{dE}\right)=\sqrt{⟨{\left(\frac{d\Phi }{dE}\right)}^2⟩-{\left(⟨\frac{d\Phi }{dE}⟩\right)}^2}$$

$$\frac{d{\Phi }_{{\nu }_e}}{dE}=R_{e\mu }\frac{d{\Phi }_{{\nu }_{\mu }}}{dE}\frac{d{\Phi }_{{\nu }_{\tau }}}{dE}=R_{\tau \mu }\frac{d{\Phi }_{{\nu }_{\mu }}}{dE}$$

$$P_{\ell {\ell }^{\prime }}=\sum _{i=1}^3|U_{\ell i}^2||U_{{\ell }^{\prime }i}^2|$$

$$R_{\ell {\ell }^{\prime }}=\frac{{\sum }_{{\ell }^{\prime \prime }}{P}_{\ell {\ell }^{\prime \prime }}{\xi }_{{\ell }^{\prime \prime }}^0}{{\sum }_{{\ell }^{\prime \prime }}{P}_{{\ell }^{\prime }{\ell }^{\prime \prime }}{\xi }_{{\ell }^{\prime \prime }}^0}\ell ,{\ell }^{\prime },{\ell }^{\prime \prime }=e,\mu ,\tau$$

$$R_{e\mu }=1.09_{-0.04}^{+0.03}{R}_{\tau \mu }=0.97_{-0.04}^{+0.03}$$

$$\frac{d{\Phi }_{\nu }(E_{\nu })}{d{E}_{\nu }}={\int }_{E_{\nu }}^{\infty }\frac{dE}{E}{\tilde{K}}_{\nu }(E_{\nu },E)\frac{d{\Phi }_{\gamma }(E)}{dE}$$

$$\frac{d{\Phi }_{\nu }(E_{\nu })}{d{E}_{\nu }}={\int }_0^1\frac{dx}{x}\tilde{K}(x)\frac{d{\Phi }_{\gamma }(x/E_{\nu })}{dE}$$

$$R_{\ell {\ell }^{\prime }}=\frac{{\zeta }_{{\nu }_{\ell }}(\gamma )}{{\zeta }_{{\nu }_{{\ell }^{\prime }}}(\gamma )}$$

$${\zeta }_{{\nu }_{\ell }}(\gamma )={\int }_0^{1}dx{x}^{\gamma -1}\left[{\tilde{K}}_{{\nu }_{\ell }}(x)+{\tilde{K}}_{{\overline{\nu }}_{\ell }}(x)\right]$$

$$R_{e\mu }=1.30\pm 0.05R_{\tau \mu }=0.92\pm 0.04$$

$$\Gamma_{\nu_{\ell}}=4\pi\times$1\text{\,}\mathrm{y}\mathrm{e}\mathrm{a}\mathrm{r}$\times\int_{$1\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$}^{$10\text{\,}\mathrm{P}\mathrm{e}\mathrm{V}$}dE\,\mathcal{A}_{\nu_{\ell}}(E)\,\frac{d\Phi_{\nu_{\ell}}}{dE}$$

$$P_{\ell {\ell }^{\prime }}=\sum _{i=1}^3|U_{\ell i}^2||U_{{\ell }^{\prime }i}^2|\ell ,{\ell }^{\prime }=e,\mu ,\tau$$

$$P_{\ell {\ell }^{\prime }}=\left(\begin{array}{ccc}1/3+2P_0& 1/3-P_0+P_1& 1/3-P_0-P_1\\ & 1/3+P_0/2-P_1+P_2& 1/3+P_0/2-P_2\\ & & 1/3+P_0/2+P_1+P_2\end{array}\right)$$

$$P_0$$

$$P_1$$

$$P_2$$

$$ϵ={\sin }^2{\theta }_{13}\alpha =\sin {\theta }_{13}\cos \delta \sin 2{\theta }_{12}\sin 2{\theta }_{23}\beta =\cos 2{\theta }_{23}\gamma =\alpha -\frac{\beta }{2}\cos 2{\theta }_{12}(1+ϵ)$$

$${\tilde{K}}_{{\nu }_{\ell }}=\sum _{{\ell }^{\prime }=e,\mu }{P}_{\ell {\ell }^{\prime }}{K}_{{\nu }_{{\ell }^{\prime }}}\ell =e,\mu ,\tau$$

$$K_{{\nu }_{\ell }}(x)={\alpha }_{\pi }\delta \left(x-(1-r_{\pi })\right)+{\alpha }_K\delta \left(x-(1-r_K)\right)+\{\begin{array}{ll}{x}^2({\beta }_0+{\beta }_{1}x)& x\le r_K\\ \sum _{n=0}^3{\chi }_n{x}^n& r_K<x<r_{\pi }\\ (1-x{)}^2({\delta }_0+{\delta }_{1}x)& x\ge r_{\pi }\end{array}$$

Full text

1 Introduction

After almost ten years of operation, the IceCube detector has provided unique and most important observations of neutrino events with energies ranging between $100\text{\,}\mathrm{G}\mathrm{e}\mathrm{V}$ and $10\text{\,}\mathrm{P}\mathrm{e}\mathrm{V}$, which resulted in the detection of an astrophysical component of neutrinos [1, 2] and in the measurement of the atmospheric components of the electron and muon neutrino spectrum [3, 4]. The prompt component of the atmospheric neutrino spectrum, which is expected to be produced in the decay of charmed mesons in the atmosphere, has not been measured yet.

On the other hand, some missing pieces of the high-energy neutrino jigsaw puzzle seem to be finally appearing: one notable astronomical coincidence [5] may hint that we are close to detecting the source of cosmic neutrinos. Moreover, two double cascade events attributable to astrophysical tau neutrinos and one candidate Glashow resonance event have been recently observed111See the talk by H. Niederhausen, on behalf of the IceCube Collaboration, at the XVIII International Workshop on Neutrino Telescopes (Venice, March 2019)..

The accepted set of assumptions on the astrophysical component, that are adopted for the interpretation of these findings, includes:

isotropy 2. 2.

standard three-flavor oscillations 3. 3.

an unbroken power law for the energy spectrum of the new component.

However, the cosmic neutrino spectrum resulting from the HESE analysis [6] is different from that obtained in the through-going muons analysis [2], as argued, for example, in [7, 8]. Moreover, if the spectrum from HESE were extrapolated down to energies lower than $100\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$, it would overshoot the gamma-ray diffuse measurements [9].

In this situation, we find it necessary to rely on theoretical guidance. In this work a primary cosmic-ray flux, fitted to the data of AMS-02 [10, 11] and KASCADE-Grande [12], is defined to numerically compute both conventional and prompt components of the atmospheric neutrino spectrum. The cosmic neutrino flux is modeled as a single-population power law, assuming $pp$ collisions in dense sources as the production mechanism. The ensuing expectations from theory are combined to the through-going muon analysis to make the astrophysical muon flux phenomenologically precise. The same muon neutrino flux is also used to predict the electron and tau cosmic neutrino fluxes. Credible regions are computed for all neutrino fluxes, which are in turn compared to the available measurements.

Finally, we comment on the results and address

the compelling issue of the low-energy softness of the HESE spectrum 2. 2.

the lack of detection of prompt neutrinos in the currently examined datasets 3. 3.

the methodological consistency of independent analyses of the various datasets.

2 The expected neutrino fluxes

In this section we describe our approach to the study of cosmic neutrinos, based upon theoretical considerations.

In figure 1 we anticipate some of our results, with the aim of emphasizing the region from $1\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ to $1\text{\,}\mathrm{P}\mathrm{e}\mathrm{V}$, where the prompt component is expected to be relevant. This is especially evident in the panels that show the expected electron (and tau) neutrino flux components. The normalization of the astrophysical component is an important ingredient of our model, and will be discussed in section 2.2.

2.1 Atmospheric neutrinos

Atmospheric neutrinos had not been studied at the highest energies before IceCube, but they can be predicted, within uncertainties linearly increasing from $\lesssim 10\%$ at $100\text{\,}\mathrm{G}\mathrm{e}\mathrm{V}$ to about 30% at $1\text{\,}\mathrm{P}\mathrm{e}\mathrm{V}$, from the observed flux of cosmic rays and the theory of strong interactions222These numbers refer to our particular computations which feature a fixed hadron interaction model, which is another source of uncertainty (around 10% at $100\text{\,}\mathrm{G}\mathrm{e}\mathrm{V}$ up to $30-40\%$ at $100\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ [13]).. Many (semi-)analytical computations [14, 15, 16, 17, 18, 19, 20, 21, 22] have been carried out, adopting different primary cosmic-ray and hadronic interaction models, in order to predict the conventional component, which results from the decay of pions, kaons, and unflavored mesons, and the prompt component, which comes from charmed meson decays – even though the latter contribution is considerably more uncertain and, moreover, undetected at present. Before proceeding, note that atmospheric neutrinos are occasionally accompanied by observable muons, that are absent in the case of cosmic neutrinos instead; the observation of these coincident atmospheric muons allows us to test the reliability of the predictions for the atmospheric neutrinos, and to veto atmospheric events when looking for astrophysical ones.

There are a few good reasons to distinguish these two components: their spectrum, angular distribution, and flavor composition are different, and thus experimentally distinguishable – in principle.

The conventional component of the atmospheric spectrum is mainly produced in the decay of pions and kaons ($\tau_{\text{rest}}\sim${10}^{-8}\text{,}\mathrm{s}), which are produced when cosmic rays collide with the atmosphere; since these mesons have time to interact and lose energy in the atmosphere before decaying, the spectrum of conventional atmospheric neutrinos is expected to be softer than that of the parent cosmic rays of a factor $E^{-1}$, i.e. $\propto E^{-3.7}$. Prompt neutrinos are also produced in hadronic collisions of cosmic rays on the atmosphere as a product of the subsequent immediate decay of charmed mesons ($\tau_{\text{rest}}\sim${10}^{-12}\text{\,}\mathrm{s}), so that their spectrum is expected to closely reflect that of the parent cosmic rays $\propto E^{-2.7}$.

The bulk of atmospheric neutrinos is due to pions and kaons; the larger the energies of the muon produced in the pion/kaon decay, the smaller the probability it decays, so that the flavor ratio of such component shifts from a $(1:2:0)$ for energies lower than $10\text{\,}\mathrm{G}\mathrm{e}\mathrm{V}$ to $(0:1:0)$ at higher energies. The conventional component is also anisotropic: the thicker the layer of the atmosphere traversed (with maximum at an azimuthal angle $\theta=\pi/2$, where $\theta=0$ refers to the upward direction and $\theta=\pi$ to the downward one), the larger the amount of possible targets for cosmic rays.

The prompt component is expected, on the contrary, to be isotropic and to follow a flavor ratio of about $(1:1:0.1)$, with the same energy distribution for all flavors333Expectations based on MCEq [23] and the branching ratios of the relevant mesons [24].. This means that the prompt contribution is the dominant one in the $\nu_{e}$ atmospheric spectrum already at few tens TeV, as we will later show, while, in the case of $\nu_{\mu}$, it reaches the level of the conventional one only at larger-than-PeV energies.

Following [25], we compute the conventional and prompt atmospheric neutrino flux with MCEq [23], adopting the most recent version of SYBILL, the 2.3c release [26], the NRLMSISE-00 [27] model of the atmosphere, and a primary CR flux defined as follows:

• •

only protons and helium nuclei are considered, because they are the most abundant elemental species, and because nuclei of mass number $A$ and energy $E$ produce neutrinos with average energy $E/(20A)$ when colliding with other nuclei;

• •

the most important part of the (Galactic) CR flux is a power-law fitted to the AMS-02 [10, 11] low energy (from 100 GeV to 10 TeV) data, and its knee is assumed to be either an “exponential-square” cutoff or a change of slope.

[TABLE]

[TABLE]

where $\alpha=2-\delta$, $\delta=1/3\equiv$ slope of the diffusion coefficient, i.e. $D(E)=D_{0}(E/E_{0})^{\delta}$. The knee is assumed to be rigidity-dependent, and its position is obtained by fitting the overall primary shape of eq. (2.4) to the KASCADE-Grande data [12];

• •

an additional, supposedly extra-galactic444In [12] this is the component which onsets after the ankle at $E\approx${10}^{17}\text{,}\mathrm{eV}$$. proton component is added to the fitting spectrum

[TABLE]

with $N_{\text{eg}}$ as the only free parameter.

[TABLE]

where the index $k$ identifies the hypothesis on the knee of cosmic rays: its value is either $k={\tt exp2}$-cut for an exponential-square cutoff or $k={\tt delta}$-slope for a change of slope.

The errors on the atmospheric neutrino flux are given by the uncertainty on the shape of the knee and by that on the normalisation and slope resulting from the fits. In table 1 we show the results of our fits.

2.2 Cosmic neutrinos

The predictions for cosmic neutrinos are much more uncertain and require more discussion. Depending upon the adopted model, the resulting neutrino spectra vary greatly in shape and normalization. One of the few stable expectations is that, according to the observed three flavor oscillation phenomena, cosmic electron, muon and tau neutrinos have to be present in similar amounts. The simplest and most popular hypothesis is that cosmic neutrinos are distributed as $E_{\nu}^{-\gamma}$ with $\gamma\sim 2$, at least in the range of energies where they become observable. Surely this case is conducive to observation, however it is important to state clearly what are its motivations, what is its extent, what are its implications; this is the aim of the present discussion.

The physical picture that we have in mind is that cosmic neutrinos are produced in collisions between the accelerated cosmic rays and the gas surrounding the accelerators. In the diffusive shock acceleration (DSA) picture, the cosmic-ray spectrum is a power-law $\propto E^{-2}$; due to the scaling associated to hadronic collisions, also the gamma-ray and neutrino spectra at the source555The neutrino spectrum is supposed to be unaltered even far from the source. will be power laws $\propto E^{-2}$. This setup can be regarded as an extension of what it is commonly supposed for the Galactic cosmic rays, where $E_{\text{max}}$ has probably a lower value than $10-$100\text{,}\mathrm{P}\mathrm{e}\mathrm{V}$$; however, the slope of the injected cosmic rays might be similar or equal to the Galactic one. The abundance of target hadrons points out instead to some dusty environment; this could the site of intense stellar formation, say, starburst and/or star-forming Galaxies.

As a specific instance, we refer to the theoretical model of Loeb and Waxman [28], according to whom the flux of cosmic neutrinos is:

[TABLE]

with

[TABLE]

in the case (as we know today) of an astrophysical neutrino spectrum extending above $100\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$. The slope of this spectrum is very close to the one we would expect in the simple case of DSA.

The assumptions above are consistent with the measurements of through-going muons [2] obtained by IceCube above $200\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ and tested by the HESE dataset [1] in the same energy region, which result in

[TABLE]

Moreover, as we will show in the following, they do not imply any clash with the measurements of the diffuse gamma-ray emission, obtained at lower energies, even if we assume that the cosmic neutrino spectrum extends well below the observed range.

Before proceeding, let us summarize our position. Our goal is to explore the hypothesis that the signal of cosmic neutrinos be produced mainly in $pp$ collisions, a hypothesis consistent with IceCube’s observations at the highest energies. Admittedly, there is no strong theoretical reason as of now to exclude a concurrent, or even dominant, $p\gamma$ origin; thus, we admit that our primary motivation for exploring the $pp$ hypothesis first is based on the Occam razor. On the other hand, the attempts of proceeding beyond this stage of the discussion, by performing accurate modelling of the astrophysical signal, have just begun and do not have yet the compelling and quantitative character that have, say, the modelling of atmospheric neutrinos discussed above.

Assuming $pp$-based sources of neutrinos, we expect pion decay as the main mechanism of neutrino production, thus resulting in a $(\nu_{e}:\nu_{\mu}:\nu_{\tau})\simeq(1:2:0)$ flavor ratio at the source; due to the impact of standard three-flavor neutrino oscillations, we expect a flavor ratio at Earth of about $(1:1:1)$.

In order to have much more precise predictions, we define our phenomenological muon neutrino flux by combining the through-going muon neutrino flux as observed in [2] and that of Loeb and Waxman of eq. (2.5). We do this in the following way:

we define the combined $\gamma_{\text{astro}}$ and $\Phi_{\text{astro}}$ - as well as their errors - as the weighted average of those from [28] and [2]; notice that, due to the large error on the normalization from [28], the weighted average on the normalization is very close to that of [2], while $\gamma_{\text{astro}}^{\text{best}}=2.14\pm 0.08$. The resulting combined muon neutrino flux is thus:

[TABLE] 2. 2.

we reproduce, using a Gaussian likelihood, the 68% CL contour in the $\gamma_{\text{astro}}$-$\Phi_{\text{astro}}$ plane in figure 6 of [2] in order to account for the correlation666Shown in figure 3 of the same paper. $\rho\sim 0.6$ between the two parameters:

[TABLE]

where

[TABLE]

and $\sigma_{\Phi}$ and $\sigma_{\gamma}$ are the errors on the flux normalization at $100\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ and the slope respectively, as taken from [2], since the correlation of the parameters can be extracted only from the data; 3. 3.

we define the best fit astrophysical neutrino flux as the flux averaged using the likelihood defined above:

[TABLE]

and its 1$\sigma$ uncertainty as:

[TABLE]

The resulting cosmic muon neutrino flux is the blue band in figure 2. We can obtain the cosmic flux of $\nu_{e}$ and $\nu_{\tau}$ simply by multiplying that of muon neutrinos by $R_{e\mu}$ and $R_{\tau\mu}$

[TABLE]

In fact, the $R_{\ell\mu}$ factors can be calculated in two different approximations. Standard three-flavor oscillations are assumed in both cases, and described by $P_{\ell\ell^{\prime}}$, the matrix of the survival/oscillation probabilities averaged over cosmic distances

[TABLE]

where $U$ is the neutrino mixing matrix. The two approaches are the following:

the “$2:1$ approximation”, in which we define

[TABLE]

where $\xi_{\ell}^{0}$ is the fraction of $\nu_{\ell}+\overline{\nu}_{\ell}$ produced at the source. We compute $R_{\ell\ell^{\prime}}$ assuming the commonly accepted flavor ratio777Known as the pion decay scenario. Choosing a generic flavor ratio with $\xi_{\mu}^{0}=x$, $\xi_{e}^{0}=1-x$, $\xi_{\tau}^{0}=0$, sampling $x$ uniformly between 0 and 1, we would have $R_{e\mu}=0.78^{+0.57}_{-0.07}$ and $R_{\tau\mu}=1.00^{+0.05}_{-0.15}$. at production of $\xi_{\mu}^{0}=2/3$ and $\xi_{\tau}^{0}=0$ and sampling $P_{\ell\ell^{\prime}}$ according to their distributions (see Appendix A) in the case of normal hierarchy of the neutrino masses. We obtain:

[TABLE] 2. 2.

the “kernel approach”, which relies on a more accurate and physically more comprehensive procedure to fully embrace the consequences of the hadronic production mechanism, namely the strict relationship between gamma rays and neutrinos. As shown in [29], and updated in Appendix A with the use of the most recent oscillation parameters from [30], the gamma-ray flux at the production site and the cosmic neutrino flux are linked by the following relation:

[TABLE]

where $\tilde{K}_{\nu}$ is a kernel which accounts also for $\nu$ oscillations (see [29] and Appendix B for more on this). Equation  (2.9) can be rewritten as

[TABLE]

and it can be easily proven that, assuming a power-law gamma-ray flux $\propto E^{-\gamma}$:

[TABLE]

where

[TABLE]

$R_{\ell\ell^{\prime}}$ as defined in eq. (2.11) depend negligibly888They change of less than $0.1\%$ in the $3\sigma$ range around $\gamma_{\text{astro}}^{\text{best}}=2.14$. on $\gamma$, and this time we obtained:

[TABLE]

We conservatively estimated the impact a 20% variation of the non-oscillated electron (anti-)neutrino kernel variation, which could be due to systematic errors, as neglecting the $K_{L}$ and $K_{S}$ contributions [29]: we found out that $R_{e\mu}$ varies of 6%, while $R_{\tau\mu}$ of 1%, which are not important considering the $30\%$ uncertainty on the astrophysical muon neutrino flux normalization. Moreover, if the electron neutrino and muon neutrino kernels were subject to the same systematic error, the factors $R_{\ell\ell^{\prime}}$ would not change.

Both procedures yield very similar tau neutrino fluxes, while the electron neutrino flux is larger in the kernel approach than in the $2:1$ approximation. This is due to the fact that in the kernel formalism the presence of neutrinos due to kaons is accounted for, the decay of which results in an electron-neutrino richer flux.

We chose to multiply the flux of muon neutrinos by $R_{\ell\ell^{\prime}}$ computed with the kernel approach to consider hadronic collisions as the neutrino production mechanism and include the effect of charged kaons, having thus a clear and precise physical picture in mind. Note that the error on the astrophysical $\nu_{\mu}$ flux normalization is 30%, much larger that those on $R_{e\mu}$ and $R_{\tau\mu}$, so that the relative error on the cosmic electron and tau neutrino fluxes will be 30% as well.

3 Results and discussion

3.1 The components of the neutrino spectra

In figure 2 we show the expectations for the various components of the muon neutrino flux, as well as the corresponding measurements by the IceCube Collaboration. We assume isotropy, but note that the IC-59 points and the HESE ones refer to two different kinds of events, namely $\nu_{\mu}$-induced tracks from the Northern sky and all-flavor High-Energy Starting Events from the whole sky respectively.

In figure 3 we show the prediction for the cosmic electron neutrino flux and the expectations for the other components of the electron neutrino spectrum, alongside with the relevant measurement by the IceCube Collaboration.

No tau neutrino data is available, also due to the fact that the atmospheric tau neutrino flux consists only in the prompt component, thus we do not show the corresponding plot as it would convey no useful information.

From figure 2 and 3 a few noteworthy features are noticeable:

the conventional atmospheric expectations obtained with MCEq, featuring SYBILL-2.3c and our primary cosmic ray spectrum eq. (2.4), agree well with the measurements from IceCube; 2. 2.

the region where the atmospheric and cosmic components of the muon neutrino spectrum cross is around $$250\text{,}\mathrm{T}\mathrm{e}\mathrm{V}$\simeq E_{\text{knee}}/20$; 3. 3.

the prompt component is always subdominant in the $\nu_{\mu}$ spectrum, so that within our model it is not surprising that no significant evidence for prompt neutrinos has been found in the through-going muons analysis; 4. 4.

the conventional atmospheric component is a factor $\sim 30$ less important for electron neutrinos compared to muon neutrinos, so that the prompt component sizably contributes to the overall flux of $\nu_{e}$ for $E_{\nu}\geq$10\text{,}\mathrm{T}\mathrm{e}\mathrm{V}$$.

3.2 Is it possible to extract the prompt neutrino signal?

Our results suggest that the best chance to detect the prompt flux of neutrinos is from electron neutrino atmospheric data for $E\gtrsim$1\text{,}\mathrm{T}\mathrm{e}\mathrm{V}$$, assuming the possibility to somehow discriminate the flavor of the events.

The cascade (or shower) event topology is the most interesting in this regard; it is one of the two kinds of events comprised in the HESE dataset, the remainder being tracks. Track-like events are produced by charged-current muon neutrino interactions, while cascades are produced in all other possible cases. It follows that the cascade sample is the one with the smallest relative contribution of muon neutrinos, which, as already seen in figure 2, are very prompt-neutrino poor. Another reason to focus on the cascade subset of the HESE dataset is that the track-like subset is compatible with being due to background events (atmospheric muons and atmospheric muon neutrinos) only [8, 31].

The contributors to the cascade dataset are:

atmospheric muons; 2. 2.

conventional atmospheric $\nu_{e}$ and $\nu_{\mu}$; 3. 3.

prompt atmospheric $\nu_{e}$, $\nu_{\mu}$ and $\nu_{\tau}$; 4. 4.

astrophysical $\nu_{e}$, $\nu_{\mu}$ and $\nu_{\tau}$.

Therefore, provided that the contamination due to muons, conventional muon neutrinos and all-flavor astrophysical neutrinos can be subtracted, the cascade sample offers us the chance to detect prompt neutrinos.

In order to quantitatively test our hypothesis we use the effective areas for cascade-like events given in figure 1 of [31] to estimate the yearly rate of cascade-like events due to neutrinos with larger-than-TeV energy. In figure 4 we show the parent distribution of cascade events, dividing them by flavor and component, in order to show the energy ranges which contribute the most to the events of table 2. The yearly rates are computed according to:

[TABLE]

where $\mathcal{A}_{\nu_{\ell}}(E)$ is the effective area for the detection of cascade-like events induced by $\nu_{\ell}$ and $\overline{\nu}_{\ell}$, and 1\text{\,}\mathrm{y}\mathrm{e}\mathrm{a}\mathrm{r}$=$\pi\text{\times}{10}^{7}\text{\,}\mathrm{s}.

The expected rate due to prompt neutrinos is small, less than 10% of the conventional contribution, which is not encouraging for the search of prompt neutrinos. However, the conventional contribution could be somewhat reduced by searching for cascades coming from below and/or by using a higher energy threshold, so as to exclude most of the conventional events.

3.3 Is there a spectral anomaly?

In [7, 8, 32] the spectral difference between the HESE and the through-going muons spectra is labelled as an anomaly, and it has been already addressed in [33, 34]. In fact, if we compare the astrophysical neutrino fits to the cascade and starting track samples, which together constitute the HESE dataset, from [31] to that from the through-going muons analysis [2] there is an evident difference, as can be appreciated from table 3 and from figure 10 of [31].

Since, as said before, the starting track sample is compatible to be due to background events only, a more accurate and interesting comparison is between the cascade and through-going muons analyses, which, however, still results in a quite evident spectral difference.

Notice that these two datasets give “complementary” indications, as cascades come from the whole sky and are due to all flavors of neutrinos, with likely a preference for electron and tau neutrinos at high energies, while through-going muons are due to muon neutrinos coming from the Northern sky.

From table 2 and figure 4 a very interesting feature emerges: the number of prompt and cosmic signals in the cascade dataset with $E_{\nu}^{\text{th}}\simeq$1\text{,}\mathrm{T}\mathrm{e}\mathrm{V}$$ are very similar to each other, both in the expected rate of events and in the range of energy in which they contribute. Taking into account also that prompt and cosmic neutrinos are expected to be isotropically distributed in the sky, it appears then difficult999It is possible when the cascade accompanying prompt neutrinos can be tagged. to disentangle the prompt component from the astrophysical one between $1\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ and $100\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ in the cascade (and thus, HESE) dataset. This is consistent with the idea that the sum of the ($\sim E^{-2.7}$) prompt and the ($\sim E^{-2.13}$) astrophysical components could produce the $\sim E^{-2.62}$ spectrum obtained in the cascade analysis. This is demonstrated in figure 5, where we show

• •

the best-fit cosmic neutrino spectrum from the HESE analysis [35]

• •

the best-fit cosmic neutrino spectrum from the cascade analysis [31]

• •

the best-fit cosmic neutrino spectrum from the through-going muons analysis [2]

• •

the sum of the prompt and cosmic components as computed with the model defined in section 2.

We show only the $\nu_{e}$ flavor contribution in figure 5 because it is the most relevant for the cascade dataset due to astrophysical neutrinos, as seen in table 2. While the resemblance of the spectral shape due to the sum of cosmic and prompt neutrino fluxes and the astrophysical best fits from [31, 35] is not perfect, the tension between the analyses below $10\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ seems somewhat alleviated. The spectral shape obviously does not change summing over the flavors, so that this result holds true, but the sum of our prompt and cosmic components would be slightly smaller than three times the best fits from [31, 35]. From this figure, it is evident that:

• •

the spectra resulting from the through-going muons and the HESE (and cascade) analyses are not compatible at low energy, which gives rise to the “spectral anomaly” of the cosmic neutrino spectrum;

• •

the theoretical expectations for the sum of prompt and cosmic neutrinos, instead, agrees within $1\sigma$ with the best-fit cosmic neutrino flux from the cascade analysis.

We conclude that the cause of the alleged spectral anomaly can be attributed to two factors;

a prompt component does contribute to the cascade dataset in the low-energy region, $\lesssim$100\text{,}\mathrm{T}\mathrm{e}\mathrm{V}$$. 2. 2.

a part of the HESE dataset is subject to background contamination due to tracks especially at the lowest energies $\lesssim$10\text{,}\mathrm{T}\mathrm{e}\mathrm{V}$$.

At this point of the discussion, it is useful to bear in mind a couple of important considerations:

1. (i)

the effectiveness of a veto system, based on the presence of muons accompanying the events with a contained vertex [36], is better in the energy range relevant for the search for a cosmic neutrino signal - namely, above several tens of TeV - rather than in the region of lower energies, which is most relevant for the search for prompt neutrinos instead; 2. (ii)

the same analysis that has obtained the cascade dataset [31] has been able to extract also a sample that is highly enriched in muons instead. Its power law description requires a slope of $2.43^{+0.28}_{-0.30}$, whose error is 3-4 times larger than for the cascade dataset and therefore is much less informative101010This is not surprising since the conventional atmospheric component, to be subtracted, is much larger for muon neutrinos - see figure 1 - and therefore, it is harder to identify new components in the track dataset..

We performed nonetheless the same exercise considering muon neutrinos, i.e. comparing the starting tracks best fit and the sum of prompt and astrophysical muon neutrino flux computed in this work. As can be understood from figure 6, no information can be extracted due to the very large uncertainties on the starting tracks sample, which could be due to difficulties in excluding atmospheric contamination.

There are other (non-exclusive) explanations of the low-energy discrepancy between the cosmic neutrino spectrum as resulting from HESE analysis and from the through-going muons analysis: a part of the low-energy soft spectrum of HESE could be due to neutrinos from the Galactic plane, and in this case one would expect a peculiar angular distribution [37, 38, 39, 40, 41] (which to date is not seen [42]); a priori, there could be other sources of extraterrestrial neutrinos which could cause such effect (but this would be at odds with the null search of prompt neutrinos [2]). Note that our proposal, concerning the role of prompt neutrinos in the cascade dataset, does not require the inclusion of hypothetical physical ingredients, and in this sense can be considered minimal.

4 Conclusions

In this work we obtained predictions for the components of the high-energy ($\gtrsim$100\text{,}\mathrm{G}\mathrm{e}\mathrm{V}$$) neutrino spectrum.

We defined a primary cosmic-ray flux, which has been tuned to the AMS-02 [10, 11] and KASCADE-Grande [12] proton and helium data, to compute the conventional and prompt components of the atmospheric neutrino spectrum.

Regarding cosmic neutrinos, we embraced hadronic collisions in dense sources, like star-forming/starburst Galaxies, as the production mechanism of astrophysical neutrinos, and combined the flux from [28] and the best fit from [2] to obtain a precise phenomenological muon neutrino flux. We accurately described neutrino oscillations updating the formalism introduced in [29] (Appendix A) with the recent oscillation parameters from [30], and predicted the astrophysical electron and tau neutrino fluxes from the muon neutrino one.

We computed the theoretical 1$\sigma$ bands for all these fluxes and compared them to the available data, obtaining a satisfactory agreement, as can be seen in figures 2 and 3.

Our results are compatible with the fact that no prompt component has been seen in the through-going muon analysis, which regards muon neutrinos only. We thus proceeded to investigate the feasibility of a measurement of prompt neutrinos from the most relevant dataset, i.e. the muon-neutrino depleted cascade dataset. The yearly rate of cascades in IceCube due to all flavors and components of the neutrino spectrum has been computed (table 2). In this framework, the extraction of the prompt-neutrino signal from the cascade dataset seems difficult, but possible.

We also addressed the so called “spectral anomaly” [7, 8] of the low-energy ($\lesssim$100\text{,}\mathrm{T}\mathrm{e}\mathrm{V}$$) part of the cosmic neutrino spectrum arising from the different IceCube analyses [2, 31, 35]. Upon visual comparison (figure 5), it is apparent how the spectral anomaly could be attributed to the presence of prompt neutrinos, which would be the reason of the discrepancy between the $\sim E^{-2.13}$ spectrum from the through-going muon analysis and the $\sim E^{-2.48}$ one from the cascade analysis.

We believe that great care should be taken when trying to disentangling atmospheric, and especially prompt, neutrinos from astrophysical ones between $1\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$ and $200\text{\,}\mathrm{T}\mathrm{e}\mathrm{V}$, and that a global analysis, adopting the same theoretical models for all datasets would be the best way to extract informations from the data.

Appendix A Neutrino oscillation/survival probabilities and their distributions

In order to account for neutrino oscillations, we followed the methodology of [43], which is based on the observation that the oscillation matrix $P_{\ell\ell^{\prime}}$ is symmetric under $\ell\leftrightarrow\ell^{\prime}$:

[TABLE]

and thus they worked out a “natural” parametrization with just three independent parameters, named $P_{0}$, $P_{1}$ and $P_{2}$. $P_{\ell\ell^{\prime}}$ is then parametrized as follows:

[TABLE]

with

[TABLE]

and

[TABLE]

As for the oscillation parameters, we considered their distributions as reported in [30] for both normal ordering (NO) and inverted ordering (IO). We performed Monte Carlo extractions of these parameters to compute the mode and the 68% confidence intervals for the parameters $P_{0}$, $P_{1}$, $P_{2}$, which are reported in table 4.

The errors are reported with the symbol $\pm$ whenever the distribution of the relevant parameter is gaussian-like, while they are reported asymmetrically (even though they could be equal) for non-gaussian distributions. All these considerations apply also to table 5, in which we report the values and confidence intervals for the oscillation/survival probabilities $P_{\ell\ell^{\prime}}$, which we obtained from the same Monte Carlo extractions used to compute the distribution of $P_{0}$, $P_{1}$, and $P_{2}$.

In figure 7 we show the distributions of $P_{i}$ and $P_{\ell\ell^{\prime}}$. We can conclude that there is compatibility with the results obtained in [43].

Appendix B The kernels linking neutrinos to photons

Due to the change over time of the oscillation parameters measurements with time, we find it useful to provide the reader with a resume of the kernels to obtain neutrinos from photons, as in eq. (2.10); in fact the oscillated kernels, which we denote with a $\sim$, are easily obtained from the non-oscillated ones like so:

[TABLE]

The generic form of the non-oscillated kernels is the following:

[TABLE]

and their parameters are listed in table 6, while $r_{i}=(m_{\mu}/m_{i})^{2}$.

Acknowledgments

The authors are thankful to Andrea Palladino, Maria Vittoria Garzelli, Paolo Lipari, Spencer Klein and the anonymous referee for valuable discussions and feedback on the text. This work was partially supported by the research grant number 2017W4HA7S “NAT-NET: Neutrino and Astroparticle Theory Network” under the program PRIN 2017 funded by the Italian Ministero dell’Istruzione, dell’Università e della Ricerca (MIUR).