Studies of Cosmic Rays at the Highest Energies with the Pierre Auger Observatory

This paper summarizes the status and the recent measurements of the Pierre Auger Observatory. The energy spectrum is described and its features discussed. Searches for anisotropy of cosmic ray arrival directions on large scales and through cor- relation with catalogues of celestial objects are reported. The first measurement of the proton-air cross section around 10 18 eV is discussed. The mass composition is addressed with measurements of the variation of the depth of shower maximum with energy and with muon density at the ground. An update on the searches for neutrinos and photons is also presented.


Introduction
Ultra high energy cosmic ray (UHECR) research probes physical laws at the highest energies, connecting the astrophysics of extreme physical systems with particle physics. In this way it explores energies and kinematic regions hitherto unexplored and continues to motivate a rich and diverse experimental program. For energies up to 10 15 eV, cosmic rays are believed to have a galactic origin and shock acceleration in supernova remnants could be the most likely source. At the highest energies, the most probable sources of UHECRs are extragalactic: jets of active galactic nuclei (AGN), radio lobes, gamma-ray bursts and colliding galaxies, among others [1]. In this paper we summarize the main experimental results from the Pierre Auger Observatory on measurements of UHECRs, the highest energy particles measured on Earth, with energy E > 0.1 EeV (1 EeV = 10 18 eV).

The Pierre Auger Observatory
The Pierre Auger Observatory is located on a plain near the town of Malargüe, in Mendoza Province, Argentina. The surface detector array (SD) [2] is composed of 1660 water-Cherenkov stations placed on a triangular grid of 1500 m that covers 3000 km 2 . Each station is filled with 12 tons of water and instrumented with three photomultipliers which detect the Cherenkov light produced by the charged particles that cross the water. The signals are read out at 40 MHz and time stamped by a GPS unit, making it possible to measure the arrival time profile of the incoming shower particles. The fluorescence detector (FD) comprises 24 optical telescopes grouped in four sites, observing the atmosphere above the array. Each telescope consists of an 11 m 2 segmented spherical mirror that focuses the fluorescence light into a camera composed of 440 photomultipliers [3]. The signal detected by them is digitized at 10 MHz providing a time profile of the shower as it deposits energy in the atmosphere.
The observatory was designed as a hybrid detector, in the sense that it combines the two different and complementary techniques of UHECR detection. This enables the energy of each event to be derived in a calorimetric way through the calibration of the shower size measured with the SD array by the energy measured with the FD telescopes on a subset of high quality hybrid events [4].
To determine the energy measured by the FD in a reliable way a few parameters have to be under control: the fluorescence yield (the proportionality factor between the number of photons emitted and the energy deposited in the atmosphere), the vertical aerosol optical depth (VAOD) that is important for estimating the attenuation of of light in the atmosphere, and finally the absolute calibration of the telescopes.
The overall up-time and efficiency of the SD is around 98% while that of FD is 13% which operates on clear moonless nights. The energy resolution of the SD alone is 12% (statistical) above 10 EeV [5] while the angular resolution is less than 1 • for this energy range [6]. The detector components are monitored second by second by hardware and software tools. As an example, the activity of each SD station is reported each second and is averaged in Fig. 1 [7].
To measure lower energy cosmic rays, a subarray of SD stations, known as "the infill", is composed of 71 stations on a denser grid of 750 m with an area of 30 km 2 [8]. This subarray is part of the AMIGA extension that will have buried muon counters at each station location. Already seven are in place [9]. Also, three High Elevation Auger Telescopes (HEAT) located at one of the fluorescence sites are dedicated to the detection of lower energy showers [10]. Other instruments for radio and microwave detection are operating onsite [11].

Energy spectrum
The energy spectrum reflects the various aspects of the propagation, production and source distribution of the cosmic rays. The measurement of the flux as a function of the energy of the incident cosmic ray exhibits two features at the highest energies: a flattening of the flux around 5 × 10 18 eV, known as "the ankle", and a strong suppression above 5 × 10 19 eV. These features can be explained by the interaction of the particles with the cosmic microwave background, but they also may be due to changes in the acceleration origin or power [12][13][14][15]. A necessary ingredient for the discrimination between theoretical models is the precise measurement of the flux of cosmic rays as a function of energy.
The Auger energy spectrum is measured with four datasets: hybrid, infill, SD vertical (zenith angle < 60 • ) and SD horizontal (62 • < zenith angle < 80 • ). For SD vertical events the signal from the particles reaching the ground is sampled with the SD water-Cherenkov detectors. The content of the

E [eV]
Auger 2013 preliminary Figure 2. Left: The energy spectra of the four Auger datasets [17]. Right: The combined energy spectrum compared to energy spectra from different astrophysical scenarios [11].
air-showers at this level depends on the amount of matter traversed by the cascade: for vertical events both the electronic and muonic component are present at the ground. S(1000), the signal at 1000 m from the axis of the air-shower, is interpolated from a lateral distribution function. A correction for the attenuation in the atmosphere is applied to S(1000) using a constant intensity method: showers with a higher zenith angle produces a signal smaller than showers with a lower zenith angle. This correction produces the final estimator for the vertical events: S38, the equivalent signal for a shower that arrives with a zenith angle of 38 degrees. The energy estimator for the other datasets is discussed in [16]. Golden hybrid events (that have independently triggered the SD array and FD telescopes) are used for the energy calibration of SD data. Fig. 2, left, shows the energy spectra obtained from the three SD datasets and the hybrid data [17]. Due to the calibration with events observed by the FD, the SD energies share the uncertainty of the FD energy scale of 14%. To produce the combined spectrum shown in Fig. 2, right, the uncertainties are taken into account in a maximum likelihood fit together with the smearing corrections due to binto-bin migrations. The unprecedented statistical accuracy allows one to identify two clear features in the spectrum: the ankle at log (E a /eV) = 18.72 ± 0.01 ± 0.02, where the spectral index changes from −3.23 ± 0.07 to −2.63 ± 0.04, and a flux suppression at log (E 1/2 /eV) = 19.63 ± 0.01 ± 0.01, where E 1/2 is the energy at which the flux has fallen to one half of the value of the power law extrapolation. When compared to a simple continuation of a power-law, the significance of the suppression is more than 20 sigma.
The combined energy spectrum is compared to fluxes from four astrophysical scenarios in Fig. 2, right [11]. Shown are spectra from models assuming pure proton or iron composition. The fluxes result from different assumptions of the spectral index β of the source injection spectrum and the source evolution parameter m. In those simulations the maximum energy of the source was set to 100 EeV (that best describes the data in the suppression region) and 300 EeV. For proton fluxes the parameter β = 2.35 and m = 5 and for iron fluxes β = 2.3 and m = 0. The model lines have been calculated using CRPropa [18] and validated with SimProp [19]. In spite of the high statistical accuracy of the measurements, it is not possible to distinguish among the various scenarios based only on the the energy spectrum. It is necessary to combine of the spectrum with other observables, such as anisotropies and mass composition. These will be discussed in the following sections.

Search for anisotropies
As mentioned above, the physical origin of the structure in the cosmic ray spectrum known as the "ankle" (at about 5 EeV) is not settled yet. There are two main candidate scenarios to explain it (see Ref. [1] for a review). In the first one, the ankle is related to the transition from a galactic component to a harder extragalactic one [13]. In the second one (so called dip scenario), cosmic rays are supposed to be extragalactic protons down to energies below 1 EeV and the concave shape is due to the effect of energy losses of protons by pair creation with cosmic microwave background photons (CMB) [12]. In order to elucidate this issue, anisotropy measurements are very important. If the ankle is due to the change from a galactic to an extragalactic component, a dipolar modulation in the cosmic ray arrival directions is expected. The amplitude of the dipole should increase with energy up to the ankle, reaching a level of a few percent, even if primaries are heavy nuclei. Predictions of the amplitude and orientation of the dipole depend on the magnitude of the intervening magnetic field and also on the source distribution. If the dip scenario is the correct explanation for the ankle, it means that UHECRs above 10 18 eV are already dominated by the extragalactic component, and their flux is expected to be highly isotropic. A dipole of less than one percent amplitude is expected due to the Compton-Getting effect. It would produce a dipolar distribution constant with energy and with an orientation related to the relative movement of the Earth with respect to the "rest frame" of the extragalactic cosmic rays [20]. Thanks to the high statistics of the SD data, a first harmonic analysis was performed in different energy ranges starting from 2.5 × 10 17 eV in a search for dipolar modulations in right ascension [20,21]. This analysis used the fact that the exposure is almost uniform in right ascension. No significant amplitude was detected. The 99% CL upper limits as a function of energy are shown in Fig. 3, right, together with predictions from two different galactic magnetic field models (A and S), for a purely galactic origin of UHECRs (Gal), and the expectations from the Compton-Getting effect  [20,21]. Right: Upper limit at 99% for the dipole amplitude as a function of energy. Some generic anisotropy expectations from stationary galactic sources distributed in the disk are also shown, for various assumptions on the cosmic ray composition [24,25].
(C-G Xgal) [20]. The upper limits shown here already exclude the model with an antisymmetric halo magnetic field (A) above energies of 0.25 EeV and the Gal model at energies of a few EeV, and are starting to become sensitive to the predictions of the model with a symmetric field [21]. In Fig. 4, left, the phase measurement as a function of the energy shows an interesting pattern: it suggests a smooth transition between a phase of ∼ 270 • (consistent with the right ascension of the galactic center) below 1 × 10 18 eV and another phase of ∼ 90 • (consistent with the right ascension of the galactic anti-center) above 5 × 10 18 eV. This is interesting since a real anisotropy would need less events to be established with high statistical confidence from the phase consistency in ordered energy intervals than by amplitude measurements [22,23]. To test the hypothesis that the phase is undergoing a smooth transition, we began to analyse new data obtained after April 2011. After 18 months, the new and independent data set is showing a similar trend [21]. When the exposure of 21,000 km 2 sr yr is attained with the independent data set, the trend will or will not be confirmed. Searches for dipole and quadruple modulations were conducted simultaneously in declination and right ascension. The upper limits presented in [24,25] are shown in Fig. 4, right. They are presented along with generic estimates of the dipole amplitudes expected from stationary galactic sources distributed in the disk considering two extreme cases of single primaries: protons and iron nuclei. While other magnetic field models, source distributions and emission assumptions must be considered, in these particular examples, the hypothesis that the light component of cosmic rays comes from stationary sources densely distributed in the Galactic disk and emitting in all directions can be excluded. This illustrates the power of the observational limits.
For the highest energy cosmic rays, possible candidates able to accelerate them up to 10 20 eV are the jets of AGN, as mentioned before. The Pierre Auger Collaboration reported [26,27] a correlation of events with the AGN directions in the VCV catalogue. The first 14 events were used for an exploratory scan that yielded the following search parameters: energy threshold (E th = 55 EeV), maximum angular separation (Ψ ≤ 3. The exposure-weighted fraction of the sky covered by the blue circles is 21%. [28]. Right: An updated data set showing the VCV correlation as an energy ordered plot [11,29]. flux were isotropic. The subsequent 13 events established a 99% confidence level for rejecting the hypothesis of isotropic cosmic ray flux. The reported fraction of correlation events was (69 +11 −13 )%. By adding data with E th = 55 EeV up to the end of 2009 (69 events in total), the amount of correlation decreased to (38 +7 −6 )% [28]. A sky map for this dataset is shown in Fig. 5, left. The VCV correlation is also shown as an energy ordered plot in Fig. 5, right. The onset of the correlation signal is visible at about 55 EeV. The current estimate of the fraction of correlating cosmic rays is (33 ± 5)%, with 21% expected under the isotropic hypothesis [29].

Mass composition
The measurement of the mass composition of UHECRs is essential to the solution of the problem of their origin, since the mass and charge (Z) distribution can give powerful constraints on their acceleration mechanisms and propagation. For UHECRs the main observable sensitive to composition is the X max , the average value of the atmospheric depth (measured in g/cm 2 ) where the shower development reaches the maximum. Proton showers reach maximum about 100 g/cm 2 deeper in the atmosphere than iron showers. In a similar way, the fluctuation of the values of X max around the mean, σ(X max ), provides another sensitive observable: iron showers fluctuate about 40 g/cm 2 less than proton showers. Those estimates are obtained from transport codes that simulate the shower development given a model for hadronic interactions. The energy evolution of the two mentioned observables for proton and iron and several hadronic models are shown in Fig. 6 together with the hybrid data [11,30]. Both plots could be interpreted as an evolution from light to heavier composition if current hadronic interaction models describe well the air shower physics.
From the measurements of X max and σ max one can obtain the first moments of the log (mass) distribution [32]. ln A and σ QGSJetII-04 interaction models in Fig. 7 [11,31]. The statistical errors are represented by error bars and systematic errors are given by the shaded bands. The overall features are similar in the plots irrespective of the interaction model assumed. The data indicate an increasing ln A as a function of energy above 2 EeV from light to intermediate masses, and a decreasing σ 2 ln A over the whole energy range. The log (mass) scale changes with the hadronic model studied but this does not prevent important properties of the mass distributions of UHECR being inferred. For some models the central predicted variance of ln A is negative but this is not the case within our systematic uncertainties. It is expected that for future analysis the systematic uncertainties will be reduced and this method can help to assess the validity of the hadronic interaction models.

Hadronic physics
The study of extensive air-showers while they develop in the atmosphere is of central importance to understanding their energy distribution and mass composition, and to probe particle interactions up to the highest energies. As shown in section 5 the uncertainties in the hadronic interaction models and the differences between them make the interpretation of the mass measurements more complicated. To improve the current models, measurements of the particle production process are very important. The Pierre Auger Observatory has two different observables associated with the particle production process: the proton-air cross-section and the number of muons at ground level.
The measurement of the proton air cross-section relies on comparisons of an appropriate air shower observable with Monte Carlo predictions. A disagreement between the measurements and the Monte Carlo predictions is attributed to an inadequate value of the proton-air cross-section. The analysis starts with a measurement of an air shower observable with high sensitivity to the crosssection that is later converted into an estimate of the proton-air cross section for particle production, σ p-air in the energy interval 10 18 to 10 18.5 eV [33]. As the primary observable, we define Λ f by means of the exponential shape dN/dX max ∝ exp (−X max /Λ f ) of the X max -distribution of the fraction f of the most deeply penetrating air showers. In this way the contribution of protons in the sample is enhanced, because the average depth at which showers reach maximum is higher in the atmo-  Figure 7. ln A and σ 2 ln A as a function of energy for different interaction models. Systematics uncertainties are represented as shaded bands [11,31].
The conversion from Λ f to the σ p-air is done with air shower simulations, which introduces systematic uncertainties due to model assumptions. The main sources of systematic uncertainty are the energy scale, the possible presence of photons as primaries, and the fraction of helium in the primary composition. The final result, taking into account the uncertainties above, is ] mb This result can be converted into the more fundamental proton-proton cross-section using the Glauber formalism. Additional systematic uncertainty is added in this model-dependent calculation and the result for the inelastic proton-proton cross-section is [33] σ inel p-p = [92 ± 7(stat) +9 −11 (syst) ± 7(Glauber)] mb The result favors a moderately slow rise of the cross-section towards higher energies, similar to results from LHC [34].
The muon content in the ground signal of individual hybrid events can be compared to the ground signal of simulated showers with matching longitudinal profiles [35]. This constitutes a test of the hadronic models, since the ground-level muonic component of ultra-high energy air showers is sensitive to hadronic particle interactions at all stages in the air shower cascade. The number of muons in the ground signal of the simulated events shows a systematically smaller number of muons than the ground signal in the data events, regardless of the hadronic model being used. To explore the Nam et al. 1975 [30] Siohan et al. 1978 [31] Baltrusaitis et al. 1984 [2] Mielke et al. 1994 [32] Knurenko et al. 1999 [19] Honda et al. 1999 [20] Belov et al. 2007 [18] Aglietta et al. 2009 [33] Aielli et al. 2009 [34] This work 0.9TeV 2.36TeV 7TeV 14TeV Figure 8. Left: Unbinned likelihood fit of Λ f to the tail of the X max distribution. Right: σ p-air compared to other measurements and model predictions. The inner error bars are statistical only, while the outer include all systematic uncertainties for a helium fraction of 25% and 10 mb photon systematics. Both plots taken from [33].  Figure 9. The best-fit values of R E and R μ for QGSJET-II-04 and EPOS-LHC, for mixed and pure proton compositions. The ellipses represent the one-sigma statistical uncertainties. The grey boxes show the estimated systematic uncertainties [35]. possible sources of the discrepancy, the ground signal has been modified in the simulated events to reproduce the ground signal in the data. Two rescaling parameters are introduced: R E , that acts as a rescaling of the shower energy, and R μ which acts as "muonic" rescaling factor because it rescales only the contribution to the ground signal of intrinsic hadronic origin. The results of fits for both rescaling parameters are shown in Fig. 9 for each model. The ellipses represent the one-sigma statistical uncertainty region in the R E -R μ plane. The grey rectangles represent the systematic uncertainties in the event reconstruction of X max , E FD and S(1000) that were propagated through the analysis by shifting the reconstructed central values by their one-sigma systematic uncertainties. A purely protonic composition is given as a benchmark. In all cases the parameter R μ is larger than one, indicating a deficit in the predictions, while R E is compatible with one for the mixed sample and for the pure proton sample only within the systematics. Independent analyses using inclined showers or using the characteristic signal shape left by muons on the SD stations also point to a deficit of muons in the simulations [36,37].

Search for photons and neutrinos
Two broad mechanisms of UHECR production models have been suggested, the "top-down" and the "bottom-up" scenarios. In the "top-down" scenario, UHECR are created in the decay of fossil Grand   References to experiments and models in [38]. Right: Differential and integrated upper limits (at 90% C.L.) from the Pierre Auger Observatory for a diffuse flux of ultra high energy neutrinos. The search period corresponds to 6 yr of a complete SD. Also shown are the integrated limits from ANITAII and RICE experiments, along with expected fluxes for several cosmogenic neutrino models as well as astrophysical sources [39].
Unification defects (topological defects, super-heavy dark matter, etc) , and no GZK cutoff is expected. In the "bottom-up" scenarios, UHECRs are supposed to be accelerated in astrophysical sources, and these should exhibit a GZK cutoff. The search for high energy photons and neutrinos is motivated by the following reasons, among others: (i) "top-down" models predict a significant fraction of photons and neutrinos at the highest energies, (ii) so called "cosmogenic" neutrinos and photons would be produced when cosmic-ray protons interact with cosmic microwave background photons and that could possibly prove the existence of the GZK effect, and (iii) they would open a new window to the most extreme Universe by possibly pointing to sources in the sky.
Air showers induced by photons are expected to develop deeper in the atmosphere (with an average value of X max larger) because of the smaller average multiplicity in electromagnetic interactions as compared to hadronic ones, even though the photon interacts first in the atmosphere. The X max for photon showers is typically a few hundred g/cm 2 larger than nuclei with the same energy. As a consequence of the deeper development and the reduced number of muons, photon showers are expected to have a smaller size at the ground, a steeper lateral distribution of secondary particles, a sparser distribution of arrival times of particles in the shower front and a larger delay with respect to a planar shower front approximation. These characteristic properties of the air showers induced by photons were explored in two independent analyses of Auger data, one using hybrid events and another using SD only events. No photon events have been identified so far with either analysis, but upper limits have been obtained at 95% confidence level on the photon flux integrated above an energy threshold, as shown in Fig. 10, left [38]. From a hybrid analysis, the upper limits are illustrated with red arrows. At higher energies, the upper limits obtained with the SD analysis are illustrated with black arrows.
Neutrinos of all flavours can interact in the atmosphere producing high energy neutrinos (downgoing) in the sub-EeV energy range and above. Tau neutrinos can also be produced by the "Earthskimming" mechanism (up-going). The SD is capable of detecting both types of neutrinos through the broad time structure of the signals induced in the SD stations. A search was performed in three angular ranges, with three different sets of identification criteria established to maximise the discrimination power [39]. The analysis starts with a "trace cleaning' procedure that removes most of the accidental signals (mainly due to atmospheric muons). PMTs not passing quality cuts are removed. After that, inclined showers are identified. The next step is to identify in the data sample containing inclined showers, those with a broad time structure (young showers). The three selection criteria were applied to data between 1 January 2004 up to 31 December 2012 in a search for neutrino candidates. This period includes a new search sample corresponding to data from 1 June 2010 until 31 December 2012 not previously unblinded under any of the three selections. After the unblinding, the compatibility of the distributions of the young shower discriminating observables in the search and training samples was tested. No candidate was found. To put the most stringent limits to the diffuse flux of UHE neutrinos, the three individual exposures were combined into a single total exposure by means of a simple procedure. The updated single-flavour 90% C.L. limit is k 90 < 1.3 × 10 −8 GeV cm −2 s −1 sr −1 and applies in the energy interval ∼ 1.0 × 10 17 − 1.0 × 10 20 eV where 90% of the event rate is expected. The result is shown in Fig. 10, right, along with the differential limit to show the energy at which the sensitivity of the search peaks. Also shown are limits of dedicated neutrino experiments and expected fluxes for several cosmogenic neutrino models as well as astrophysical sources.

Conclusion and outlook
The Pierre Auger Observatory has detected high quality events and has made key measurements of the highest-energy cosmic rays. In spite of that, many issues remain open. The spectrum suppression is well established. The measurements of the large scale anisotropy provides upper limits in the dipole amplitude that constrain theoretical models in the ankle region. The Auger data provide evidence for a weak correlation between arrival directions of cosmic rays above 55 EeV and the positions of AGNs with z < 0.018. The fraction of correlating cosmic rays is about 33% compared with 21% expected for isotropy. The astrophysical interpretation of the data is subject to many unknown parameters in the source properties and distribution, extra-galactic magnetic fields, and other uncertainties. However, the data seems to indicate that the cut-off region represents more a consequence of the source maximal acceleration energy than a propagation effect as expected from the GZK scenario, together with a proton sub-dominant fraction. Still in the cut-off region, another possible interpretation is a change in the hadronic interactions of protons at highest energies [11]. Detailed studies of hadronic interactions in the atmosphere together with a much larger data sample may provide new information that will help to answer the remaining questions about the UHECRs. The Auger observatory will continue to take data and a detector upgrade is under way to improve the muon detection capability.

Acknowledgments
I would like to thank the organizers for the invitation for this very stimulating conference and to my colleagues of the Pierre Auger Collaboration for all the work described in this contribution. I also would like to thank the Universidade Federal do Rio de Janeiro, CAPES and FAPERJ for the partial funding of my scientific activities and CNPq for a research grant.