Phenomenology of atmospheric neutrinos

The detection of astrophysical neutrinos, certainly a break-through result, introduced new experimental challenges and fundamental questions about acceleration mechanisms of cosmic rays. On one hand IceCube succeeded in finding an unambiguous proof for the existence of a diffuse astrophysical neutrino flux, on the other hand the precise determination of its spectral index and normalization requires a better knowledge about the atmospheric background at hundreds of TeV and PeV energies. Atmospheric neutrinos in this energy range originate mostly from decays of heavy-flavor mesons, which production in the phase space relevant for prompt leptons is uncertain. Current acceleratorbased experiments are limited by detector acceptance and not so much by the collision energy. This paper recaps phenomenological aspects of atmospheric leptons and calculation methods, linking recent progress in flux predictions with particle physics at colliders, in particular the Large Hadron Collider.


Introduction
Hadronic interactions of high energy cosmic rays (CR) with the atmosphere produce extensive cascades of secondary particles.The electromagnetic component of electrons and -rays dissipates its energy in form of bremsstrahlung and pair production.In a particle cascade, hadrons can either interact, producing a sub-cascade similar to the original one but at lower energy, or decay into leptons or other hadrons.Since leptons interact very weakly with the atmosphere, they reflect the history of the air-shower development and easily reach down to the surface.
A Monte Carlo program seems to be the most reasonable calculation method, because by using random numbers it can generate an infinite number of cascade topologies and also employ detailed physical models of sub-leading processes, such as multiple scattering.The high degree of realism is, on the other hand, a drawback of the Monte Carlo calculation, since the history of every particle has to be traced, and there are billions of particles in a high energy particle shower.Due to the requirements on computing time and storage, one has to deal with a limited number of repetitions constraining the achievable precision.Rare processes or subtle signatures might therefore be smeared out by the statistical error.The flux of long-lived mesons and the resulting conventional lepton flux component attenuate quickly at very high energy due to re-interactions with air nuclei.Biasing techniques in a e-mail: anatoli.fedynitch@cern.chThis is an Open Access article distributed under the terms of the Creative Commons Attribution License 4.0, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
EPJ Web of Conferences simulation, as implemented in FLUKA [1,2] or in upcoming Corsika [3] versions, overcome this problem to some extent.
Additional challenges stem from the inclusion of very rare processes such as the decay of heavy flavor, charm or bottom, hadrons inside air-showers.Heavy mesons, mostly D ± and D 0 , decay immediately and produce prompt muons and neutrinos.This component does not suffer from reinteractions with air up to energies of several tens of PeV.But a drop in calculation efficiency comes from difference in inclusive production cross-sections.The ratio of production cross-sections in a p-air collision at 100 PeV and high x = E secondary E CR ∼ 0.5 is ∼ ± / D ± ∼ 10 −4 − 10 −5 , i.e. 10 4 simulated showers yield on average a few prompt leptons.At this point one can start introducing additional approximations or optimizations to further reduce computational time in exchange for precision, loosing at some point the justification to run complex Monte Carlo codes at all.
So, why should we be interested in a very efficient calculation?First of all, because the models used as ingredients in an atmospheric flux calculation carry sizable uncertainties coming from models of primary cosmic ray spectrum and composition, from the seasonal and geographic variability of the atmosphere and from the modeling of hadronic interactions.By having the possibility to easily re-run calculations with varying model assumptions we can propagate the effect of each model uncertainty on the input side to the flux prediction.
In the next section a numerical method is introduced, which competes with full Monte Carlo calculations in terms of precision, but requires only a tiny fraction of computational time.

Matrix cascade equations
The flux of particles in air is described by the coupled cascade integro-differential equations Particle interactions and decays are the two competing processes and it is convenient to write the equation as a function of the slant depth X [g/cm 2 ].Continuous energy losses arising from the passage of charged particles through air can be neglected at energies above several tens of GeV.Index h determines the type of particle, like , K, , D ± , , ν, etc.The sums over all cascade members l define the couplings between the evolution equations for each particle type.The inclusive particle production spectra or redistribution functions dN dE describe the spectrum of secondaries created in particle interactions or decays.The process rates are defined in terms of interaction and decay lengths in units of X.More details about the equation and solutions can be found in [4,5].
In a typical semi-analytical approach one tries to decouple the system by neglecting rare processes, for example, the production of nucleons in secondary pion interactions or the production of K − in interactions of K + .Often these approximations are valid under certain restrictions, such as under the assumption of steep power-law primary spectra or Feynman-scaling of particle production spectra.
Very Large Volume Neutrino Telescope (VLVnT-2015) One can criticize that additional uncertainties due to the calculation method are introduced (and need quantification), on top of those from the physical models themselves.
An approach aimed to avoid these calculation uncertainties is based on numerical integration of the coupled system by discretization on an energy grid.Eq. (1) becomes where the interaction coefficients absorb the finite interval size and the inclusive production cross-sections.The energy grid extends from 50 to 10 10 GeV using 8-9 bins per decade.Analogously to Eq. ( 3) the inclusive decay coefficients are obtained by summing over all spectra of exclusive branching channels i that contain the secondary particle h.By combining the solutions for each particle type h and each energy bin E i in a vector the coupled system of equations can be expressed in matrix notation The interaction and decay lengths are arranged in diagonal rate (∼ 1/ ) matrices , where the air density dependence has been factorized out.The interaction or decay coefficients from Eqs. ( 3) and ( 4) are arranged in sparse matrices C and D.
By numerically integrating Eq. ( 6), the coupled cascade equations are evolved collectively for all contributing particle types.Numerical issues arise when short-lived hadrons are included, such as D mesons or C baryons, since these particles decay on a much shorter timescale compared to conventional hadrons.By introducing a special type of steady-state approximation, the so-called resonance approximation, these stiffness problems can be resolved.The full system of the order 6000 × 6000 can be solved on an ordinary PC in a few seconds.More details about this approach and general results can be found in [6,7].The code "MCEq" is available online1 .

Connection with LHC
The Large Hadron Collider (LHC) at CERN with its center of mass energies between 900 GeV and 13 TeV covers an energy regime which is very relevant for inclusive lepton fluxes.Figure 1 shows that the LHC reaches cosmic ray interaction energies at which 1 PeV to 10 PeV atmospheric neutrinos are produced.However, the data do not constrain significantly atmospheric flux calculations.One of the reasons is that the LHC is running in pp, pP b or P bP b mode.Lead's high mass produces an extreme energy density emphasizing collective effects of Quark-Gluon-Plasma rather than the role of identified high-x mesons which are of higher importance for cosmic ray applications.These data is therefore of limited use for cosmic ray interaction model developers who aim to constrain nuclear effects up to the iron mass.This situation can be significantly improved by running the LHC with carbon or nitrogen nuclei.
A more serious constraint is the phase space coverage of current LHC experiments.Except a few dedicated roman pots or zero-degree calorimeters, none of the general purpose detectors covers the xrange relevant for inclusive fluxes.Typical Feynman-x F of secondary particles rarely exceed 0.05 before leaving the acceptance region.Using the matrix code it is possible to study the contribution of ranges in x ≈ x F (at high energies) to the total inclusive flux.Bands in Fig. 2 illustrate the contribution of different kinematical regions to the inclusive lepton flux by using the ability of the numerical integrator to apply cuts on ranges in x.Less than 10% of conventional muons and muon neutrinos originate from secondaries created at an x value < 0.09, giving an idea to what extent an ordinary LHC detector with particle identification capabilities would experimentally constrain the flux prediction.Muon neutrinos exhibit a pronounced shift towards lower x above 100 TeV, where prompt neutrinos dominate.They are sensitive to the production properties of charmed mesons which are predicted to be distributed more centrally in Sibyll-2.3[8,9].These results are model dependent and by exchanging the interaction model some shifts in the importance of x-ranges can be observed.The central result is that the interaction models are weakly constrained by LHC observations in the phase space relevant for atmospheric leptons.

Comments on prompt neutrino calculations
The prompt component is one of the dominant backgrounds in the determination of the spectral index of the astrophysical neutrinos with IceCube [10].For an unambiguous identification of an up-going astrophysical flux a precise determination of the contribution from prompt neutrinos is even of crucial importance.So far, IceCube has not claimed any a signal from prompt leptons.However, we should be certain that this flux exists, since heavy flavor mesons are created in all particle collisions of sufficient energy.
We have estimated this flux using a preliminary version of the sibyll 2.3 Monte Carlo generator, which has been cross-checked against LHC data, and a new, preliminary, version of dpmjet-III [11].In the left panel of Fig. 3 is shown an updated version of our plot from [6].It compares more popular calculations ERS [12], TIG [13] and MRS [14] with the more recent BERSS [15] and GMS [16] calculations.The band around the GMS line represents uncertainties from the inputs of the nextto-leading order (NLO) perturbative calculation, i.e. renormalization and factorization scale, parton distribution function, charm mass and parton shower matching.A similar calculation [17], aiming for extraction of state-of-the-art knowledge about forward charm production in accordance with LHCb results, arrives at comparable errors.On one hand, the high uncertainty is expected given the limited acceptance of LHC experiments with respect to the phase space relevant for atmospheric leptons.On the other hand, none of the listed calculations, although carried out with different methods and models, contradicts the accelerator-based estimation.This is encouraging, since a measurement of the flux normalization by neutrino telescopes, like IceCube or ANTARES, will give insight into small-x physics inaccessible to the current collider generation.11010-p.5

EPJ Web of Conferences
However, the flux variation due to the cosmic-ray composition is of similar magnitude, as it can be seen by comparing with the right panel of Fig. 3. H3a [18] and GST3 [19] are mixed composition scenarios, while H4a and GST4 are model versions with only protons at energies above the ankle of cosmic ray energy spectrum.With more data and higher precision from the upgraded versions of the Pierre Auger Observatory and Telescope Array it should be possible to further constrain this particular uncertainty.
I am grateful to Ralph Engel, Felix Riehn, Thomas K. Gaisser and Todor Stanev for inspiring discussions.

Figure 1 .
Figure 1.Energy spectrum for priamry cosmic rays that lead to the production of mono-energetic secondary leptons.Line styles represent different leptons and colors the lepton energy.

Figure 2 .
Figure 2. Contribution of (Feynman-)x ranges in the calculation of the inclusive atmospheric lepton flux for muons (left) and electron neutrinos (right).

Figure 3 .
Figure 3. Prompt muon or electron neutrino flux calculated for a simple broken power-law proton spectrum vs. other calculations described in the text (left).Muon neutrino flux calculated using different composition scenarios of ultra-high energy CR (right).Dashed lines represent the conventional flux.