Determining the muon charge ratio using an experimental measurements and the CORSIKA simulation code

Nuclear Inst. and Methods in Physics Research, A 945 (2019) 162635

Contents lists available at ScienceDirect

Nuclear Inst. and Methods in Physics Research, A journal homepage: www.elsevier.com/locate/nima

Determining the muon charge ratio using an experimental measurements and the CORSIKA simulation code M. Bahmanabadi ∗ Department of Physics, Sharif University of Technology, P.O. Box 11155-9161, Tehran, Iran ALBORZ Observatory, Sharif University of Technology, Tehran, Iran

ARTICLE

INFO

Keywords: Cosmic rays telescope Atmospheric muon charge ratio

ABSTRACT The muon charge ratio contains important information about the flux of atmospheric neutrinos and the hadronic interactions. Using a cosmic ray telescope, the atmospheric muon charge ratio has been studied. The result of this experiment is compared with simulation results using the CORSIKA code.

1. Introduction Muon is a fundamental particle with a rest energy of ∼106 MeV, which have a mean lifetime of 2.2 μs in its rest reference frame. The ratio of the number of positive to negative charge atmospheric muons reaching the Earth’s surface is known as the muon charge ratio, 𝑅𝜇 . The charge ratio of the positive and negative muons, 𝑅𝜇 = 𝑛+ ∕𝑛− , yields important information on the flux of atmospheric neutrinos and the hadron interactions. In the range of low energy muons (∼1 GeV), the muon flux is an important source for the production of electron and muon neutrinos and anti-neutrinos, and the ratio of the number of positive to negative muons directly relates to the number of electron neutrinos to anti-neutrinos (𝑛+ ∕𝑛− = 𝑛𝜈𝑒 ∕𝑛𝜈̄𝑒 ). Since the interaction of positive and negative muons with the matter is not the same, finding the ratio of the positive to negative muons can be important for the response of neutrino detectors. Typically, positive and negative muons and relevant neutrinos and anti-neutrinos are produced from the following decays in the atmosphere as, ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩

𝑒± + 𝜈𝑒 (𝜈̄𝑒 ) + 𝜈̄𝜇 (𝜈𝜇 ) ↗ 𝜇 ± + 𝜈𝜇 (𝜈̄𝜇 ) ; ∼ 100%, 𝜏0 = 26 ns ↖ 𝐴𝐶𝑅 + 𝐴𝐴𝑖𝑟 → 𝜋 ± + 𝜋 0 + K ± + other hadrons. ↙ 𝜇 ± + 𝜈𝜇 (𝜈̄𝜇 ) ; ∼ 63.5%, 𝜏0 = 12 ns ↘ 𝑒± + 𝜈𝑒 (𝜈̄𝑒 ) + 𝜈̄𝜇 (𝜈𝜇 )

where 𝐴𝐶𝑅 represents a cosmic ray (such as proton, alpha, carbon, etc.), 𝐴𝐴𝑖𝑟 an air nucleus (such as nitrogen, oxygen), and 𝜋 ± , 𝜋 0 , and K± pions and kaons, respectively. It is expected from these decays that the ratio of muon neutrinos to electron neutrinos, (𝑛𝜈𝜇 + 𝑛𝜈̄𝜇 )∕(𝑛𝜈𝑒 + 𝑛𝜈̄𝑒 ), is equal

2, but Super-Kamiokande experiments showed a significantly smaller value for this ratio. This anomalous behavior was described based on the neutrino oscillations [1]. Different methods are used to measure the muon charge ratio. In some experiments, spectrometers with a magnetic field are used to separate positive and negative muons from each other. But in low energies, these spectrometers are not suitable due to the presence of many electrons, and the different acceptance for positive and negative particles [2]. Here, we obtain the muon charge ratio of atmospheric muons with a delayed coincidence method, based on the interaction between the positive and negative muons with the material. In fact, positive muons decay with an average lifetime of 𝜏0 ≃ 2.2 μs by interacting with the matter after stopping, while when negative muons stop in matter, they can be bound to the nucleus of the atoms of the matter just the same way as electrons do. The mass of the muon is about 200 times heavier than the mass of the electron, so its orbit is closer to the nucleus and can interact with the proton and produces neutron and neutrino (𝜇 − + 𝑝 → 𝑛 + 𝜈𝜇 ). Since there are two ways to disappear a negative muon, the lifetime of negative muons in matter is somewhat less than the lifetime of the positive muons, which do not interact with the second mechanism. In this way, the negative muon lifetime is as 𝜏1 = 𝛬𝑐 + 𝜏1 , where 𝛬𝑐 is the probability of capturing − ◦ negative muons by the nucleus of the atom and 𝜏0 is its lifetime in the free space. The muon capture rate is proportional to 𝑍 4 , where 𝑍 is the atomic number of the nuclei, so with increasing 𝑧 the muon negative lifetime is decreased. In this paper, using a cosmic ray telescope, the muon charge ratio in the energy range between 0.8 GeV and 1.6 GeV is obtained. Meanwhile, using the CORSIKA simulation code and various models for Hadron interactions, the charge ratio is obtained.

∗ Corresponding author at: Department of Physics, Sharif University of Technology, P.O. Box 11155-9161, Tehran, Iran E-mail address: [email protected].

https://doi.org/10.1016/j.nima.2019.162635 Received 31 July 2019; Received in revised form 19 August 2019; Accepted 25 August 2019 Available online 28 August 2019 0168-9002/© 2019 Elsevier B.V. All rights reserved.

M. Bahmanabadi

Nuclear Inst. and Methods in Physics Research, A 945 (2019) 162635

Fig. 1. A schematic of experimental layout and electronic system used to record data.

2. Description of the experimental setup

that these muons are formed at the highest possible height from the Tehran’s level (897 g/cm2 ), according to the ionization loss relation, the energy of the muons recorded in the detector (0.8 ≤ 𝐸 (GeV) ≤ 1.6) should be between 2.6 GeV and 3.4 GeV at the top of atmosphere. Regarding the relativistic effects, a muon with energy between 2.6 ≤ 𝐸 (GeV) ≤ 3.4 (i.e., the Lorentz factor, 𝛾 = 𝐸𝜇 /(𝑚𝜇 𝑐 2 ), between 24.5 and 32.1) has a lifetime, 𝛾𝜏0 , between 53.9 and 70.62 μs in the Earth’s frame. Hence, they can reach Earth’s surface from a height of more than ℎ = (𝛾 2 − 1)0.5 𝑐𝜏0 ≃ 16 km, but due to the ionization loss, this height is reduced to about 10 km. Thus, the particles recorded in the telescope are muons that are produced at altitudes more than 10 km above Tehran’s surface and with energies between 2.6 and 3.4 GeV, and then decay at the Tehran’s level with energy between 0.8 GeV and 1.6 GeV in the telescope.

A cosmic ray telescope, containing 2 plastic scintillator blocks (30 × 10×1 cm3 ) and 100 cm apart from each other, has been used to find the muon charge ratio. Fig. 1 shows a view of the telescope and the electronic device used. With this arrangement, we record the muons decay. The opening angle of the telescope is about 17◦ relative to the vertical axis. When a muon passes through the top scintillation detector and, of course, passes through the opening angle of the telescope and then decays into the bottom scintillation detector, an event is recorded. The signals produced in photomultiplier tubes (PMT, 9813B with diameters of 5 cm), after passing through a fast amplifier(×10, CAEN N412), are connected to channels 1 and 2 of a 8-channel fast discriminator (CAEN N413A). The threshold level of the top scintillation detector was set at −50 mV and of the bottom one at −23 mV. Each channel of the discriminator has two outputs. One of them connects to a logic unit (CAEN N405) with a gate width of 150 ns. The output of the logical unit was connected via a 10 m coaxial and LEMO cable (equivalent to a time delay of 50 ns) to the START input of a time to amplitude converter (TAC, ORTEC 566), which was set to a full scale of 10 μs. If the muon decays in the bottom scintillator, the electron or positron signal was also sent to the STOP input of the TAC after passing through the amplifier (×10) and the discriminator. Eventually, the TAC output is fed into a multichannel analyzer (MCA) via an analog to digital converter (ADC) unit. In this way, the interval time between the signals produced by the stopping of the muon in the bottom scintillator and the electron generated from decay is recorded by a computer. These experiments were carried out at the Tehran’s level (35◦ 43′ N, 51◦ 20′ E, 1200 m a.s.l= 897 g/cm2 ).

4. Analysis of experimental data With the arrangement shown in Fig. 1, the muon decay experiment was carried out for a one-week period and repeated six times. In all experiments, the telescope was in the right position (𝜃 = 0). A sample of distribution of the Muon decay time is shown in Fig. 2. The time distribution function used for the arrangement is, in fact, a superposition of the decay of two types of muons (𝜇+ + 𝜇 − ) and two types of noises, one with a fast decay and with a time constant 𝜏𝑏𝑔 , and another with a slow decay, which is shown with a constant value 𝑆𝑏𝑔 as follows [2]: 𝑡+𝑡 𝑡+𝑡 𝑑) −( 𝑑 ) 𝑛+ −( 𝑡+𝑡 𝑛 −( 𝑑 ) 𝑑𝑁 = 𝑒 𝜏+ + − 𝑝𝑑𝑒𝑐𝑎𝑦 𝑒 𝜏− + 𝐹𝑏𝑔 𝑒 𝜏𝑏𝑔 + 𝑆𝑏𝑔 . 𝑑𝑡 𝜏+ 𝜏−

3. Investigation of the energy spectrum

(1)

where 𝜏+ and 𝜏− are the mean lifetime of positive and negative muon, respectively, 𝑝𝑑𝑒𝑐𝑎𝑦 is the decay probability of a negative muon in the scintillator and 𝑡𝑑 is the delay time of the arrival of the signal to the START input of the TAC, which is caused by a 10 m cable in the circuit (Fig. 1). To reduce the parameters of the fit, we set the values of these three parameters, 𝜏+ = 2.197 μs, 𝜏− = 2.026 μs, and 𝑝𝑑𝑒𝑐𝑎𝑦 = 92.15% [4], in Eq. (1) and we obtained other parameters by fitting to the data points. So the muon ratio, 𝑅𝜇 = 𝑛+ ∕𝑛− was determined. Table 1 shows the values of 𝑅𝜇 for intervals of one to six weeks. The muon charge ratio value does not change after 3 weeks, which indicates that this period is sufficient for the experiment.

In a decay experiment, only particles that stop inside the detector and decay in it are considered. In the present experiment, muons cross the top scintillator and stop in the bottom scintillator and then decay in it. So, the minimum energy of the particles is related to the particles that enter the detectors vertically, and after passing through the top scintillator, equal to 1 cm, they eventually decay at the beginning of entering the bottom scintillator. The maximum energy also refers to particles that enter from one of the two ends of the top scintillator and pass through the just against end of the bottom scintillator and stop at the longest possible passage in it. The minimum and maximum energies are calculated by using the Bethe–Bloch formula [3]. In this way, the detection energy range is between ∼0.8 GeV and ∼1.6 GeV. The muons entered into the detector are actually produced at high altitudes of the atmosphere. Muons lose their energy by 𝑑𝐸𝜇 ∕𝑑𝑋 = −2 MeV cm2 /g due to ionization in the atmosphere. If we assume

Our experimental results show that the muon charge ratio at Tehran’s level, between energy range 0.8 GeV and 1.6 GeV, is 𝑅𝜇 = 1.18 ± 0.02. 2

M. Bahmanabadi

Nuclear Inst. and Methods in Physics Research, A 945 (2019) 162635

Fig. 3. A description of vertical height of atmosphere, 𝑑, zenith angle, 𝜃, length of a slant trajectory, 𝑙, the radius of the Earth, 𝑅⊕ , and other notations used in the text. Fig. 2. A sample of distribution of the muon decay time. The solid line shows the fitting of Eq. (1).

Table 1 The values of the muon charge ratio in the period of 1 to 6 weeks. Period of time (per week)

1

2

3

4

5

6

𝑅𝜇

1.15 ±0.03

1.17 ±0.03

1.18 ±0.03

1.18 ±0.02

1.18 ±0.02

1.18 ±0.02

5. Positive and negative muons produced from CORSIKA code simulations The simulations of extensive air showers caused by different atomic nuclei, such as protons, alpha, etc., are usually done with the CORSIKA code [5]. In this code, different models are available for interactions between Hadrons. The GHEISHA [6] and UrQMD [7] models for hadronic interactions below 𝐸𝑙𝑎𝑏 = 80 GeV, and the QGSJET (qgsjet01.f package) model [8] for hadronic interactions above 𝐸𝑙𝑎𝑏 = 80 GeV, are used. The observation level of the simulated extensive air showers is considered for the Tehran’s level. The geomagnetic field components for the latitude and longitude of Tehran are: 𝐵𝑥 = 27.97 μT, 𝐵𝑧 = 39.45 μT [9]. Air showers are produced at different zenith angles. So the air showers pass through various geomagnetic field, until eventually reach the surface of Tehran. But despite this, we assume that the geomagnetic field is constant. To prove this assumption, the following equation can be easily obtained in accordance with Fig. 3: 𝑑 ≃ 𝑙 cos 𝜃 +

𝑙2 sin2 𝜃, 2𝑅⊕

(𝑙∕𝑅⊕ ≪ 1).

Fig. 4. Energy distributions of the muons at Tehran’s level, and solid line shows the fitting of Eq. (4).

6 GeV to 100 GeV were chosen for their energy, and in each state, 2 × 104 air showers were simulated. The zenith angles of the primary particles from 0◦ to 60◦ , and all azimuth angles between 0◦ and 360◦ were uniformly taken. The number of positive and negative muons generated from each simulation were counted and the energy of each muon is obtained using the momentum components from the simulation output data as 𝐸𝜇 = (𝑝2 𝑐 2 +𝑚20 𝑐 4 )0.5 , with 𝑝2 = 𝑝2𝑥 +𝑝2𝑦 +𝑝2𝑧 . With a differential flux given by 𝑑𝑁∕𝑑𝐸 ∝ 𝐸 −2.7 for primary particles, we can obtain the energy spectrum of muons as Fig. 4. This distribution is fitted with a function as [10],

(2)

where 𝑑, 𝑙, 𝜃, and 𝑅⊕ are vertical height of atmosphere, length of a slant trajectory, zenith angle, and the radius of the Earth, respectively. For 𝜃 less than 60◦ , the first term is important and can be ignored the second one. If the zenith angle of a muon be 𝜃, the variation in its latitude value from the start point to the end of the path, 𝛥𝜆, is given as: 𝑙 sin 𝛥𝜆 = sin 𝜃. (3) 𝑅⊕ + 𝑑

𝑑𝑁(𝐸𝜇 ) 𝑑𝐸𝜇

=𝑘

(𝐸0 + 𝐸𝜇 )−𝛼 1+

𝐸𝜇

.

(4)

𝜖

where the values of the fit parameters for the Tehran’s level are 𝑘 = 166.8 ± 37.9 (GeV)𝛼−1 , 𝐸0 = 9.8 ± 0.3 GeV, 𝛼 = 2.5 ± 0.1 , and 𝜖 = (8 ± 2) × 103 . This spectrum is flat in energies less than 1 GeV, but it gradually drops.

Since the distribution of muons per solid angle is as 𝑑𝐼∕𝑑𝛺 = 𝐼0 cos2 𝜃, 60◦ the mean value of the zenith angle is ⟨𝜃⟩ = ∫0◦ 𝜃 cos2 𝜃 sin 𝜃𝑑𝜃∕ ◦ 60 ∫0◦ cos2 𝜃 sin 𝜃𝑑𝜃 ≃ 34◦ 34◦ . So for values 𝜃 = 34◦ , 𝑑 = 16 km, 𝑙 = 𝑑∕ cos 𝜃 ≃ 19 km, and 𝑅⊕ = 6400 km, we calculate 𝛥𝜆 = 0.1◦ . With this latitude variation, the particles are exposed different geomagnetic fields. But with the assumption of a dipole approximation for the geomagnetic field, that is 𝐵𝑥 = 𝑟𝑏3 cos 𝜆, and 𝐵𝑧 = − 2𝑏 sin 𝜆, where 𝑟3 for the Earth’s dipole 𝑏 ≃ 8 × 1015 Tm3 , we obtain |𝛥𝐵𝑥 ∕𝐵𝑥 | = 0.001, and |𝛥𝐵𝑧 ∕𝐵𝑧 | = 0.002. Therefore, the variation in total geomagnetic field is: 𝛥𝐵 = [(𝐵𝑥 𝛥𝐵𝑥 )2 + (𝐵𝑧 𝛥𝐵𝑧 )2 ]0.5 ∕𝐵 = 0.07 μT. So the constant geomagnetic field approximation creates a very small error that can be ignored. Since about 90% of the primary particles are protons, only this particle is considered in these simulations. 14 discrete states between

6. The muon charge ratio obtained from simulation data Positive and negative muons are separated based on the particle code in CORSIKA and using the three momentum components of the particle, the muon energy, 𝐸𝜇 , and zenith angle, 𝜃, are calculated as, 𝐸𝜇 = (𝑝2 𝑐 2 + 𝑚20 𝑐 4 )0.5 , cos 𝜃 = 𝑝𝑧 / (𝑝2𝑥 + 𝑝2𝑦 )0.5 . If the muon energy was within the energy range of the telescope’s detection, 0.8 ≤ 𝐸 (GeV) ≤ 1.6, as well as its zenith angle in the range of the telescope’s opening, 𝜃 ≤ 17◦ , that muon would be recorded as a visible event in the telescope. The simulation results for the two sets of models used in this work, QGSJET- GHEISHA and QGSJET-UrQMD, are presented in Tables 2 and 3, respectively. The standard deviation for each muon 3

M. Bahmanabadi

Nuclear Inst. and Methods in Physics Research, A 945 (2019) 162635

Table 2 The number of positive and negative muons with energy range 0.8 ≤ 𝐸 (GeV) ≤ 1.6, and muon charge ratio in different primary particle energies using the models GHEISHA and QGSJET. Primary particle energy (GeV)

𝑛+

𝑛−

𝑅𝜇

6 7 8 9 10 20 30 40 50 60 70 80 90 100

31 59 135 170 306 1 874 4 067 6 520 8 891 11 180 13 635 16 030 20 029 22 275

18 52 93 162 242 1 572 3 233 5 447 7 641 9 692 11 968 13 949 17 429 19 255

1.72 1.13 1.45 1.05 1.26 1.19 1.26 1.20 1.16 1.15 1.14 1.15 1.15 1.16

± ± ± ± ± ± ± ± ± ± ± ± ± ±

UrQMD. The simulation results show that model GHEISHA for low energy hadronic interaction is more consistent with our experimental results than model UrQMD. 7. Conclusion

0.51 0.22 0.20 0.12 0.11 0.04 0.03 0.02 0.02 0.02 0.01 0.01 0.01 0.01

Using a cosmic ray telescope consisting of two scintillators, above each other and one meter apart, the number of positive and negative muons in the muon energy range between 0.8 GeV and 1.6 GeV was obtained by a delayed coincidence method in Tehran’s level. The muon charge ratio value was obtained 𝑅𝜇 = 1.18 ± 0.02 by this experiment. We also obtained the muon charge ratio using the CORSIKA code for the same energy range for muons, and by using two different models GHEISHA and UrQMD for hadronic interactions in low energy , and model QGSJET for the hadronic interactions of high energy. The mean values in the energy range 0.8 ≤ 𝐸 (GeV) ≤ 1.6, are 1.23 ± 0.13 and 1.25 ± 0.16, respectively for two sets of models QGSJET-GHEISHA and QGSJET-UrQMD, which shows that the first model has better compatibility with the result of the experiment.

Table 3 The number of positive and negative muons with energy range 0.8 ≤ 𝐸 (GeV) ≤ 1.6, and muon charge ratio in different primary particle energies using the models UrQMD and QGSJET. Primary particle energy (GeV)

𝑛+

𝑛−

𝑅𝜇

6 7 8 9 10 20 30 40 50 60 70 80 90 100

6 26 50 111 167 1 562 3 568 6 178 9 336 11 979 14 747 17 485 20 008 22 482

7 16 35 83 126 1 268 2 962 5 035 7 687 10 073 12 453 14 692 17 136 19 341

0.86 1.63 1.43 1.34 1.33 1.23 1.20 1.23 1.21 1.19 1.18 1.19 1.17 1.16

± ± ± ± ± ± ± ± ± ± ± ± ± ±

Acknowledgments I would like to thank Miss S. Abdollahi for all her meritorious helps. This work was supported by the office of vice president for science, research and technology of Sharif University of Technology.

0.48 0.52 0.31 0.19 0.16 0.05 0.03 0.02 0.02 0.02 0.01 0.01 0.01 0.01

References [1] [2] [3] [4] [5] [6] [7] [8]

0.5

charge ratio is calculated as 𝜎 = (𝑛2+ 𝑛− + 𝑛2− 𝑛+ ) ∕𝑛2− . Considering the differential energy spectrum proportional to 𝐸 −2.7 for primary particles, the mean number of positive and negative muons is as follows.: 𝑛̄ ± =

14 (𝑛 ) 𝐸 −2.7 𝛴𝑖=1 ± 𝑖 𝑖 14 𝐸 −2.7 𝛴𝑖=1 𝑖

.

(5)

[9] [10]

Hence, the mean values of the muon charge ratios, in the energy range between 0.8 GeV and 1.6 GeV, are 1.23 ± 0.13 and 1.25 ± 0.16, respectively for two sets of models QGSJET-GHEISHA and QGSJET-

4

T.K. Gaisser, T. Stanev, Phys. Rev. D 38 (1998) 85. B. Vulpescu, et al., Nucl. Instrum. Methods A 414 (1998) 205. C. Grupen, Astroparticle Physics, Springer-Verlag, Berlin, Heidelberg, 2005. T. Suzuki, D.F. Measday, J.P. Roalsvig, Phys. Rev. C 35 (1987) 2212. D. Heck, et al., Report FZKA6019, Forschungszentrum Karlsruhe, 1998. H. Fesefeldt, (GHEISHA), RWTH Aachen, Germany Report No. PITHA-85/02, 1985 (unpublished). S.A. Bass, et al., Prog. Part. Nucl. Phys. 41 (1998) 225. N.N. Kalmykov, S.S. Ostapchenko, Yad. Fiz 56 (1993) 105; Phys. At. Nucl. 56 N3 (1993) 346. N.N. Kalmykov, S.S. Ostapchenko, A.I. Pavlov, Izv. RAN Ser. Fiz. 58 (N12) (1994) 21; N.N. Kalmykov, S.S. Ostapchenko, A.I. Pavlov, Bull. Russ. Acad. Sci. (Phys.) 58 (1994) 1966; N.N. Kalmykov, S.S. Ostapchenko, A.I. Pavlov, Nucl. Phys. B (Proc. Suppl.) 52B (1997) 17. https://www.ngdc.noaa.gov/geomag-web/. M. Bahmanabadi, L. Rafezi, Phys. Rev. D 98 (2018) 103003.