Comparison of the Measurement and Simulation with KM2A Prototype Array

. The Large High Altitude Air Shower Observatory (LHAASO) is a new hybrid array for very high energy gamma ray astronomy and cosmic ray physics. The KM2A array, one of the main parts of LHAASO, covers an area of 1.3 km 2 to observe gamma rays above 10 TeV up to 1 PeV for many sources. A prototype with 1% the size of the whole KM2A has been in stable operation for more than two years. A Monte Carlo simulation program named G4KM2A was developed; based on this work, the trigger rate, hit multiplicity, angular and core reconstruction are compared with KM2A prototype data. Finally, the moon shadow with -6.5 signiﬁcance was obtained. The G4KM2A simulation results are consistent with KM2A prototype data and can be used for the whole KM2A array in future.


Introduction
Very high energy (VHE>30 GeV) gamma-rays, a neutral radiation carrying the information of primary sources without any deflection, play an important role to solve the origin and acceleration mechanism of cosmic rays. The flux of VHE gamma-rays is so low that detectors with large effective area are very essential for high energy gamma-ray observations. LHAASO is a large scale hybrid extensive air shower (EAS) array for VHE gamma-ray astronomy and cosmic ray physics, and consists of one square kilometer extensive air shower array (KM2A) covering an area about 1.3 km 2 , 78000 m 2 water Cherenkov detector array (WCDA) and 12 wide-filed air Cherenkov/fluorescence telescopes (WFCTA) [1].
The KM2A array consists of 5195 electromagnetic detectors (EDs, 1m 2 each) and 1171 muon detectors (MDs, 36m 2 each) to detect secondary particles (γ,e ± ,µ ± ) in the EAS with a full duty cycle, and is dedicated to the observation of gamma-ray sources from 10 TeV to 1 PeV [2,3]. The whole construction of KM2A array will be finished by 2021, and a prototype with 1% the size of KM2A was constructed at the YangBaJing cosmic ray observatory, and operated from 2014 September to 2016 October. The time calibration for all data has been completed, and some data has been collected to analyze muonic and electromagnetic components in cosmic ray air showers [4].
In order to achieve the original energy, component, flux, direction and other properties of primary events, a Monte Carlo simulation program G4KM2A based on GEANT4 was developed [5]. G4KM2A has some advantages such as the array distribution can be changed conveniently for different requirements. In this work, the KM2A prototype data were used to check G4KM2A to ensure its validity for the whole array. The moon shadow with two years observation of prototype was also analyzed.

The KM2A prototype
The KM2A prototype consists of 39 EDs and 2 MDs, as shown in Figure 1. The EDs are designed mainly for the electromagnetic particles in EAS, with a time resolution better than 2 ns. The MD is a kind of water Cherenkov detector to detect µ ± with a threshold of about 1 GeV; each MD is located under a shielding of 2.5 m thick overburden soil to prevent low energy electromagnetic particles reaching the water. A detailed description of ED and MD is given in [1,6], and the performance of the KM2A prototype array met the original design [7].
A real shower event is only recorded when more than 5 detectors are fired (one fire is denoted as one "hit") within a time window of 200 ns, and thus the array is recognized to be triggered. All hits in a time window of 5000 ns before and after this trigger are recorded as a complete event time, so the noise distribution is mainly before 4500 ns. The noise rate can be obtained by the prototype data with about 2000Hz for all 39 EDs, and was used to set the noise rate value in G4KM2A.
Some information about the shower and the triggered detector are recorded to reconstruct the origin energy and direction of the primary event, such as the arrival time of electromagnetic particles in the air shower. The hit detector's ID (with fixed coordinate) and integral electronic count are recorded in a data file. The long time performance of the KM2A prototype detectors was reported in [4], and the mean value of the trigger rate is about 50.4 Hz.

KM2A simulation and reconstruction
The CORSIKA package is used to simulate the propagation of EAS through atmosphere. The composition and energy spectrum of P, He, CNO, MgAlSi and Fe was simulated following the H3a model according to Gaisser [8] from 100 GeV to 10 PeV. The interaction models in COR-SIKA are QGSJETII to describe high-energy hadronic interactions and GHEISHA to describe hadronic collisions up to some 100 GeV. The zenith angles are from 0 • to 60 • , and azimuth angles from 0 • to 360 • . The G4KM2A based on the GEANT4 tool package is used to simulate the detector response to secondary particles in EAS [5].
The effective area of the detector is different according to the shower size for lower and higher energy primary events. In the simulation, the shower core position is randomly sampled within a radius of 300m for energies lower than 100TeV, 600m for 100TeV-1000TeV and 1200m for 1PeV-10PeV around the detector center. The areas are large enough to cover all locations of primary events.
Secondary electromagnetic particles of one cosmic ray event form a conical plane in front of the shower caused by different interaction lengths with the atmosphere. So the direction of primary cosmic rays can be reconstructed using this time information. A conical fitting procedure by minimizing following χ 2 function is commonly used to reconstruct the direction, where l=sinθcosφ, m=sinθsinφ, θ is the zenith angle, φ is the azimuth angle, α is cone factor, 0.03 is used in this work, c=29.98cm/ns, (x i ,y i ) is the coordinate of the i th detector, t i is the recorded time of one hit, D i is the distance between detector and shower axis calculated using the reconstructed core location. The shower core location can be reconstructed by fitting the NKG function.

Comparison with experiment
The trigger mode in the simulation program is the same as in the experiment, i.e. one shower should have fired at least five detectors simultaneously (N hit ≥5) within 200ns. The noise rate input in G4KM2A is equal to that in the KM2A prototype, about 2000Hz. Trigger rate, hit multiplicity, shower direction reconstruction and core position reconstruction will be compared between experiment and simulation in this section.

Trigger rate
The simulated energy distributions for all particles are shown in Figure 2 and normalized to the spectrum given by Gaisser [8]. From this simulation result, the mode-energy of the KM2A prototype for all particles is about 19-31TeV. The differences of mode-energy, are mainly caused by the different cross section of each component interacting with nuclei in the atmosphere. Another aspect, the trigger result can also be influenced by simulation statistical numbers of origin shower in a light extent. Trigger rates for each component are listed in Table 1. The main nuclei are protons and helium with 29.1Hz and 16.0Hz respectively. Another three components are much smaller with 6.5Hz in total. So the main components in cosmic rays detected by the KM2A prototype array are H and He nuclei which occupied about 87.47% in these five components. The total trigger rate in the simulation is 44.03 Hz which is a little lower than the mean trigger rate of two years experiment data as shown before.

Hit multiplicity
Higher energy primary cosmic rays can produce more secondary particles, and more detectors will be hit, so the distribution of hit number can reflect the primary energy of cosmic ray in a preliminary way. Because the thickness of the shower front can make the event detection continue for some time, and one ED can be fired many times when the arrival time period between two hits is larger than 8ns, therefore the hit above 39 are produced, as shown in Figure 3. The simulated result matches experiment much better when the hit number is less than 50, the disagreement in higher hit numbers may be caused by the array's detection character which is under consideration to improve these differences. The relationship between hit numbers and primary energy is obtained as shown in Figure 4.

Reconstruction result
The times recorded in the KM2A prototype data have been calibrated for reconstruction. The reconstruction result of  the shower core and direction are finally obtained by fitting the NKG function and equation 1. The distance (R) between the reconstructed core location and the array center are displayed in Figure 5 (left), and a comparison of the χ 2 value is shown in Figure 5 (right). These two comparison results imply that simulation matched experiment data well, and give very close reconstruction results.
Since the atmosphere thickness increases with increasing zenith angle, the number of EAS events is strongly related to θ, as shown in Figure 6. The reconstructed direction of data and MC sample coincide with the primary direction. The azimuth angle distribution shows a modulation effect in different directions, which may be caused by the array structure, and this modulation is greater at larger zenith angles as shown in Figure 7.   The direction differences between primary and reconstructed data can be used to illustrate the angular resolution of the detector array. For this KM2A prototype array, the space angle ( θ) has been fitted by the following equation where θ r is the reconstructed zenith angle, φ r is the reconstructed azimuth angle, σ illustrates the angular resolution and concerned to point spread level. The fitting result is shown in Figure 8, and the angular resolution with all cases is about 1.32 • for this prototype array, which meets the design requirement for future KM2A physical studies.  Figure 9. The moon shadow obtained by KM2A prototype data.

Moon shadow
The Moon blocks galactic cosmic rays casting a shadow that appears as a deficit region in the relative region, which is called moon shadow [9]. The data of the KM2A prototype array with two years operation period was taken into account to analysis the moon shadow, as shown in Figure 9. This map is based on a 2D histogram on a grid of 0.1 • ×0.1 • with a smoothing angular radius 1.32 • . The shadow position, deficit event, angular spread of measurement are shown in Table 2 by fitting the shadow area with a 2-D Gaussian function. The most significant point with -6.5σ of the moon shadow was obtained using the Li-Ma formula [10]. The deficit extension of the moon shadow is determined by fitting the distribution of Ω 2 as shown in Figure  10, where Ω is the angular distance to the position of the Moon shadow. The measured deficit events for each nest ring is the sum of (N on − N b ) on this scale, and the expected deficit event is the whole moon shadow that should be distributed in relatively ring scale caused by the angular spread effect. N on is the on-source observation data, N b is off-source observation data. The measured deficit event and angular spread is consistant with expectation value. This result confirmed that the KM2A prototype operation and its time calibration are successful. The G4KM2A program performance and reconstruction method are generally correct for the KM2A prototype array.

Conclusion
A G4KM2A program based on the GEANT4 package has been developed to simulate the KM2A array. The G4KM2A can get a relatively much close detector response to particles in EAS as compared with experiment data, such as trigger rate, hit multiplicity, reconstruction direction and core location. The angular resolution for the KM2A prototype array is about 1.32 • . Finally the moon shadow with two years observation of KM2A prototype array was obtained with significance about -6.5σ. The shadow position, deficit event number and elongation angle distribution show a consistancy between experiment and expectation. This result confirmed the time calibration of the KM2A prototype array and reconstruction method are valid. In future, this G4KM2A program will be further modified and used for the whole KM2A array.
This work is supported in China by NSFC (No. 41504080, No. 11575203 and No. 11505028).