Investigating cosmic rays and air shower physics with IceCube / IceTop

IceCube is a cubic-kilometer detector in the deep ice at South Pole. Its square-kilometer surface array, IceTop, is located at 2800 m altitude. IceTop is large and dense enough to cover the cosmic-ray energy spectrum from PeV to EeV energies with a remarkably small systematic uncertainty, thanks to being close to the shower maximum. The experiment offers new insights into hadronic physics of air showers by observing three components: the electromagnetic signal at the surface, GeV muons in the periphery of the showers, and TeV muons in the deep ice. The cosmic-ray flux is measured with the surface signal. The mass composition is extracted from the energy loss of TeV muons observed in the deep ice in coincidence with signals at the surface. The muon lateral distribution is obtained from GeV muons identified in surface signals in the periphery of the shower. The energy spectrum of the most energetic TeV muons is also under study, as well as special events with laterally separated TeV muon tracks which originate from high-pT TeV muons. A combination of all these measurements opens the possibility to perform powerful new tests of hadronic interaction models used to simulate air showers. The latest results will be reviewed from this perspective.


Introduction
Cosmic-ray research is a conundrum.Despite high-quality data and sophisticated theoretical models, many open questions remain.We can identifiy two reasons why this is so, one concerns the propagation of cosmic rays in space, the other the indirect detection of high-energy cosmic rays.The second issue is connected to theoretical uncertainties in the description of soft quantum chromodynamics, and therefore to an active field of fundamental physics research.
The first obstacle is that cosmic rays are electrically charged.They are bent onto complex paths by the various magnetic fields present in outer space [1].The arrival direction of a cosmic ray at Earth can be accurately measured, but it cannot be traced back to its source.Therefore, we cannot select between several plausible theories of origin and propagation of cosmic rays [2].Searches for sources are reduced to searches for anisotropies in the cosmic ray arrival directions.Several experiments have observed such anisotropies (a recent summary can be found in Ref. [3]).
The two remaining identifying properties of a cosmic ray are its energy and mass.Measurements of the flux of cosmic rays combined with the mass composition are able to discriminate between different origin and propagation scenarios [2], but here we face the second obstacle.
Cosmic rays above about 0.1 PeV are too rare and too energetic to be directly measured in detectors mounted on balloons or satellites.One retreats to indirect observation through air showers initiated in Earth's atmosphere.The a e-mail: hdembins@mpi-hd.mpg.deb http://icecube.wisc.edudownside of this approach is that the inference of energy and mass from air shower data is now model-dependent and requires accurate simulation of several cascades of hadronic particle interactions in the regime of soft quantum chromodynamics.
Since these simulations carry uncertainties themselves, they have to be tested in experiments.This requires both dedicated analyses of data from collider experiments which measure particles going into the forward direction [4,5], and complementary analyses of data from air shower experiments themselves.The latter can test models by measuring several observables for the same set of air showers simultaneously.A model passes the test if it is able to describe all observables consistently.To avoid degeneracy, at least three observables are needed; two to determine the primary energy and mass, and a third to test the hadronic model.
Muons in air showers play an important role here.They carry information about the last hadronic interaction that created their parent pion.Due to the near lightlike propagation of muons, this information is reflected in the lateral density profile and in the arrival times of muons when an air shower reaches the ground.Several experiments have measured the lateral density profile and found a consistent deficit in the simulated numbers of muons in the lateral profile compared to experimental data [6][7][8].This discrepancy is called the muon puzzle.
The IceCube Neutrino Observatory is an excellent detector for measuring cosmic rays and muons in air showers.It has unique capability to shine new light on the muon puzzle: the ability to separately measure the GeV and TeV muon component in air showers (see Fig. 1).Measurements further benefit from the high altitude of IceTop, its surface detector array.The proximity of an array to the shower maximum significantly reduces the systematic uncertainties in the inferred cosmic-ray energy.
Predictions for muon densities vary in the two energy regimes for different hadronic interaction models, as shown in Fig. 2. For determining the cosmic-ray mass, measurements of TeV muons are better suited, because the model-uncertainty is smaller compared to the effect from the primary mass.A simultaneous measurement of the GeV and TeV muon abundance allows one to test hadronic interaction models.
Although not discussed in this proceeding, we also emphasize that the absence of TeV-scale muons in an otherwise normal-looking air shower is strong evidence for a shower produced by a PeV-scale gamma ray.Limits on the flux of PeV gamma rays have been published [14], and studies with more recent data are ongoing.
In the following we will summarize results obtained with IceCube on the flux of cosmic rays, their mass composition, and on air shower physics.

IceCube Neutrino Observatory
The IceCube Neutrino Observatory, shown in Fig. 3, is a cubic-kilometer Cherenkov detector located at the South Shown are 2 PeV air showers with a zenith angle of 60 • simulated using CORSIKA [9] with several different hadronic interaction models [10][11][12][13] to a slant depth of 1000 g.cm −2 .The plot was kindly contributed by F. Riehn.Pole [15].The main detector, IceCube, is embedded in the Antarctic ice shelf.Its surface component, IceTop [16], is the primary detector for cosmic-ray induced air showers.IceCube detects charged particles via Cherenkov light generated in the ice, which is collected and digitized by Digital Optical Modules (DOMs).The bulk of the DOMs are attached to 86 strings in the deep ice, with a horizontal spacing of 125 m and a vertical spacing of 17 m.The DeepCore sub-detector has a denser instrumentation [17].IceTop uses the same DOMs embedded in light-tight tanks filled with ice.
IceCube's buried array is protected from low-energy cosmic radiation by a 1.5 km-thick ice shield.The protective layer imposes an energy threshold of about 0.3 TeV for particles with vertical incidence.At the surface, 81 strings are equipped with an IceTop station.Each station consists of two IceTop tanks separated by 11 m.Each tank contains two DOMs with different gain settings to cover a large dynamic signal range.

Cosmic-ray energy spectrum and mass composition
Air shower properties are primarily reconstructed from the size and time of signals collected by IceTop.The air shower direction is reconstructed from the arrival times of signals in at least three stations.The achieved angular resolution is better than one degree above a PeV [16].
The energy of an air shower is inferred from the shower size parameter S 125 , the equivalent signal at a lateral distance of 125 m from the shower axis.The equivalent signal is obtained by fitting a model of the lateral density to the measured signals.The energy resolution obtained is better than 25% above 2 PeV, improving with energy.
Thanks to the high altitude of the IceTop array, the conversion of S 125 to energy has a very small systematic uncertainty of only 10%, which is illustrated by the overlap of the proton and iron showers in Fig. 4. Another important systematic effect are the varying snow layers on top of the IceTop tanks.While there is very limited snow fall at the South Pole, there is considerable snow drift, which leads to snow buildup of up to 1 m per year.The snow buildup is precisely monitored by shifters at the South Pole and is corrected for using an exponential attenuation model [16].The snow-corrected signals show a negligible remaining dependency on the snow height, except for sub-PeV showers near the trigger threshold of the array, which are not used in standard analyses.
The small systematic uncertainty of the energy scale leads to a corresponding small systematic uncertainty in the flux of cosmic rays, shown in Fig. 5.The latest spectrum was derived from three years of data and shows two peaks, the so-called knee around 5 PeV and a second "knee" around 100 PeV.The discovery of the second knee structure is relatively recent, but it is now established by at least three different experiments [19,20].Another estimate of the total flux was derived from coincident events seen in IceTop and IceCube, which does not depend on the H4a model and is in perfect agreement with the main result [18].
IceCube was not built with intrinsic particle identification capabilities, but it turns out that muons can nevertheless be separated from other charged particles.GeV muons are detected in IceTop on a statistical basis at large lateral distance from the shower axis, where they leave a characteristic signal (see next section).TeV muons are the only charged particles able to penetrate the 1.5 kmthick ice shield to reach IceCube.
TeV muons are primarily used to estimate the cosmicray mass.The number of TeV muons can be roughly deduced from the reconstructed energy loss per path length dE/dX in the ice.For the mass composition analysis, IceCube uses the median energy loss for a reference path length of 1500 m, shown in Fig. 4. The median energy loss is the most important mass-sensitive observable.Since TeV muons are not minimum-ionizing, they also undergo stochastic losses.These show up as spikes in the energy loss profile.The number of such stochastic losses above two different thresholds is used to further refine the discrimination of different masses.
The shower size S 125 , the zenith angle of the shower, the median energy loss, and the two counts of stochastic energy losses are fed into a neutral network trained on simulations to produce estimates of the energy and mass of the cosmic rays.
The output of the network is binned in the energy and mass estimator for proton, helium, nitrogen, and iron showers simulated with CORSIKA and SIBYLL2.1 [9,10], and likewise for real events.Binning the output for simulated known masses produces templates, which are fitted to the binned histograms of the mass estimator obtained for real events.This allows us to determine which fraction of the showers in each energy bin are similar in mass to proton, helium, and so on.
With this method, the component spectra shown in Fig. 6 are produced.We observe an initially light composition that becomes increasingly heavy as the energy increases, then stabilizes around 100 PeV [18].Our measurement indicates a heavier composition around 1 EeV than measurements from the Pierre Auger Observatory based on the depth of shower maximum [21].This apparent inconsistency could originate from a deteriorating description of TeV muons in SIBYLL2.1 at such high energies (hints for such a deterioration are found in GeV muons and discussed in the next section).We emphasize, however, that the results were crosschecked using QGSJet-II.03[22], which gave a consistent interpretation in terms of cosmic-ray mass.
A sudden drop in the helium and iron spectra around 6 PeV is observed, with corresponding elevated levels of proton and oxygen.This feature is under intense study.So far the most likely explanation is a statistical fluctuation.We note that the plotted standard error bars cannot be interpreted in the usual way, because they do not show the correlations.The four component spectra are not statistically independent, since they sum up to the smooth total flux by construction.

GeV muons in IceTop
IceTop tanks are ice-Cherenkov detectors.A minimumionizing particle which passes vertically through the detector produces a characteristic signal, which produces a stable peak in a histogram of photo-electron counts.IceTop uses this very feature for its continuous signal calibration.Signals expressed in units of Vertical Equivalent Muons (VEMs) are comparable from tank to tank and between real and simulated events, while the number of photoelectrons collected per vertical muon hit may vary from tank to tank and over time.
It was determined by IceTop that this effect can be used to statistically separate muon hits from hits of other charged particles in tanks far from the shower axis.If the tank is sufficiently far from the shower axis, the probability of an individual tank to be hit by a particle is much smaller than one.Then, the typical charge deposited by a single muon (for vertical showers it is around 1 VEM by definition of the unit) stands out statistically over the continuum of charges from electrons and photons.
By counting how many tanks in a narrow ring around the shower axis have a muon hit, the average muon density for air showers within a narrow interval in cosmic-ray energy can be derived using Poisson statistics [23,24] and the known cross-sectional area of an IceTop tank.For near vertical showers, the average lateral muon density profiles derived in this way are shown in Fig. 7.These profiles are raw results, not corrected for biases in the method.
A detailed study of simulated proton and iron showers showed that the muon density obtained in this way is slightly biased.The bias is corrected using these simulations, by taking the average effect for proton and iron showers into account.This leads to a small systematic uncertainty due to the mass composition.Another systematic uncertainty arises from a variation of the empirical statistical model of the non-muon hits.
For two lateral distances at 600 m and 800 m the muon density is shown as a function of the air shower energy in Fig. 8, and compared to simulations for proton and iron ISVHECRI 2016  showers.The measurement is dominated by the systematic uncertainty, which however remains below 15%.
For the hadronic interaction model SIBYLL2.1, we find a similar behavior in GeV muons and in TeV muons.The results are compatible with a light mass composition around 1 PeV and then indicate an increase in the average mass.However, around 100 PeV the muon density rises even beyond the value for iron showers, which is astrophysically not plausible.The more likely explanation is a deterioration of the description of muons in SIBYLL2.1 as the air shower energy increases.This interpretation is consistent with another measurement by the HiRes-Mia group [6], which measured the muon density in air showers from 100 PeV to 1 EeV and found a similar deficit of muons in simulations compared to real showers.

High-energy muon flux
As previously mentioned, TeV muons in the deep ice occasionally have stochastic energy losses which generate large localized amounts of light.This experimental signature can be used to identify air showers events with high-energy muons (HE-muons).HE-muons are typically produced by prompt decays of charmed hadrons, which also give rise to prompt neutrinos.
The IceCube collaboration is very interested in reliable predictions of the flux of PeV neutrinos that originate from air showers, because they form the main background to astrophysical neutrinos.Measuring the flux of HEmuons allows one to constrain the flux of atmospheric PeV neutrinos very well.
In a recent study [25], machine learning tools were applied to discriminate air showers with a HE-muon from other events, which was defined as a single muon that carries more than 50% of the energy of the whole muon bundle that reaches the deep detector.The classifier was trained on simulated air showers, using 30 variables which quantify various aspects of the air shower event in IceCube.
After training, a signal-to-background ratio of better than 10 is achieved, with about 40% of the signal surviving.
Estimating the energy of the HE-muon based on its energy loss in the deep ice is still difficult, so the energy resolution is not very high.Unfolding methods are used to compensate distortions in the estimated flux due to low energy resolution, and also correct for the detector acceptance and selection losses.The final result is the unfolded flux of HE-muons at the South Pole surface, shown in Fig. 9.The measurement is very promising, but so far still compatible with both a flux predicted by SIBYLL2.1, in which charmed hadrons are not simulated, and an alternative flux where decays of such charmed hadrons have been added.
The situation will improve in the near future.The current measurement is limited by the statistical accuracy of the simulated response matrices which are used in the flux unfolding.More simulations are in preparation to significantly reduce the flux uncertainty between 0.1 and 1 PeV.

Laterally-separated TeV muons
The last muon analysis presented here targets a special class of air shower events discovered by IceCube.These  events feature a coincidence of two parallel tracks in the deep ice, which are interpreted as a muon bundle and a single laterally separated muon (LS-muon), see Fig. 10.An early measurement of the distribution of the lateral separation d T , shown in Fig. 11, has been already published already in 2013 [26].
Because of the light-like propagation of muons, the lateral separation d T of the LS-muon is geometrically connected to the transverse momentum p T of the meson which decayed into the LS-muon, roughly described as where H int is the height of the primary interaction of the cosmic ray, E µ is the energy of the LS-muon and θ is the zenith angle of the air shower.Using simulations, it was shown that hadrons which produce LS-muons are themselves produced predominantly in the primary interaction of the cosmic  + x/400) n , which describes a transition from an exponential, expected for a soft process, and a power law, expected for a hard process.The P-value for this fit is 4%.For a purely exponential function it is 0.001% [26].ray.LS-muons are therefore direct probes of the primary interaction, which is very sensitive to the mass of the cosmic ray.Simulated distributions of transverse momentum of hadrons produced in the first interaction vary greatly from model to model and with the mass of the cosmic ray, as shown in Fig. 12.
A new study is underway, which benefits from a larger data set and a new sophisticated simulation approach for laterally separated muons.Muons with large transverse momentum are separately generated and then embedded in regular CORSIKA air showers [27,28].In the new approach, all details in the production of LS-muons can be controlled.Also, the Monte-Carlo production of such special events is greatly accelerated, which rarely occur in standard CORSIKA simulations.
The analysis of LS-muons is particularly interesting.If the high-p T tail of the LS-muon distribution is dominated by hard processes, the tail could be predicted to necessary accuracy based on pertubative quantum chromodynamics.This opens a possibility to avoid the theoretical uncertainties which are always present when conventional hadronic interaction models are used to infer the cosmic-ray mass [29].The change of slope in Fig. 11 is a hint that the tail is indeed dominated by hard processes.If this approach turns out to be feasible, it would allow us to infer the mass composition of PeV cosmic rays to unprecedented theoretical accuracy.

Conclusion and outlook
The IceCube Neutrino Observatory is a powerful and versatile cosmic-ray detector.It has the unique ability to detect TeV muons, which along with neutrinos are the only air shower particles able to reach the deep detector under the ice shield.For the first time, predictions from hadronic interaction models for TeV muons can be experimentally tested.New measurements, which have not been possible before, include the flux of muons with hundreds of TeV and the lateral distribution of high-p T muons.
Thanks to the high altitude of the IceTop array, the energy of air showers can be measured with systematic uncertainties smaller than 10%, which results in a very precise measurement of the cosmic-ray flux.IceTop is among the first experiments which observed the so-called second knee structure around 100 PeV.
Finally, the ice-Cherenkov tanks of the IceTop array allow us to statistically identify GeV muon hits based on the deposited signal alone, and correspondingly measure the lateral density of GeV muons far from the shower axis without additional detector equipment.The preliminary results indicate a deficit of GeV muons in air showers simulated with SIBYLL2.1 around 100 PeV, with is consistent with similar measurements done by other experiments at higher cosmic-ray energies.
Studies are currently ongoing for a proposed extension of the IceCube Neutrino Observatory, shown schematically in Fig. 13.The goal is to increase IceCube to 10 km 3   and likewise IceTop to 10 km 2 .This would extend the accessible range of cosmic-ray energies to 3 EeV and the rate of coincident events observed in both IceTop and IceCube by a factor of 50, due to the increased acceptance in the zenith angle.Also planned is an even larger veto array of about 75 km 2 made of simpler detectors to veto air showers for neutrino studies [30].
Both the near and far future of cosmic-ray research with IceCube/IceTop are very promising, with many interesting results in preparation.

Figure 1 .
Figure 1.Sketch of an air shower which develops over the IceCube Neutrino Observatory.The IceTop array is close to the shower maximum.Signals near the shower axis are dominated by contributions from electrons and photons, which allows us to infer the shower energy with comparably small systematic uncertainty.Signals far from the shower axis are dominated by hits of GeV muons.TeV muons form collimated particle bundles, and are the only air shower particles apart from atmospheric neutrinos which penetrate the ice shield and reach the deep detector.

Figure 2 .
Figure 2. Average ratio of muon spectra in air showers initiated by proton and iron showers as a function of the muon energy.Shown are 2 PeV air showers with a zenith angle of 60 • simulated using CORSIKA[9] with several different hadronic interaction models[10][11][12][13] to a slant depth of 1000 g.cm −2 .The plot was kindly contributed by F. Riehn.

Figure 3 .
Figure 3.The IceCube Neutrino Observatory.Charged particles produce Cherenkov light in the clear deep ice at the South Pole.Cherenkov light is detected by Digital Optical Modules (DOMs) placed on 86 strings between roughly 1450 m to 2450 m below the surface.The strings in the center are more densely instrumented (DeepCore).The IceTop array is formed by 81 stations with two ice-filled tanks each.

Figure 4 .
Figure 4. Top: S 125 as a function of primary energy in simulated air showers, which has minimal sensitivity to composition.Bottom: Muon energy loss at a path length of 1500 m in the ice for the same air showers, which is highly composition-sensitive[18].

Figure 5 .
Figure 5. Unfolded all-particle flux of cosmic rays obtained from 3 years of IceTop data, derived using the H4a mass composition assumption.The gray band indicates the systematic uncertainty of the flux [18].

Figure 6 .
Figure 6.Flux of cosmic rays split into components (colored dots), based on similarity to simulated proton, helium, oxygen, and iron showers.The sum of the components is shown in each plot for reference (black dots) [18].

Figure 7 .
Figure 7. Raw lateral profiles of muon density observed in IceTop for near vertical showers (average zenith angle is 13 • ) obtained for air showers in several energy bins (average energy is shown in rotated text).Dotted and dashed lines indicate estimates based on alternative models used to fit the background of non-muon signals.The impact of non-muon background on the analysis is larger for bins closer to the shower axis up to a point where the systematic uncertainty becomes too large and results are discarded (open dots).The profiles are raw results obtained from fits of signal distributions, before applying a composition-dependent correction derived from air shower simulations.

Figure 8 .
Figure 8. Top: muon density at 600 m and 800 m as a function of primary energy for near vertical showers.Lines show expected values for proton and iron showers simulated with CORSIKA and SIBYLL2.1.Bottom: same as above, but with all densities normalized to the expected density for proton showers.Brackets indicate systematic uncertainties.The densities shown here were corrected for a small composition-dependent bias found when applying the method to air shower simulations.

ISVHECRI 2016 Figure 9 .
Figure 9. Unfolded spectrum of TeV muons measured with IceCube (black dots).Shown for comparison are a previous result (gray line and error band), a prediction computed with SIBYLL2.1 without prompt hadrons (dashed line), the same with charmed hadrons (ERS) added (solid line) [25].

Figure 10 .
Figure10.Event in IceCube with a muon bundle and a laterallyseparated muon track that arrives in coincidence and with nearparallel direction.The experimental observable for these type of events is lateral distance d T between bundle and track.

Figure 11 .
Figure 11.Distribution of lateral separation muon bundle and isolated muon, after subtracting a background contributed by random coincidences.Also shown are the best-fit parameters for the empirical model N = exp( A + Bx) + 10 C (1+x/400) n , which describes a transition from an exponential, expected for a soft process, and a power law, expected for a hard process.The P-value for this fit is 4%.For a purely exponential function it is 0.001%[26].

Figure 12 .
Figure 12.Distributions of transverse momentum for hadrons produced in the first interaction of a 0.1 PeV proton or iron cosmic ray with nitrogen, as predicted by several hadronic interaction models [27].

Figure 13 .
Figure 13.Proposed for the Generation-Two IceCube Neutrino Observatory.IceCube is extended by the High-Energy Array to 10 km 3 , similarly IceTop is extended to 10 km 2 .A larger surface of about 75 km 2 is covered by the Surface Veto Array, which consists of simpler detectors and used to veto air shower events for neutrino searches.