Air shower simulation for background estimation in muon tomography of volcanoes

One of the main sources of background for the radiography of volcanoes using atmospheric muons comes from the accidental coincidences produced in the muon telescopes by charged particles belonging to the air shower generated by the primary cosmic ray. In order to quantify this background effect, Monte Carlo simulations of the showers and of the detector are developed by the T OMUVOL collaboration. As a first step, the atmospheric showers were simulated and investigated using two Monte Carlo packages, CORSIKA and GEANT4. We compared the results provided by the two programs for the muonic component of vertical proton-induced showers at three energies: 1, 10 and 100 TeV. We found that the spatial distribution and energy spectrum of the muons were in good agreement for the two codes.


Introduction
Two aspects of muon tomography of volcanoes make relevant, if not necessary, the use of Monte Carlo simulations to accurately evaluate the attenuation of the flux of atmospheric muons crossing the volcano.The first aspect is related to the fact that the incident muon flux has to be described with sufficient precision and analytical approximations are not accurate enough, in particular at low muon energies or for nearly horizontal muons.The atmospheric muon flux depends not only on the altitude and on the zenith angle of the muons, but also on their energy, which is generally not measurable with the experimental setups of the present tomography experiments.Also, for low energy muons (below a few GeV), the Earth's magnetic field induces a dependence of the flux with the azimuth.
The second argument in favour of developing Monte Carlo simulations comes from the intense atmospheric background the detector is exposed to while performing the tomography.In particular, various particles produced in air showers can hit the detector coherently, creating accidental coincidences and faking the passage of a muon through it.Due to the randomness of those events, the frequency at which they occur and the characteristics of the induced tracks are very hard to determine analytically.
In this context, a detailed Monte Carlo study of the atmospheric muon flux is being conducted in TOMUVOL (Cârloganu et al., 2012).As a first step, a simulation of high energy protons interacting with the atmosphere and creating showers was performed using GEANT4 (Agostinelli et al., 2003).After briefly reviewing here how the atmospheric muons are created during the development of air showers, we present the simulation code in development within the TO-MUVOL collaboration, in the framework of GEANT4.Then, the muonic component of the vertical showers induced by protons of 1, 10 and 100 TeV, simulated with GEANT4, is compared with the muonic component of similar showers simulated with CORSIKA (Heck et al., 1998).CORSIKA is a dedicated software for the simulation of atmospheric showers over a very large energy spectrum, widely used in cosmic ray experiments.

The physics of air showers
Cosmic rays are particles, mostly protons (∼ 88 %) and helium nuclei (∼ 10 %), originating from various (extra)galactic sources that interact in our atmosphere.High energy cosmic rays initiate events called air showers, i.e. extensive cascades of particles collimated in the direction of the incident cosmic particle.The interaction of the incoming cosmic particle with an atmospheric nucleus is characterised by the fragmentation of the nucleus and the creation of a number of secondary hadrons.These can, if energetic enough, experience in turn inelastic interactions with others atoms, or decay into hadrons, photons and leptons.In this scheme, the atmospheric muons (Cecchini and Spurio, 2012) used in tomography are mostly coming from the decays of two kind of hadrons: charged pions and kaons.As pions constitute the most important part of the hadronic component of the showers, and since the charged ones decay almost exclusively into muons, their contribution is dominant.It is worth mention-ing that in addition to the hadronic decays described above, a tiny fraction of atmospheric muons is created by the conversion of photons into a pair of muons.The charm contribution through the decay of η, D's and c particles (the so-called prompt contribution) does not become important until very high energies and will be neglected here (Costa , 2001).
The processes contributing to the muon flux are given below with the corresponding life-times and branching ratios: The development of the hadronic showers depends not only on the energy and on the type of the incoming cosmic ray, but also on the amount of matter the cosmic particle will cross when propagating in the atmosphere.The secondary particles created in the shower lose part of their energy when propagating.As an example, high energy pions tend to interact before decaying, leading to less higher energy muons being created and therefore to a steeper energy spectrum.
Although the muons will interact only relatively weakly with the atmosphere, they also lose a small fraction of their initial energy until they reach the ground, where they are detected.In the case of horizontal muons, the very long path they have to travel in the atmosphere until ground leads to significant energy losses, and an important fraction of the low energy horizontal muons decay before reaching the ground.This is the reason why the energy spectrum of the horizontal muons is much flatter than in the case of the vertical muons.
We see that a multitude of physical processes are involved between the arrival of a cosmic particle at the top of the atmosphere and the detection of a muon at ground level.These processes are of stochastic nature, and will occur many times on average for each event.Thus, to extrapolate the expected muon flux at detection level from the (relatively well-known) primary cosmic ray flux, the use of Monte Carlo simulation tools seems to be a convenient way to go.In the next sections, the efforts made in TOMUVOL to develop such a Monte Carlo simulation are presented.

Simulation of air showers using GEANT4
GEANT4 is a widely used simulation software, designed to simulate the interaction of particles with matter in a fully customizable environment.The advantage of this program in the context of muon tomography is the possibility to describe in the same simulation framework the development of air showers, the interactions of the muons in the volcano, and finally  their detection.The full program chain is currently in development within the TOMUVOL collaboration, and we focus here on the part dealing with the development of air showers.The present simulation setup consists of an atmosphere modelised with 1960 spherical layers of air, with temperature and density set according to the mid-latitude winter model (Kneizys et al., 1996).This data model is also currently used in CORSIKA.At the top of the atmosphere, vertical protons of predefined energies are shot and let to interact with the air molecules.The used interaction model is QGSP BERT (GEANT4 Collaboration, 2012;Apostolakis et al., 2010;Dotti, 2011), which is standard in several high energy physics experiments, completed with the γ → µ + µ − process.The output of the code is a file containing the information about the particle flux at an altitude of 870 m, which is the current TOMUVOL detector operation altitude in front of the Puy-de-Dôme.
One of the features of GEANT4 is that each particle carries information about its creation.In particular, it is possible to know the nature of the particles which created the muons crossing the plane z = 870 m.As an example, the energy spectrum of the muons created in showers initiated by 10 TeV incident protons is given in Fig. 1 for all the different production channels active in our physics list.As a very general observation, one can see that the relative contributions from the different mother particles are compatible with the description of the formation of the atmospheric muons flux made in the previous section.Before extending the simulation to the generation of a cosmic ray flux which includes all the nuclear species and covering a wide energy range and the full sky aperture, a more accurate check of the generated muon fluxes is needed.In order to do so, the code was tested by comparing its results with those from another simulation software, CORSIKA.

Comparison with CORSIKA results
CORSIKA is a program specifically designed for the simulation of air showers, commonly used in cosmic ray physics, maintained and improved for decades as new experimental data on air showers was obtained.Therefore, the results provided by this program will be considered as a reference in this study.In order to validate our GEANT4 code with respect to the simulations of CORSIKA, the atmospheric muons fluxes recorded at an altitude of 870 m (a.s.l.) in both programs were compared.The GEANT4 data which have been used corresponds to vertical proton primary particles of energies 1, 10 and 100 TeV.In CORSIKA we simulated a realistic cosmic ray spectrum, over the full energy and zenith angle ranges and including also Helium nuclei.To have a meaningful comparison with the GEANT4 simulation, we restricted the CORSIKA sample to proton-initiated showers with zenith angles less than 10 • and energies within 10 % of the corresponding energy of the GEANT4 showers.Note that a cut on the energy of the muons, E µ > 10 GeV, was applied, so the effects of the Earth's magnetic field, not yet included in our GEANT4 code but present in CORSIKA, could be neglected.The muons at z = 870 m were characterised by their energy and their spatial and arrival time distributions within the shower.The energy spectra of the muons (the number of entries in the histogram is normalised to the total number of muons present in the showers) are shown in Fig. 2 for both GEANT4 simulation and CORSIKA, and the agreement is good.
For the spatial distributions, the mean centre of gravity of all the muons in the shower at the z = 870 m plane was subtracted from the position of each individual muon, such that the normalised spatial distributions shown in Fig. 3 are centred on zero.Again, there is a reasonable agreement between the GEANT4 and CORSIKA simulations, though CORSIKA displays a slight tendency towards wider high energy showers.
The analysis concerning the arrival time profiles of the muons obtained from the two softwares is still ongoing.Indeed, the distributions obtained in CORSIKA and GEANT4 after performing our analysis are not in agreement with each other and are both showing unexpected features.This effect needs to be investigated because the number of accidental coincidences induced by the showers in the muon telescope depends on the time distribution of the particles reaching the detection level.

Conclusions
The first step toward the Monte Carlo study of the atmospheric background in muon tomography is the simulation of the air showers.A choice has been made in TOMUVOL to develop such a program in the framework of the GEANT4 simulation code, since the next stage of this work should involve the simulation of the detection of the particles in the telescope, which GEANT4 is designed for.In order to validate the part of the code related to the shower simulation, the muon fluxes at an altitude of 870 m produced by vertical 1, 10 and 100 TeV primary protons were compared to those from the well-tested CORSIKA air shower simulation code.We found that the lateral and energy profiles of the muons were in relatively good agreement in both codes.The arrival time distributions of the muons that we obtained from GEANT4 and CORSIKA were not in agreement, and the analysis to determine the origin of this effect is still ongoing.As soon as the GEANT4 results for vertical showers can be validated against CORSIKA, an inclusive cosmic ray flux will be simulated.Together with the detector simulation, this should allow us to estimate the contamination from the accidental coincidences induced by the air showers that could be expected in our detector configuration.

Fig. 1 .
Fig. 1.Energy distributions, at an altitude of 870 m, of the muons (E µ >1 GeV) from showers initiated by 10 TeV vertical protons, as simulated in GEANT4.The contributions from the different mother particles are dissociated: pions (solid line), charged kaons (dashes), neutral kaons (dot dashes), gammas (dots).The vertical axis indicates the mean number of particles per shower.

Fig. 2 .
Fig. 2. Energy distributions of the muons at z = 870 m as obtained with CORSIKA (dashes) and GEANT4 (solid line), normalised to the total number of muons produced in each simulation code.

Fig. 3 .
Fig.3.Lateral distributions of the muons at z = 870 m as obtained with CORSIKA (dashes) and GEANT4 (solid line), normalised to the total number of muons produced in each simulation code.r is the distance to the centre of gravity of the muons in each shower, computed in the detection plane.
et al.: Air shower simulation for background estimation in muon tomography of volcanoes