Auger

Contains different functions needed for analysing data of the Pierre Auger Observatory:

It makes use of the following reference

  1. Manlio De Domenico et al., JCAP07(2013)050, doi:10.1088/1475-7516/2013/07/050

    1. Adeyemi and M.O. Ojo, Kragujevac J. Math. 25 (2003) 19-29

  3. JCAP 1302 (2013) 026, Interpretation of the Depths of Maximum of Extensive Air Showers Measured by the Pierre Auger Observatory, DOI:10.1088/1475-7516/2013/02/026

  4. Auger ICRC’15, data file received from Ines Valino on 2015-07-30

  5. GAP-2014-083, Update of the parameterisations given in “Interpretation of the Depths …” in the energy range 10^17 - 10^20 eV

  6. Long Xmax paper

Tools related to the Pierre Auger Observatory

auger.charge_fit_from_auger(log10e)[source]

Gives charge fractions of H/He/CNO/Fe for highest energies oriented on Auger data (see GAP2016_047).

Parameters

log10e – Input energies (in log10(E / eV)), array-like. Output charges have same size of first dim.

Returns

charges in size (4, shape(log10e)), last dimension contains fractions in order H/He/CNO/Fe

auger.convert_spectrum(data_17)[source]

Converts units of ICRC17 spectrum to 2015 form: units before are: Flux J*E*10**(-27) [ m^-2 s^-1 sr^-1 ] units after are: Flux J [ eV^-1 km^-2 sr^-1 yr^-1 ]

Parameters

data_17 – numpy array with 2017 spectrum data

Returns

converted array, same form

auger.event_rate(log10e_min, log10e_max=21, zmax=60, area=3000, year=17)[source]

Cosmic ray event rate in specified energy range assuming a detector with area ‘area’ and maximum zenith angle cut ‘zmax’. Uses AUGERs energy spectrum.

Parameters

  • log10e_min – lower energy for energy range, in units log10(energy/eV)

  • log10e_max – upper energy for energy range, in units log10(energy/eV)

  • zmax – maximum zenith angle in degree (default: 60)

  • area – detection area in square kilometer (default: Auger, 3000 km^2)

  • year – take ICRC 15 or 17 data

Returns

event rate in units (1 / year)

auger.geometrical_exposure(zmax=60, area=3000)[source]

Geometrical exposure with simple maximum zenith angle cut and certain area.

Parameters

  • zmax – maximum zenith angle in degree (default: 60)

  • area – detection area in square kilometer (default: Auger, 3000 km^2)

Returns

geometrical exposure, in units sr km^2

auger.gumbel(x, log10e, mass, model='EPOS-LHC', scale=(1, 1, 1))[source]

Gumbel Xmax distribution from [1], equation 2.3.

Parameters

  • x – Xmax in [g/cm^2]

  • log10e – energy log10(E/eV)

  • mass – mass number

  • model – hadronic interaction model

  • scale – scale parameters (mu, sigma, lambda) to evaluate the impact of systematical uncertainties

Returns

G(xmax) : value of the Gumbel distribution at xmax.

auger.gumbel_cdf(x, log10e, mass, model='EPOS-LHC', scale=(1, 1, 1))[source]

Integrated Gumbel Xmax distribution from [2]

Parameters

  • x – upper limit Xmax in [g/cm^2]

  • log10e – energy log10(E/eV)

  • mass – mass number

  • model – hadronic interaction model

  • scale – scale parameters (mu, sigma, lambda) to evaluate the impact of systematical uncertainties

Returns

integral -inf, x of G(xmax) : value of the Gumbel distribution

auger.gumbel_parameters(log10e, mass, model='EPOS-LHC')[source]

Location, scale and shape parameter of the Gumbel Xmax distribution from [1], equations 3.1 - 3.6. parameters from [10]

Parameters

  • log10e (array_like) – energy log10(E/eV)

  • mass (array_like) – mass number

  • model (string) – hadronic interaction model

Returns

mu (array_like, location paramater [g/cm^2]), sigma (array_like, scale parameter [g/cm^2]), lamda (array_like, shape parameter)

Return type

tuple

auger.gumbel_sf(x, log10e, mass, model='EPOS-LHC', scale=(1, 1, 1))[source]

Integrated Gumbel Xmax distribution from [2]

Parameters

  • x – lower limit Xmax in [g/cm^2]

  • log10e – energy log10(E/eV)

  • mass – mass number

  • model – hadronic interaction model

  • scale – scale parameters (mu, sigma, lambda) to evaluate the impact of systematical uncertainties

Returns

integral x, inf of G(xmax) : value of the Gumbel distribution

auger.ln_a_moments(log10e, mass, weights=None, bins=None)[source]

Energy binned <lnA> and sigma^2(lnA) distribution

Parameters

  • log10e (array_like) – energies in [EeV]

  • mass (array_like) – mass numbers

  • weights (array_like) – optional weights

  • bins (array_like) – energies bins in log10(E/eV)

Returns

tuple consisting of:

  • energy bin centers in log10(E/eV), array_like

  • <ln(A)>, mean of ln(A), array_like

  • sigma^2(ln(A)), variance of ln(A) including shower to shower fluctuations, array_like

auger.ln_a_moments2xmax_moments(log10e, m_ln_mass, v_ln_mass, model='EPOS-LHC')[source]

Translate <lnA> & Var(lnA) into <Xmax> & Var(Xmax) according to [3,4].

Parameters

  • log10e – Array of energy bin centers in log10(E/eV)

  • m_ln_mass – <ln(A)> in the corresponding energy bins

  • v_ln_mass – Var(ln(A)) in the corresponding energy bins

  • model – Hadronic interaction model

Returns

Tuple consisting of

  • mXmax : <Xmax> in the corresponding energy bins

  • vXmax : Var(Xmax) in the corresponding energy bins

auger.mean_xmax(log10e, mass, model='EPOS-LHC')[source]

<Xmax> values for given energies log10e(E / eV), mass numbers A and hadronic interaction model, according to [3,4].

Parameters

  • log10e (array_like) – energy log10(E/eV)

  • mass (array_like) – mass number

  • model (string) – hadronic interaction model

Returns

mean Xmax value in [g/cm^2]

auger.plot_spectrum(ax=None, scale=3, with_scale_uncertainty=False, year=17)[source]

Plot the Auger spectrum. 2017 spectrum is in arbitrary units because J0 is not given in publication

auger.plot_var_ln_a(ax=None, model='EPOS-LHC', with_legend=True)[source]

Plot the Auger Var(lnA) distribution.

auger.rand_charge_from_auger(log10e, model='EPOS-LHC', smoothed=None, year=17)[source]

Samples random energy dependent charges from Auger’s Xmax measurements (arXiv: 1409.5083).

Parameters

  • log10e – Input energies (in log10(E / eV)), array-like. Output charges have same size.

  • model – Hadronic interaction model [‘EPOS-LHC’, ‘QGSJetII-04’, ‘Sibyll2.1’]

  • smoothed – if True, smoothes the charge number (instead binned into [1, 2, 7, 26])

  • year – year of publication (15 or 17 available)

Returns

charges in same size as log10e

auger.rand_charge_from_exponential(log10e)[source]

Samples random energy dependent charges from charge_fit_from_auger().

Parameters

log10e – Input energies (in log10(E / eV)), array-like. Output charges have same size.

Returns

charges in same size as log10e

auger.rand_energy_from_auger(n, log10e_min=17.5, log10e_max=None, ebin=0.001, year=19)[source]

Returns energies from the analytic parametrization of the Auger energy spectrum units are 1/(eV km^2 sr yr)

Parameters

  • n – size of the sample

  • log10e_min – minimal log10(energy) of the sample

  • log10e_max – maximal log10(energy) of the sample: e<emax

  • ebin – binning of the sampled energies

  • year – take 15 or 17 ICRC data

Returns

array of energies (in log10(E / eV))

auger.rand_gumbel(log10e, mass, size=None, model='EPOS-LHC')[source]

Deprecated funcion -> See rand_xmax()

auger.rand_xmax(log10e, mass, size=None, model='EPOS-LHC')[source]

Random Xmax values for given energy E [EeV] and mass number A, cf. [1].

Parameters

  • log10e (array_like) – energy log10(E/eV)

  • mass (array_like) – mass number

  • model (string) – hadronic interaction model

  • size (int or None) – number of xmax values to create

Returns

random Xmax values in [g/cm^2]

Return type

array_like

auger.spectrum(log10e, weights=None, bins=None, normalize2bin=None, year=17)[source]

Differential spectrum for given energies [log10(E / eV)] and optional weights. Optionally normalize to Auger spectrum in given bin.

Parameters

  • log10e – Input energies (in log10(E / eV))

  • weights – Weight the individual events for the spectrum

  • normalize2bin – bin number to normalize

  • year – take data from ICRC 15/17

Returns

differential spectrum

auger.spectrum_analytic(log10e, year=17)[source]

Returns a analytic parametrization of the Auger energy spectrum units are 1/(eV km^2 sr yr) for cosmic-ray energy in log10(E / eV)

Parameters

  • log10e – Input energies (in log10(E / eV))

  • year – take ICRC 15, 17 or 19 data

Returns

analytic parametrization of spectrum for given input energies

auger.var_xmax(log10e, mass, model='EPOS-LHC')[source]

Shower to shower fluctuations sigma^2_sh(Xmax) values for given energies log10e(E / eV), mass numbers A and hadronic interaction model, according to [3,4].

Parameters

  • log10e (array_like) – energy log10(E/eV)

  • mass (array_like) – mass number

  • model (string) – hadronic interaction model

Returns

variance of Xmax in [g/cm^2], due to shower to shower fluctuations

auger.xmax_acceptance(x, log10e, zsys=0, fov_cut=True)[source]
Xmax acceptance from [4] parametrized as a constant with exponential tails
exp(+ (Xmax - x1) / lambda1) Xmax < x1
eps(Xmax) = | 1 for x1 < Xmax < x2
exp(- (Xmax - x2) / lambda2) Xmax > x2
Parameters

  • x – Xmax,true in [g/cm^2]

  • log10e – log10(E/eV)

  • zsys – standard score of systematical deviation

  • fov_cut – was the fiducial volume cut applied?

Returns

Relative acceptance between 0 - 1

auger.xmax_energy_bin(log10e)[source]

Get xmax energy bin for given log10(energy) log10

Parameters

log10e – energy log10(E/eV)

Returns

Xmax energy bin number

auger.xmax_moments(log10e, mass, weights=None, model='EPOS-LHC', bins=None)[source]

Energy binned <Xmax>, sigma^2(Xmax), cf. arXiv:1301.6637

Parameters

  • log10e – Array of energies in [eV]

  • mass – Array of mass numbers

  • weights – Array of weights (optional)

  • model – Hadronic interaction model

  • bins – Array of energies in log10(E/eV) defining the bin boundaries

Returns

tuple consisting of

  • lEc : Array of energy bin centers in log10(E/eV)

  • mXmax : Array of <Xmax> in the energy bins of lEc

  • vXmax : Array of sigma^2(Xmax) in the energy bins of lEc

auger.xmax_moments2ln_a_moments(log10e, m_xmax, v_xmax, model='EPOS-LHC')[source]

Translate <Xmax> & Var(Xmax) into <lnA> & Var(lnA) according to [3,4].

Parameters

  • log10e – Array of energy bin centers in log10(E/eV)

  • m_xmax – <Xmax> in the corresponding energy bins

  • v_xmax – Var(Xmax) in the corresponding energy bins

  • model – Hadronic interaction model

Returns

tuple consisting of

  • mlnA : <ln(A)> in the corresponding energy bins

  • vlnA : Var(ln(A)) in the corresponding energy bins

auger.xmax_resolution(x, log10e, zsys=0, fov_cut=True)[source]

Xmax resolution from [4] parametrized as a double Gaussian R(Xmax^rec - Xmax) = f*N(sigma1) + (1-f)*N(sigma2)

Parameters

  • x – Xmax,rec in [g/cm^2]

  • log10e – log10(E/eV)

  • zsys – standard score of systematical deviation

  • fov_cut – was the fiducial volume cut applied?

Returns

Resolution pdf

auger.xmax_scale(log10e, zsys)[source]

Systematic uncertainty dX on the Xmax scale from [4] Xmax,true is estimated to be within [sigma-, sigma+] of the measured value.

Parameters

  • log10e – log10(E/eV)

  • zsys – standard score of systematical deviation

Returns

Systematical deviation dX = Xmax,true - Xmax,measured