Documentation for crflux

Coming soon….

Module documentation

models — models of the cosmic ray flux at the top of Earth’s atmosphere

This module is a collection of models of the high-energy primary cosmic ray flux, found in the literature of the last decades. The base class PrimaryFlux has various methods, which can be employed in lepton flux calculations, cosmic ray physics, radiation physics etc.

The numbering scheme for nuclei is adapted from the Cosmic Ray Air-Shower Monte Carlo CORSIKA. Protons have the ID 14. The composite ID for nuclei follows the formula ID=100 \cdot A + Z, where A is the mass number and Z the charge. With this scheme one can easily obtain charge and mass from the ID and vice-versa (see PrimaryFlux.Z_A()).

The physics underlying each of the models can be found following the references in this documentation.

Note

As always, if you use one of the models for your work, please cite the corresponding publication!

Example:

To generate the plots from below, run:

from crflux.models import test
test()

(Source code)

class models.PrimaryFlux[source]

Base class for primary cosmic ray flux models.

nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
total_flux(E)[source]

Returns total flux of nuclei, the “all-particle-flux”.

Args:
E (float): laboratory energy of particles in GeV
Returns:
(float): particle flux in \Phi_{particles} in (\text{m}^2 \text{s sr GeV})^{-1}
tot_nucleon_flux(E)[source]

Returns total flux of nucleons, the “all-nucleon-flux”.

Args:
E (float): laboratory energy of nucleons in GeV
Returns:
(float): nucleon flux \Phi_{nucleons} in (\text{m}^2 \text{s sr GeV})^{-1}
nucleon_gamma(E, rel_delta=0.01)[source]

Returns differential spectral index of all-nucleon-flux obtained from a numerical derivative.

Args:
E (float): laboratory energy of nucleons in GeV rel_delta (float): range of derivative relative to log10(E)
Returns:
(float): spectral index \gamma of nucleons
nucleus_gamma(E, corsika_id, rel_delta=0.01)[source]

Returns differential spectral index of nuclei obtained from a numerical derivative.

Args:
E (float): laboratory energy of nuclei in GeV corsika_id (int): corsika id of nucleus/mass group rel_delta (float): range of derivative relative to log10(E)
Returns:
(float): spectral index \gamma of nuclei
delta_0(E)[source]

Returns proton excess.

The proton excess is defined as \delta_0 = \frac{\Phi_p - \Phi_n}{\Phi_p + \Phi_n}.

Args:
E (float): laboratory energy of nucleons in GeV
Returns:
(float): proton excess \delta_0
p_and_n_flux(E)[source]

Returns tuple with proton fraction, proton flux and neutron flux.

The proton fraction is defined as \frac{\Phi_p}{\Phi_p + \Phi_n}.

Args:
E (float): laboratory energy of nucleons in GeV
Returns:
(float,float,float): proton fraction, proton flux, neutron flux
lnA(E)[source]

Returns mean logarithmic mass <ln A>/

Args:
E (float): laboratory energy of particles in GeV
Returns:
(float): mean (natural) logarithmic mass
Z_A(corsika_id)[source]

Returns mass number A and charge Z corresponding to corsika_id

Args:
corsika_id (int): corsika id of nucleus/mass group
Returns:
(int,int): (Z,A) tuple
class models.PolyGonato(constdelta=False)[source]
    1. Hoerandel, Astroparticle Physics 19, 193 (2003).
nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
class models.HillasGaisser2012(model='H4a')[source]

Gaisser, T.K., Astroparticle Physics 35, 801 (2012).

Model is based on Hillas ideas and eye-ball fits by T.K. Gaisser. H3a is a 3-‘peters cycle’ 5 mass group model with mixed composition above the ankle. H4a has protons only in the 3. component.

Args:
model (str): can be either H3a or H4a.
nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
class models.H3a_polygonato(model='H3a')[source]

Modified version of Gaisser, T.K., Astroparticle Physics 35, 801 (2012).

Model is based on Hillas ideas and eye-ball fits by T.K. Gaisser. H3a is a 3-‘peters cycle’ 5 mass group model with mixed composition above the ankle. H4a has protons only in the 3. component.

This version has a five component poly-gonato at lower energies and resembles H3/4a at energies above the knee.

Args:
model (str): can be either H3a or H4a.
class models.GaisserStanevTilav(model='3-gen')[source]
    1. Gaisser, T. Stanev, and S. Tilav, arXiv:1303.3565, (2013).
Args:
model (str): 3-gen or 4-gen
Raises:
Exception: if model not properly specified.
nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
class models.CombinedGHandHG(model='H3a')[source]
  1. Fedynitch, J. Becker Tjus, and P. Desiati, Phys. Rev. D 86, 114024 (2012).

In the paper the high energy models were called cHGm for GH+H3a and cHGp for GH+H4a. The names have been change to use the quite unintuitive names H3a and H4a in ongoing literature.

nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
class models.ZatsepinSokolskaya(model='pamela')[source]

The model was first released in V. I. Zatsepin and N. V. Sokolskaya, Astronomy and Astrophysics 458, 1 (2006).

Later, the PAMELA experiment has fitted the parameters of this model to their data in PAMELA Collaboration, O. Adriani et al., Science 332, 69 (2011). Both versions of parameters can be accessed here.

The model does not describe the flux above the knee. Therefore, the highest energies should not exceed 1-10 PeV.

Args:
model (str): ‘default’ for original or ‘pamela’ for PAMELA parameters
nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
class models.GaisserHonda(*args, **kwargs)[source]

5 mass group single power-law model from T.K. Gaisser and M. Honda, Annual Review of Nuclear and Particle Science 52, 153 (2002)

This model was tuned to lower energy baloon data. It fails to describe the flux at and above the knee. A safe range for using it is < 100 TeV/nucleon.

nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
class models.Thunman(*args, **kwargs)[source]

Popular broken power-law flux model.

The parameters of this model are taken from the prompt flux calculation paper by M. Thunman, G. Ingelman, and P. Gondolo, Astroparticle Physics 5, 309 (1996). The model contians only protons with a power-law index of -2.7 below the knee, located at 5 PeV, and -3.0 for energies higher than that.

nucleus_flux(corsika_id, E)[source]

Broken power law spectrum for protons.

class models.SimplePowerlaw27(*args, **kwargs)[source]

Simple E**-2.7 parametrization based on values from Thunman below knee.

nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
class models.GlobalSplineFit(*args, **kwargs)[source]

Data-driven fit of direct and indirect measurements of the cosmic ray flux and composition. Tracks the mass composition using four leading elements (p, He, O, Fe), whose flux is modeled by shaped spline functions. Assumes fixed flux ratios in rigidity for the flux of subleading elements. Covers the whole rigidity range from 10 GV to 10^11 GeV.

nucleus_flux(corsika_id, E)[source]

Returns the flux of nuclei corresponding to the corsika_id at energy E.

Args:
corsika_id (int): see crflux for description. E (float): laboratory energy of nucleus in GeV
Returns:
(float): flux of single nucleus type \Phi_{nucleus} in (\text{m}^2 \text{s sr GeV})^{-1}
p_and_n_flux(E)[source]

Returns tuple with proton fraction, proton flux and neutron flux.

The proton fraction is defined as \frac{\Phi_p}{\Phi_p + \Phi_n}.

Args:
E (float): laboratory energy of nucleons in GeV
Returns:
(float,float,float): proton fraction, proton flux, neutron flux
total_flux(E)[source]

Returns total flux of nuclei, the “all-particle-flux”.

This version of the method does not vectorize the nucleus_flux method.

Args:
E (float): laboratory energy of particles in GeV
Returns:
(float): particle flux in \Phi_{particles} in (\text{m}^2 \text{s sr GeV})^{-1}
lnA(E)[source]

Returns mean logarithmic mass <ln A>/

This version of the method does not vectorize the nucleus_flux method.

Args:
E (float): laboratory energy of particles in GeV
Returns:
(float): mean (natural) logarithmic mass
dump_nucleon_flux_splines(emin=1.0, emax=1000000000000.0, nbins=1000)[source]

Dumps a nucleon flux splines of the full GSF model to a pickled file.

Energy and flux coordiante are interpolated as a natural logarithm of the values.

Args:
emin (float): minimal energy for the spline range emax (float): maximal energy for the spline range nbins (int): number of energy steps for interpolation
class models.GlobalSplineFitBeta(spl_fname=None)[source]

Data-driven fit of direct and indirect measurements of the cosmic ray flux and composition. Tracks the mass composition using four leading elements (p, He, O, Fe), whose flux is modeled by shaped spline functions. Assumes fixed flux ratios in rigidity for the flux of subleading elements. Covers the whole rigidity range from 10 GV to 10^11 GeV.

Tabulated ICRC 2017 version. The full interface will become available, when the model is published. The class picks the most recent spline file in the current directory.

Use for reference: Dembinski et al., PoS ICRC2017 533 https://inspirehep.net/literature/1639832

p_and_n_flux(E)[source]

Returns tuple with proton fraction, proton flux and neutron flux.

The proton fraction is defined as \frac{\Phi_p}{\Phi_p + \Phi_n}.

Args:
E (float): laboratory energy of nucleons in GeV
Returns:
(float,float,float): proton fraction, proton flux, neutron flux
tot_nucleon_flux(E)[source]

Returns total flux of nucleons, the “all-nucleon-flux”.

Args:
E (float): laboratory energy of nucleons in GeV
Returns:
(float): nucleon flux \Phi_{nucleons} in (\text{m}^2 \text{s sr GeV})^{-1}
nucleus_flux(corsika_id, E)[source]

Dummy function, since for particle fluxes are not supported in the spline interface version

Indices and tables