High-energy atmospheric neutrinos
Seminar at Chiba University, Japan, 28th July 2017

Anatoli Fedynitch
DESY ZeuthenAtmospheric leptons> For high precision calculations allphenomena need accurate modeling> Not so well known “ingredients”:§ Cosmic ray spectrum and composition§ Hadronic interactions§ Atmosphere (dynamic, depends on use case)§ (Rare) decays (solved)> Energy range MeV – EeV!

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 2

Production channels

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 3

General shape of high energy lepton spectra> 30 km
Zenith angle changes competition
between decay and (re-)interaction
Depth X(g/cm2)

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 4

Transport equations (hadronic cascade equations)
System of non-linear PDE for each particle species h (~62 x #E-bins) :cosmic ray physics
Interactions with air

atmospheric physics

DecaysEnergy losses (radiative)
Re-injection frominteractions
Re-injection fromdecaysparticle physics

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 5Inclusive secondary particle distributions (as in fixed-target experiment)
Measurements at low energy,
no clear extrapolation recipe

No kaon data forproton-carbon!

F. Riehn, R. Engel, A. Fedynitch, T. K. Gaisser and T. Stanev.
arXiv:1510.00568, PoS ICRC 2015, 1129 (2015)
SIBYLL 2.3c: Riehn et al. PoS(ICRC2017)301

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 6

MCEq: open-source Python code

> Simultaneous solution of up to 8000 kinetic equations> Energy range 1 (30) GeV – 1011 GeV> All models included> High optimization: multi-core, GPU, … (BLAS, MKL,CUDA) (~milli-seconds!)> MIT licensed @ https://github.com/afedynitch/MCEq

CORSIKA: A. Fedynitch, J. Becker Tjus and P. Desiati, PRD 2012
MCEq: A. Fedynitch, R. Engel, T. K. Gaisser, F. Riehn and S. Todor. PoS ICRC 2015, 1129
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 7

Relevant hadrons for neutrino production
Muon neutrinos

pion decay kaon decay

charm decay

Electron neutrinos

muon decaykaon decay

charm decay

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 8

Relevant particles for muon production
Muons..from unflavored mesons

103

Fµ (E/GeV)3 [cm2 s sr GeV] 1101
100101

K
D± , D0Ds
Lc

Fµ (E/GeV)3 (cm2 s sr GeV) 1

total
total conv.
total prompt
other prompt
other conv.
p102

10 2

unflavored

10 2
10 3
10 410h

fh0w

r0
J/Y

310 4

10 5

10 610 510 6102

103

104

105

106

107

108

Eµ [GeV]pions + kaonscharm + unflv.

unflavored109

10 7102103104105
106
Eµ [GeV]

These mesons decay into muonpairs with small branching ratios
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 9107

108109Relation between lepton and cosmic ray energy

Energy range covered byparticle accelerators.

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 10

If energy is not a problem…
Beam pipe!

Kinematic variables

For atmospheric leptons

Theory control
The LHC is great!

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 11

Secondary hadron (pion/kaon/charm) energy fractions
xlab regions contributing to inclusive
muon neutrino flux

> Atmospheric leptons are sensitive mostlyto xlab ~ 0.2 and above> Reason: steepness of primary CRspectrum

NA49 pC @ 158 GeV:
x2.5dN/dxL

CR Knee, change ofspectral index
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 12

Role of secondary interactions in cascade
Inclusive leptons

It’s not 90% from first interaction(as often assumed)

100 PeV p air shower

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 13

Comic ray spectrum at the top of the atmosphere> None of the shown models is a real fit(with errors and covariance matrix)> GST-X and HXa are quite extreme
assumptions for UHECR> No error estimates!
HXa:
T. K. Gaisser,
Astropart. Phys. 35 (2012) 801–806
GH:
T. K. Gaisser and M. Honda,
Ann. Rev. Nucl. Part. Sci. 52 (2002) 153–199
GST-X:T. K. Gaisser, T. Stanev, and S. Tilav,
Front. Phys.(Beijing) 8 (2013) 748–758
Update of GH (not shown):
J. Evans, D. Porzio et al., 1612.03219

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 14

New GSF all-nucleon flux (input for atm. neutrino calculations)> Novel global fit to entire CR spectrum: Dembinskiet al., PoS(ICRC2017) 533> Wavy feature at lower energies are due to thehardening of proton and He spectra> Increase of error around 10 TeV because of thegap between direct and indirect exp.
> Higher flux at the knee and harder spectral index
between knee and ankle> Latter effect comes from the lighter
composition at the knee as in other models> Mainly driven by KASKADE Gr. and latest datafrom IceTop and TUNKA

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 15

Hadronic uncertainties: current state of the art> “Uncertainties in atmospheric neutrino fluxes”, G. D. Barr, S.
Robbins, T. K. Gaisser, and T. Stanev, Phys. Rev. D 74,094009 (2006)> Cut phase-space in regions/slices in Elab and xlab and assign
uncertainty to each slice (uncorrelated)> Draw-back 1: Uncertainty assigned by hand and judged only
from availability of experimental data (not how well TARGETdescribed it)> Draw-back 2: The “central value” is assumed to be TARGET.
Scheme doesn’t tell anything about ”best estimate”.

Up to 10 TeV neutrino energy, no muons.
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 16

Implementation of the “Barr scheme” in MCEq
The regions> Compute partial derivatives wrt. phase-space
regions, i.e.
> No correlations between phase-space regions(as in Barr et al.) or add. correlations
Elements of Jacobian (numerical)

CH_A
10 %Charm+ 25 %70CH_B% + max.
Energydependence
dependence
EnergyError propagation

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 17

Uncertainties of lepton fluxes

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 18

Uncertainties of lepton ratios

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 19

Hadronic uncertainties for high energy leptons

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 20

Prompt neutrinos from decays of charmed mesons> SIBYLL 2.3c is the only full MC model> Compatible with LHC data and IceCube limit> New CR flux model (GSF) changes situation a bit> Uncertainties from QCD very large and calculations
are compatible> Uncertainties from pQCD calculations are presumably
overestimated, since LHC measurements of charm are
not taken into account

IceCube: Astrophys.J. 833 (2016)
GMS: Garzelli et al., JHEP 1510 (2015) 115
BERSS: Bhattacharya et al. JHEP 2015: 110
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 21

Prompt neutrinos from decays of charmed mesons> SIBYLL 2.3c is the only full MC model> Compatible with LHC data and IceCube limit> New CR flux model (GSF) changes situation a bit> Uncertainties from QCD very large and calculations
are compatible> Uncertainties from pQCD calculations are presumably
overestimated, since LHC measurements of charm are
not taken into account

IceCube: Astrophys.J. 833 (2016)
GMS: Garzelli et al., JHEP 1510 (2015) 115
BERSS: Bhattacharya et al. JHEP 2015: 110
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 22Conclusions> Relevant cosmic ray energies, contributing particles and their phase-space
understood in detail.
> Progress transparent to the community through an open source code, support, etc.
> Cosmic ray uncertainties are quantified to our best knowledge. The numerical model
including the covariance matrix will be open source, as well.
> Uncertainties, at the level of Barr et al. can be reproduced with MCEq andcorresponding tools will be published, asap.
> High-quality atmospheric muon and neutrino data (incl. systematics, error correlation
matrix, etc.) would definitely help

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 23

Hadronic uncertainties below 200 GeV

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 24

Hadronic uncertainties computed with different interaction models

1. If models with known problems are excluded, then errors are~universal = approx. independent on the choice of interaction model
2. “Bracketing” with different models can be a good estimate, sometimes

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 25Leading particles

Many “slow”
particles

Few “fast”
particles

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 26

Main changes in SIBYLL 2.3c are related to scaling and leading Kaons> Problems in SIBYLL 2.3 related to remnant excitation model> No changes to charm in 2.3c
> Small changes for air-shower physics

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 27

Contribution of the leading particle effectModeling of leading particle effect contributes to large uncertainties

Larger effect in SIBYLL
compared to other models
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 28

Charm in SIBYLL“non-perturbative/leading” charm“pQCD/minijet” charm

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 29

Charm in SIBYLLPROSA,
Garzelli et al. 2017

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 30

Muon yield in single air-shower

High uncertainty
for air-showers

Scaling

Ratios of secondary spectra

Large uncertainties among models

High uncertainty
for inclusive fluxes

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 31

A few words on 3D calculations

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 32

Classical solutions
1D MC calculation

Semi-analytical
T. K. Gaisser, R. Engel and E. Resconi,
Cosmic Rays and Particle Physics, 2nd edition

Very slow> 3 – 6 particle
speciesdetail andLimited

precision> De-coupling of equations> Numerous assumptions about spectral
shape, atmospheric model, etc.
> Only asymptotical low- and high-energysolutions

M. Honda et al. 1995, Agrawal et al. 1995, Fedynitch et al. 2012

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 33

MCEq: Matrix Cascade Equations

A. Fedynitch, R. Engel, T. K. Gaisser, F. Riehn and S. Todor
PoS ICRC 2015, 1129 (2015), EPJ Web Conf. 99, 08001 (2015)
and EPJ Web Conf. 116, 11010 (2016). Main paper in prep.

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 34

Competition between decay and interactions> Above critical energy, probability for interaction is higher
than decay probability> Flux from one particular hadron species gets attenuatedCritical energy

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 35

Hadrons in MCEq

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 36

Matrix forms

Decay matrix D

Interaction matrix C

matrices are sparsehigh performance

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 37

Stiffness is tricky

Eigenvalues of matrix equation

Idea: transport particles only when needed

Fastest eigenvalue controls integration step

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 38

A new take on CR flux uncertainties> Combine datasets from direct andindirect observations> Use only “golden” datasets withsystematic errors (incl. energy scale)> Fit whole energy range 10 – 1011 GeV> Fit composition:§ Direct: elements§ Indirect: mass groups

PoS(ICRC2017)533:
Hans Dembinski, AF, Ralph Engel, Tom K. Gaisser,
Felix Riehn, Todor Stanev

Idea: propagate uncertainties into fluxes
fluxes using partial derivatives from MCEq
Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 39

How to create mass groups (indirect) from element fluxes (direct)
Mass group ⟷ Element

Fit of direct observations> Ratio dN/dR roughly equalfor neighboring elements ofa mass group§ The only assumption (at allenergies)> Fit 4 mass groups with
distinct sets of B-splines

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 40

Reproduction of data from air shower based detectors

Dotted: H3a
Dashed: H4a

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 41

GSF: Global Spline Fit> Data-driven representation of thecosmic ray flux
> Full covariance matrix for allparameters> Serves as input for flux
calculations and error propagation

Energy scale

Anatoli Fedynitch | Discussion and seminar, Chiba, Japan | 2017/07/28 | Page 42