Average shape of longitudinal shower profiles measured at the Pierre Auger Observatory

The average profiles of cosmic ray shower development as a function of atmospheric depth are measured for the first time with the Fluorescence Detectors at the Pierre Auger Observatory. The profile shapes are well reproduced by the Gaisser-Hillas parametrization at the 1% level in a 500 g/cm2 interval around the shower maximum, for cosmic rays with log(E/eV)>17.8. The results are quantified with two shape parameters, measured as a function of energy. The average profiles carry information on the primary cosmic ray and its high energy hadronic interactions. The shape parameters predicted by the commonly used models are compatible with the measured ones within experimental uncertainties. Those uncertainties are dominated by systematics which, at present, prevent a detailed composition analysis.


Introduction
The Fluorescence Detector (FD) of the Pierre Auger Observatory [1] has collected unprecedented statistics of high quality data, imaging the longitudinal profile of the electromagnetic component of showers induced by ultra-high energy cosmic rays.The integral of the shower profiles gives direct calorimetric measurements of the primary particle energy, with small corrections due to the energy carried away by muons and neutrinos [2].The depth at which the profile maximum occurs, X max , is the main observable for analysis of the mass composition of primary cosmic rays.
This work follows the hybrid data reconstruction procedure and event selection developed earlier for a composition unbiased X max measurement [3].The determination of each shower geometry uses the timing of the triggered FD pixels and one station in the Surface Detector, which should be at less than 1.5 km from the shower core, and have an expected trigger probability greater than 95%, for proton and iron showers.The energy deposit profile is obtained from the photons detected at the telescope, taking into account the different emission mechanisms and the light absorption and scattering in the atmosphere; only events with viewing angles above 20 • from the shower axis are accepted, to avoid direct Cherenkov light affecting the determination of X max .In addition, a fiducial field of view is defined as a function of energy to ensure a uniform acceptance of most of the X max distribution observed in the data; finally, for each shower track, at least 300 g/cm 2 must be observed.This should include the X max position, for which the expected resolution, given the geometry, must be under 40 g/cm 2 .
The main goal of this contribution is to present a precise measurement of the shape of the longitudinal electromagnetic profile of cosmic ray showers [4], and of the remaining information that it can carry about the interactions of primary cosmic rays.

Measuring the profile shape
Because individual events are subjected to high statistical uncertainties, we construct a high accuracy average profile that can be analysed in detail and provide information about the sample.The profile of each event is normalized to its maximum energy deposit (dE/dX) max , and centred at X max .For each 10 g/cm 2 bin of X = X − X max , the normalized energy deposit is averaged, taking into account the respective uncertainties of each contributing event.
Figure 1(top) shows a reconstructed average profile in data, and the fractions of light components contributing at each X .After the geometry of each event has been obtained, the light arriving at the FD at a given time is converted to its emission point at the shower axis, taking into account the attenuation by Rayleigh scattering in air and Mie scattering on aerosols.Most of the detected photons are from isotropic fluorescence emission, but a large component of the forward Cherenkov beam -integrated along the shower axis -can reach the telescope through scattering.The later part of the profile has a lower contribution of fluorescence light (the one directly proportional to the energy deposition) -and so will be more subject to corrections of atmospheric effects and assumptions used in the reconstruction.Figure 1(bottom) compares the average profiles of simulated energy deposition to the ones obtained after a full detector simulation and applying the same reconstruction procedure as used for data.The full chain is applied separately to proton and iron showers to show that the early part of the profile is the one that keeps information on the interaction of the different primary particles, and that this information is preserved by the reconstruction of the high accuracy average profile.
Figure 2 shows again the reconstructed average profile in data, in the region of [−300, +200] g/cm 2 around the maximum, for which the detailed quantitative analysis is done.Also shown are the bin-by-bin residuals to the fitted Gaisser-Hillas function [5], as an experimental test that it is a good representation of the shower longitudinal energy deposit profiles.The Gaisser-Hillas shape is confirmed at a 1% level in each of six energy bins with log(E/eV) > 17.8.The Gaisser-Hillas function is written as [6]: where R = λ/|X 0 |, L = |X 0 |λ and X 0 ≡ X 0 −X max .L can be seen as a Gaussian width, and R an asymmetry parameter, with smaller correlations than the λ and X 0 parameters that are more commonly used.These parameters keep information even when applied to average profiles [7].

Systematic uncertainties
Table 1 shows the maximum deviation caused by several categories of systematic uncertainties, compared to the much higher statistical accuracy.Most of the effects on the shape parameters are asymmetric and energy dependent.In the table only the highest values are presented.The small differences found between the simulated and reconstructed profiles (fig. 1) do not change significantly the values of R and L obtained in the chosen fit region, except for the first energy bin (up to 1 EeV).A small bias correction is determined by averaging the effect over proton and iron to return to the values as obtained with the generated dE/dX values.Half of this correction is added as a systematic uncertainty, presented together with the effect of the energy scale uncertainty [8].
The FD consists of 24 telescopes viewing elevations from 2 to 32 degrees, overlooking the 3000 km 2 Surface Detector (SD) which samples the particles arriving at ground with an array of stations separated by 1.5 km 1 .One of the telescopes was excluded from the start of the analysis, since it showed large time residuals in the geometry fit.The stability of the average profile across the remaining 23 telescopes was used to check for detector effects, and the differences accounted for as a small detector systematic uncertainty.The effect of the hybrid geometry reconstruction is studied directly by varying within their uncertainties five independent parameters, representing the shower detector plane and the shower axis.Possible variations on the average shower profiles obtained with different selections in shower zenith angle or distance to the telescope were searched for, but found to be contained within the previously derived geometric systematic uncertainties in L and R.
The determination of single event profiles, and their X max and (dE/dX) max , is clearly a very fundamental piece of the analysis.In particular, the Gaisser-Hillas fit of single event profiles uses constraints to guarantee the convergence when short tracks are observed, affecting mainly low energy events.These were derived earlier from a low statistics sample of events with long tracks, and are expressed as λ = L • R = (61 ± 13) g/cm 2 and X 0 = X max + L/R = (121 ± 172) g/cm 2 ; an additional constraint on the normalized integral was derived from simulations, compatible with all available models and compositions.Changing the constraints coherently by one standard deviation in each parameter results in systematic deviations of up to 3 g/cm 2 in L and 0.01 in R.
The fluorescence and Cherenkov yields were also varied within uncertainties and the reconstruction was repeated with and without multiple atmospheric scattering corrections, but a larger effect comes from separating events according to the percentage of light associated with direct fluorescence light: when selecting only events with more than 10% of non-direct fluorescence light (the mean value in the analysed sample), the resulting average profile is wider, with L increased by up to 2 g/cm 2 (while the change in R is negligible).
The atmospheric effects play an important role in the profile reconstruction and the Pierre Auger Observatory also has a large set of instruments dedicated to monitoring the atmospheric conditions [1].The data was separated according to the seasons of the year, and the impact of cloud patches in the sky was also studied.The main effect comes, however, from the overall aerosol content and height profile.The atmospheric aerosol attenuation, τ A (h), is obtained hourly by firing a laser from the centre of the SD array [9,10] and comparing the observed number of photons in the FD with the ones observed in a clear reference night.The correlated uncertainties in the definition of the clear night and the uncorrelated uncertainties from the hourly measurements were propagated by reconstructing again the full sample, and they lead to the largest systematic effect in the values of L, that changes by ± 5 g/cm 2 , and R, that changes by ± 0.02.

Two new observables
Figure 3 shows the measured values of L and R, compared to the predictions of three hadronic interaction models for proton and for iron primaries.We conclude that the measured shower shapes are compatible with model predictions within uncertainties and so that the models provide a fair description of the electromagnetic shower component, measurable by the Fluorescence Detectors.The predictions of the models for the L variable are different, with Sibyll 2.3c [12] predicting higher values than either QGSJetII-04 [13] or EPOS-LHC [14].The data is compatible with both proton or iron values in Sibyll 2.3c, while it is above the prediction for pure iron from the other two models.On the other hand, the three models predict similar values for R, higher for iron than for proton showers, and decreasing with energy for each pure composition.In this case, the data indicates a slower energy evolution, but the present systematic uncertainty prevents an analysis in terms of primary mass composition.
Figure 4 puts together the two observables in two energy bins, highlighting the correlations of statistical and systematic uncertainties.It also shows the prediction for all possible mass compositions (obtained by combining different fractions of H, He, Ni and Fe), according to each of the high energy hadronic interaction models.The combination of the two variables allows a clearer separation for different predictions.For example, in the energy bin of log (E/eV) ∈ [18, 18.2], the profile shape parameters predicted by QGSJet-II04 for any composition are at more than one standard deviation from the values allowed by data, while the shape predicted by Sibyll 2.3c is fully compatible with the measurements for both energies.A measurement of the average shape with lower systematic uncertainties would provide constraints on the hadronic models, specially when put together with other composition sensitive variables.

Summary
We present the first measurement of the shape of the average shower profiles as a function of atmospheric depth, with high statistical precision, and a 1% level check that these follow the Gaisser-Hillas parametrization in the region of [-300, 200] g/cm 2 around the shower maximum.
The analysis provides a systematic check of the shower reconstruction procedures used in the Pierre Auger Observatory, and it is concluded that the present measurement of the profile shape parameters is mostly affected by the (E/eV)  L vs R for the energy bin 10 18 to 10 18.2 eV (left) and from 10 18.8 to 10 19.2 eV (right).The inner dark grey ellipse shows the fitted value for data and its statistical uncertainty, and the outer light grey area the systematic uncertainty.For each hadronic model proton, helium, nitrogen and iron showers were simulated and average profiles were built making all possible combinations.Each of the points represents the value of L and R for a given model and composition combination, so the phase space spanned by each model is contained in its respective coloured area.Pure proton is, for each model, on the upper left side (high L and low R) and the transition to iron goes gradually to the lower right side.
uncertainties in the description of the height profiles of atmospheric aerosols.
The predictions of the hadronic interaction models for the longitudinal profile of the electromagnetic shower component are compatible with the measurements, within present experimental uncertainties.However, it is also shown that the two new observables (L and R) that describe the average profile shape have the potential of further discriminating between models and mass composition scenarios.

Figure 1 .
Figure 1.Average longitudinal shower profile in data (top) and simulation (bottom) for events with energies between 10 18.8 and 10 19.2 eV.Top: The measured normalized energy deposit is shown in black, and the coloured regions (detailed in the legend) represent the average fraction of direct and scattered Cherenkov light in the photons measured at the telescope aperture, computed in the individual shower reconstruction.Bottom: Both generated and reconstructed profiles are shown (and their ratio in the inset plot) in blue for proton and red for iron simulated showers.

Figure 2 .
Figure 2. Measured average longitudinal shower profiles for energies between 10 18.8 and 10 19.2 eV.Data is shown together with the Gaisser-Hillas fit to the profile.The residuals of the fit are shown in the bottom inset.

Figure 4 .
Figure 4. L vs R for the energy bin 10 18 to 10 18.2 eV (left) and from 10 18.8 to 10 19.2 eV (right).The inner dark grey ellipse shows the fitted value for data and its statistical uncertainty, and the outer light grey area the systematic uncertainty.For each hadronic model proton, helium, nitrogen and iron showers were simulated and average profiles were built making all possible combinations.Each of the points represents the value of L and R for a given model and composition combination, so the phase space spanned by each model is contained in its respective coloured area.Pure proton is, for each model, on the upper left side (high L and low R) and the transition to iron goes gradually to the lower right side.

Table 1 .
Breakdown of systematic uncertainties for R and L. Uncertainties are energy dependent and asymmetric, so only the largest value is reported.
[11]re 3. R (left) and L (right) as a function of energy.The data are shown in black, with the vertical line representing the statistical error and the brackets the systematic uncertainty.Air showers were simulated with CORSIKA[11]using different interaction models (see legend).The predictions from simulations are shown with full lines for proton and dashed lines for iron primaries.