Reconstruction of a muon bundle in the JUNO central detector

Reconstruction of a muon bundle in the JUNO central detector 1 first-author 1 corresp Yong-Bo Huang, huangyb@gxu.edu.cn Yong-Bo Huang, huangyb@gxu.edu.cn huangyb@gxu.edu.cn https://orcid.org/0000-0002-3801-9749 2 corresp Ji-Lei Xu, xujl@ihep.ac.cn Ji-Lei Xu, xujl@ihep.ac.cn xujl@ihep.ac.cn 2 3 2 2 https://orcid.org/0000-0003-0137-1797 2 1 1 School of Physical Science and Technology, Guangxi University, Nanning 530004, China School of Physical Science and Technology, Guangxi University, Nanning 530004, China Institute of High Energy Physics, Beijing 100049, China Institute of High Energy Physics, Beijing 100049, China University of Chinese Academy of Sciences, Beijing 100049, China University of Chinese Academy of Sciences, Beijing 100049, China The Jiangmen Underground Neutrino Observatory (JUNO) is a multi-purpose neutrino experiment. One of the main goals is to determine the neutrino mass ordering by precisely measuring the energy spectrum of reactor antineutrinos. For the detection of reactor antineutrinos, cosmogenic backgrounds, such as 9 Li/ 8 He and fast neutrons induced by cosmic muons, should be rejected carefully by applying muon veto cuts, which require good muon track reconstruction. With a 20-kton liquid scintillator detector, the simulation shows the proportion of muon bundles (muon multiplicity ≥ 2) to be approximately 8% in JUNO, whereas its reconstruction has been rarely discussed in previous experiments. This study proposes an efficient algorithm for muon track reconstruction based on the charge response of a photomultiplier tube array. This is the first reconstruction of muon bundles in a large-volume liquid scintillator detector. In addition, the algorithm shows good performance and potential for reconstruction for both a single muon and double muons (muon multiplicity = 2). The spatial resolution of a single-muon reconstruction was 20 cm, and the angular resolution was 0.5°. For double-muon reconstruction, the spatial and angular resolutions could be 30 cm and 1.0°, respectively. Moreover, this paper also discusses muon classification and the veto strategy. JUNO Liquid scintillator detector Muon reconstruction Muon bundle Veto strategy Introduction 1 Long-lived radioactive isotopes produced by energetic cosmic muons are one of the main sources of large-volume neutrino and dark matter detectors [ 1 1 - 9 9 1-9 ]. Locating the detector deep underground to provide an essential overburden is an effective approach for suppressing muon-induced backgrounds. The Jiangmen Underground Neutrino Observatory (JUNO) [ 3 3 , 10 10 , 11 11 ] is a multi-purpose neutrino experiment with a detector approximately 700 m underground. One of the main goals is to determine the neutrino mass ordering by precisely measuring the energy spectrum of reactor antineutrinos at a site approximately 53 km from the reactors of the Yangjiang and Taishan nuclear power plants. Figure 1 Figure 1 shows a schematic of the JUNO detector. The central detector (CD) is a 35.4-m-diameter sphere with a 20-kton liquid scintillator detector inside that uses 17,612 20-inch photomultiplier tubes (PMTs) and 25,600 3-inch PMTs as photosensors. The water Cerenkov detector is a cylindrical water pool with 2400 20-inch PMTs filled with 35 kton of ultrapure water as a muon veto and radioactive shielding. The top-tracker detector can be used to tag muon events. In the muon simulation using the JUNO offline software [ 12 12 , 13 13 ], the cosmic air shower simulation at the JUNO site, local geological map, and details of the detector structure were considered, resulting in a muon rate of approximately 0.004 Hz/m 2 with an average muon energy of approximately 207 GeV [ 10 10 ]. The simulation shows that there are approximately 92% single-muon events, approximately 6% double-muon events, and less than 2% of other events with muon multiplicity of three or more. The muon energy deposit ranged from 0 to approximately 100 GeV. For a giant volume detector such as JUNO, the muon rate ( ∼ 3.4 × 10 5 /day) is still much higher than the reactor antineutrino rate ( ∼ 83 inverse beta decay (IBD) events caused by the reactor per day before the selection cut). According to a preliminary antineutrino selection criterion [ 3 3 ], there are approximately 71 residual cosmogenic 9 Li/ 8 He events per day, which is comparable to IBD events (73 per day) after the antineutrino selection cut. Thus, the β-n decay from 9 Li/ 8 He is the most seriously correlated background to reactor antineutrinos and needs to be further rejected with an efficient muon veto strategy. Vetoing the entire volume for some time after a muon event is both uneconomical and would result in too much exposure. Because most cosmogenic isotopes are generated near the muon track, a strategy that vetoes a sufficient cylindrical volume along the muon trajectory has been adopted, which has been proven to be a good approach [ 3 3 ]. For low-background experiments, such as neutrino and dark matter experiments, muon-induced background suppression is the main motivation for muon reconstruction. There are several types of muon tracks in the JUNO spherical detector. According to the position relative to the detector and the behavior of the muon tracks, the tracks can be classified into the following categories, as shown in Fig. 2 Fig. 2 . (1) Through-going muon: This corresponds to a muon that goes through the CD and the distance between the muon track and the center of the CD is less than 16 m. (2) Clipping muon: This corresponds to a muon that passes through the CD but leaves only a short track in the detector (in other words, the distance between the muon track and the center of the CD is larger than 16 m). (3) Stopping muon: This corresponds to a muon track that has stopped in the CD. In addition, based on muon multiplicity, muons can be classified into the following categories: single, double, and multiple ( Fig. 2 Fig. 2 ). Muon bundles [ 14 14 - 18 18 14-18 ], corresponding to muons with multiplicity of two or more, are remnants of air showers produced in the atmosphere by high-energy cosmic-ray nuclei. When multiple muons enter an underground detector simultaneously, their tracks are approximately parallel because they are produced in the same air shower, which typically occurs kilometers away from the detector [ 16 16 ]. In previous experiments, the reconstruction of single muons has been well studied in both water/ice Cherenkov detectors [ 19 19 - 21 21 19-21 ] and liquid scintillator detectors [ 22 22 , 23 23 ]; however, the reconstruction of muon bundles has only been studied in large water/ice Cherenkov detectors [ 24 24 - 26 26 24-26 ]. In JUNO, three reconstruction algorithms have been developed for single muons, such as a method with a geometrical model which utilizes the geometrical shape of the fastest light [ 27 27 ], a method with the fastest light model which utilizes the minimization of the first hit time (FHT) [ 28 28 ], and new technology using deep learning and GPU acceleration [ 30 30 ]. However, reconstruction of muon bundles was not included. The time needed for the first two algorithms is at the level of a few seconds. For experiments with detectors of tens of meters, such as JUNO and LENA [ 31 31 ], the proportion of muon bundles is higher. Hence, it is necessary to develop an algorithm to reconstruct multiple muon tracks. This paper proposes an efficient algorithm for the track reconstruction of both single and double muons which also de monstrates good reconstruction abilities for clipping muons based on experiments at JUNO. However, a few muon tracks are accompanied by electromagnetic or hadronic showers, called shower muons, and it is difficult to estimate where the showers take place. A dedicated algorithm for shower muons is being developed; however, it is not the subject of this study. Charge pattern of the PMT hit array 1 When a cosmic muon passes through the JUNO detector, the muon leaves track information in the detector, including time and charge information. For example, we can observe a PMT charge cluster at the incident point and another charge cluster at the exit point. To deal more easily with the charge clusters, the PMT hit array was projected onto a two-dimensional (2D) plane to show the charge pattern ( Fig. 3 Fig. 3 ). In the charge pattern, the number of PMTs in each pixel varied from zero to eight depending on their locations. The fastest light [ 27 27 - 29 29 27-29 ] information around the muon track is widely used for single-muon reconstruction, but it is more complicated to handle in the case of multiple muon tracks. To reconstruct muon bundles, the charge pattern of the PMT hit array was investigated, and different types of muons can provide different charge clusters ( Fig. 3 Fig. 3 ). The 20-inch PMTs near the muon incident point and the exit point receive several thousands of photoelectrons (PEs); they are expected to be saturated and are not applied in this analysis. The 3-inch PMTs, with a photocathode area approximately 40 times smaller than that of the large PMTs, work in a larger dynamic range and are suitable for this study [ 13 13 , 32 32 ]. Different types of muons have different charge pattern cluster features. The cluster number is related to the muon track number and track behavior, which are listed in Table 1 Table 1 . A single through-going muon event ( Fig. 3 Fig. 3 ) has two clusters but a single clipping-muon event ( Fig. 3 Fig. 3 ) has only one cluster, because its incident and exit points are very close. A single stopping-muon event ( Fig. 3 Fig. 3 ) has only one cluster because it has only one incident point, but not an exit point. A through-going muon event ( Fig. 3 Fig. 3 ) has four clusters that correspond to two incident and two exit points. There might be one or both tracks clipping in double muons, so the number of clusters is three or two. Similarly, a stopping double-muon event includes three clusters (one muon stops) or two clusters (both muon stop). For muons with multiplicity of three or more (multiple muons), their reconstruction corresponds to a low reconstruction efficiency and is not discussed further in this paper. 1 1   1 1 Through-going 1 1 Clipping 1 1 Stopping 1 1 Single muon 1 1 2 clusters 1 1 1 cluster 1 1 1 cluster 1 1 Double muons 1 1 4 clusters 1 1 2–3 clusters 1 1 2–3 clusters 1 1 Multiple muons 1 1 >4 clusters 1 1 >2 clusters 1 1 >2 clusters Reconstruction algorithm 1 Using information from the PMT charge pattern, a reconstruction algorithm was developed for both single and double muons. The following is a brief step list of the algorithm flow ( Fig. 4 Fig. 4 ). More details are introduced later. (1) Project the PMT charge pattern onto the θ – ϕ plane. (2) Discard PMTs with charges smaller than 75 PEs, further details are given in Fig. 5 Fig. 5 . (3) Apply the charge smoothing method (Sect. 3.2) to reduce the number of fake clusters found when applying the fitting in the next step. (4) Fit with the ROOT [ 33 33 ] tool “TSpectrum2” [ 34 34 ] to determine the center positions of clusters (cluster seeds) for all possible cluster candidates. (5) Rotate the PMT ball by 90° so that the PMTs close to the north and south poles are moved to the equator (Sect. 3.2) and then repeat steps (1)–(4). When cluster seeds are found, their spatial coordinates are converted back to real coordinates before rotation. This rotation operation is applied to avoid distortions when projecting PMTs from a spherical surface to a 2D plane, especially in the area near the north and south pole. (6) Charge-weighted calculation to correct each cluster seed and discard fake seeds. (7) Adjacent candidate seeds are merged if they are too close. (8) We matched two candidate seeds as a muon track. An additional strategy is used to handle cases where there are more than two candidate seeds. More details can be found in Sect. 3.3. Methodologies 2 For the JUNO detector, the PMTs closely packed around the LS sphere form a PMT shell, which is a 2D surface instead of a 3D surface, so the PMT charge hit pattern can be projected onto a 2D plane to perform convenient fitting. Similar to the world map of the Earth, there are several projection methods to obtain a 2D map, such as the longitude–latitude projection method and the spherical projection method, but every method introduces distortion through the projection. For reconstruction, different projection methods yielded similar reconstruction results. Finally, we selected the longitude–latitude projection method. The ϕ angle of the PMT position is for the X -axis on the 2D plane and the θ angle is for the Y -axis. When a muon passes through the LS, more than 90% of the photons are emitted from the scintillation process, which emits isotropic light. The closest PMTs to the muon incident point or exit point will collect the most PEs, which will create two clusters on the PMT charge-hit pattern ( Fig. 3 Fig. 3 ). Figure 5 Figure 5 shows an example of the distribution of the number of photoelectrons (nPE) for a single through-going-muon event and a double through-going-muon event. PMTs with a charge greater than 75 PEs generally correspond to clusters caused by muon tracks. To highlight the clusters that were caused by muon incident and exit points, we needed to eliminate the PMT noise ( ∼ 1 PE), radioactivity (several PEs), and PMTs far from clusters (tens of PEs). As a result, PMTs with charges smaller than 75 PEs are discarded. In a simulation study using the JUNO offline software [ 12 12 , 13 13 ], when the maximum discard charge varies from 70 to 80 PEs, the reconstruction result is not changed significantly; more details can be found later (Section 4.2 and Fig. 14 Fig. 14 ). For the PMT charge pattern, a 2D peak searching and fitting algorithm named “TSpectrum2” was investigated and applied. “TSpectrum2” is a ROOT tool [ 34 34 ], which is based on a 2D Gaussian function fitting and it searches for peaks in source spectrum with a deconvolution method. “TSpectrum2” can obtain the maximum local peak by calculating the spectrum of the standard deviations of the smoothed second derivatives in two dimensions [ 35 35 ]. If the fit is successful, the number of identified peaks and their locations are returned. The identified peaks were arranged according to their heights in the spectrum. In our study, the maximum number of clusters was set to seven to list the highest seven peaks only, which is sufficient for searching for single muons (2 peaks), double muons (4 peaks), and triple muons (6 peaks). The tracking precision of triple muons is much lower than that of double muons, so it is not described in this paper, but we still set the maximum to seven to retain the reconstruction potential for triple muons when using this method, which will be studied further in the future. “TSpectrum2” runs very fast (Sect. 4.3) and it can mathematically find the cluster seed candidates. However, the “TSpectrum2” fitting usually finds more seeds than at the muon incident and exit points. The unmatched seeds were caused by the statistical fluctuation of the PMT charge. In the simulation, the collected charges of adjacent PMTs surrounding the largest charge collection PMT do not continuously drop from the cluster center to the cluster edge. However, some muons have a relatively large energy deposit at a certain position along the track. All of these may result in an incorrect cluster number and fake cluster seeds when matching with the true muon hit points. Figure 6a Figure 6a shows an example of the PMT charge hit pattern of a single through-going-muon event before charge smoothing above the PMT charge threshold. This shows that the fitting finds three cluster seeds near the muon exit point, where we need only one cluster seed. The fitting results for the cluster seeds must be carefully checked and analyzed. To discard fake cluster seeds and obtain true hit points, a charge-weighted algorithm was developed. First, based on each cluster seed, we calculated the charge-weighted position of the fired PMTs in a 4-m radius (corresponding to the size of the normal cluster) using the following formula: <math display="block" id="M1"><mrow><mover accent="true"><mi>x</mi><mo>¯</mo></mover><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mi>n</mi></msubsup><mrow><msub><mi>q</mi><mi>i</mi></msub><msub><mi>x</mi><mi>i</mi></msub></mrow></mstyle></mrow><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mi>n</mi></msubsup><mrow><msub><mi>q</mi><mi>i</mi></msub></mrow></mstyle></mrow></mfrac></mrow></math> where x i is the spatial position of PMT i , n is the number of fired PMTs in a 4-m radius around the cluster seed, and q i is the charge collected by PMT i . If the total charge in a 4-m radius around the cluster seed is less than 1000 PEs, the cluster seed is discarded. Second, after the calculation, replace the cluster seed with the charge-weighted position ( <math id="M2"><mover accent="true"><mi>x</mi><mo>¯</mo></mover></math> ). Then, if two new cluster seeds are too close ( < 4 m), they are merged into one point according to their charge weights. This step can filter out most of the fake seeds given by “TSpectrum2” fitting and give a correct position. Third, if the gathered cluster charge within a radius of 4 m around the cluster seed is too small ( < 10% of the total charge of all clusters), the seed is also discarded. This can eliminate the small-charge clusters found by “TSpectrum2” fitting which cannot match with the true incident and exit points of the muons. Finally, the charge-weighted positions of the remaining clusters were considered as muon incident and exit points. Distortion elimination and performance optimization 2 To make the muon track reconstruction more precise, we eliminated distortion using the rotation method, removed the statistical fluctuation of the charge cluster using the smoothing method, and transferred the cluster seeds on the PMT charge pattern onto the true muon track. As mentioned previously, the projection method (normal projection) causes non-negligible distortion, particularly near the polar zone. To avoid projection distortion, we rotated the PMT positions 90° around the X -axis to allow the north and south poles to be placed at the equator ( Fig. 7 Fig. 7 ), and then, the new PMT positions are again projected onto θ – ϕ 2D plane (called rotation- X projection). The cluster seeds in the θ range (-50°, 50°) were selected in both the normal projection (before the rotation) and rotation- X projection (after the rotation). Some repeated seeds in the overlap area of the two types of projection are merged during the charge-weighted step. Sometimes, some clusters are split from one cluster into two clusters when the incident or exit points are located at the edge of the 2D plane ( ϕ > 160° or ϕ < -160°). First, the PMTs were rotated 90° around the Z -axis (rotation- Z projection) to obtain the precise position of the clusters that were located at the edge of the 2D plane. Then, combined with the normal projection and rotation- Z projection, we check whether some fitted points are located at the edge. If a cluster is split into two parts, the fitting result from the rotation- Z projection is used. In addition, apart from the peaks induced by the incident and exit points of the muon, there were also some relatively small peaks. This is because of the statistical fluctuation in the PE number collected by the PMTs or a relatively large-energy deposit at a certain location along the muon track. Figure 6a Figure 6a shows all the reconstructed peaks. Meanwhile, the “TSpectrum2” tool may supply fitted points at the non-hit position. To handle this effect, each pixel's charge was smoothed ( q s ) based on the adjacent pixel charges using the smoothing formula <math display="block" id="M3"><mrow><msup><mi>q</mi><mi>s</mi></msup><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mi>N</mi></msubsup><mrow><mrow><mo>(</mo><mrow><mfrac><mrow><msub><mi>w</mi><mn>0</mn></msub></mrow><mrow><msub><mi>w</mi><mn>1</mn></msub><mo>+</mo><msub><mi>d</mi><mi>i</mi></msub></mrow></mfrac><msub><mi>q</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow></mrow></mstyle></mrow><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mi>N</mi></msubsup><mrow><mrow><mo>(</mo><mrow><mfrac><mrow><msub><mi>w</mi><mn>0</mn></msub></mrow><mrow><msub><mi>w</mi><mn>1</mn></msub><mo>+</mo><msub><mi>d</mi><mi>i</mi></msub></mrow></mfrac></mrow><mo>)</mo></mrow></mrow></mstyle></mrow></mfrac></mrow></math> where q i is the charge of the i th adjacent pixel, w 0 and w 1 are the weight factors, with the best tuning result w 0 / w 1 = 4/9, and d i is the distance factor between the smoothing pixel and the i th adjacent pixel. The distance factor of the nearest pixel was one. The further the pixel distance, the larger d i is and the lower the weight of q i . Here N is the number of adjacent pixels that are used for charge smoothing and surround the specified pixel, which is set to 24. A comparison of Figs. 6a and b Figs. 6a and b shows that the charge smoothing method can make the cluster peaks protrude further, reduce the number of false fitting peaks, and can also improve the reconstruction efficiency. In fact, the charge cluster center is not exactly equal to the true incident and exit points of the muon. This is because there is a water buffer with a thickness of 1.8 m between the PMT spherical shell and the LS sphere ( Fig. 8 Fig. 8 ). When a muon enters the PMT shell, most photons originate from the LS, and the PMT charge cluster is close to the muon hit point of the LS. As a result, the track reconstructed by the PMT charge cluster center does not correspond to the true muon track. For this reason, the cluster charge centers (yellow clusters in Fig. 8 Fig. 8 ) must be translated into true muon hit points on the LS sphere (red clusters in Fig. 8 Fig. 8 ), according to their spatial positions. Finally, we reconstructed the muon track (green dotted line) using the reconstructed incident and exit points. Double-muon event reconstruction strategies 2 As listed in Table 1 Table 1 , a double-muon event typically corresponds to four hit points. Every two points comprise a track candidate. There are three combinations of four hit points, and two of them are in parallel; thus, they can be considered as muon bundle candidates, as shown in Fig. 7(a) Fig. 7(a) . The PMT timing information was used to further identify the true muon tracks. With each reconstructed point in the center, the FHTs of the surrounding PMTs within 4 m were averaged as the cluster time. The two earlier times were treated as two incident points, and the other two points were treated as two exit points. This time strategy can help to determine the correct combination. If the two tracks are too close, then only two cluster seeds can be found ( Fig. 9(b) Fig. 9(b) ). In this case, the two cluster seeds are connected and reconstructed into a single track. Although the two tracks are close, the overlapped cluster is larger than the normal cluster. The distance between two muon tracks can be deduced from the overlapped cluster information. The primary requirement of muon track reconstruction is to veto muon-induced isotopes through a cylindrical volume cut along the muon track. Although only one track is reconstructed, a cylindrical cut with a larger radius will be considered and can provide the same background reduction ability. If double muons have only three charge clusters, it is not easy to distinguish between the two muon tracks. In this study, we identified three cases in a three-cluster situation ( Figs. 9(c)–(e) Figs. 9(c)–(e) ). One case occurred when two incident points overlapped ( Fig. 9(c) Fig. 9(c) ). We calculated the average FHT ( <math id="M4"><mover accent="true"><mi>T</mi><mo>¯</mo></mover></math> ) of the surrounding PMTs within 4 m of each cluster. Because two tracks in the double muons shoot in the LS at the same time, the <math id="M5"><mover accent="true"><mi>T</mi><mo>¯</mo></mover></math> of the incident or exit points are close. In three values of <math id="M6"><mover accent="true"><mi>T</mi><mo>¯</mo></mover></math> , if the two <math id="M7"><mover accent="true"><mi>T</mi><mo>¯</mo></mover></math> differences are less than 5 ns, both clusters are recognized as exit points and the other is the incident point. Then, the two muons can be correctly tracked using the incident point twice. Although the incident point is not very precise, the larger radius cylindrical cut according to the incident cluster size along each track can still provide the same background reduction of the muon-induced isotopes. The other two cases, clipping ( Fig. 9(d) Fig. 9(d) ) or stops in the CD ( Fig. 9(e) Fig. 9(e) ), also produced three clusters. The projection on the 2D charge histogram is shown in Fig. 10 Fig. 10 , respectively. Although the two clusters' time ( <math id="M8"><mover accent="true"><mi>T</mi><mo>¯</mo></mover></math> ) difference is larger than 5 ns, we still need to select two seeds to make a track. A direction determination method was developed for the 2D cluster plot. When a muon track in the LS is further from the CD center, such as a clipping muon, the PMT charge cluster is shaped like an ellipsoid. The direction of the elliptical major axis was along the muon track direction. The major axis directions of the ellipse can be expressed by the vectors <math id="M9"><mrow><mover accent="true"><mrow><msub><mi>V</mi><mi>A</mi></msub></mrow><mo stretchy="true">←</mo></mover></mrow></math> , <math id="M10"><mrow><mover accent="true"><mrow><msub><mi>V</mi><mi>B</mi></msub></mrow><mo stretchy="true">←</mo></mover></mrow></math> , and <math id="M11"><mrow><mover accent="true"><mrow><msub><mi>V</mi><mi>C</mi></msub></mrow><mo stretchy="true">←</mo></mover></mrow></math> for the three clusters ( Fig. 10 Fig. 10 ). The angles of the candidate track ( L AB , L BC , or L AC ) with its corresponding cluster vectors can be calculated using the following equations: <math display="block" id="M12"><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow><mo>∠</mo><mi> </mi><mn>1</mn><mo>=</mo><mo>∠</mo><mi> </mi><mover accent="true"><mrow><msub><mi>V</mi><mi>A</mi></msub></mrow><mo stretchy="true">←</mo></mover><mo>&</mo><msub><mi>L</mi><mrow><mi>A</mi><mi>B</mi></mrow></msub><mo>+</mo><mo>∠</mo><mi> </mi><mover accent="true"><mrow><msub><mi>V</mi><mi>B</mi></msub></mrow><mo stretchy="true">←</mo></mover><mo>&</mo><msub><mi>L</mi><mrow><mi>A</mi><mi>B</mi></mrow></msub></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow><mo>∠</mo><mi> </mi><mn>2</mn><mo>=</mo><mo>∠</mo><mi> </mi><mover accent="true"><mrow><msub><mi>V</mi><mi>B</mi></msub></mrow><mo stretchy="true">←</mo></mover><mo>&</mo><msub><mi>L</mi><mrow><mi>B</mi><mi>C</mi></mrow></msub><mo>+</mo><mo>∠</mo><mi> </mi><mover accent="true"><mrow><msub><mi>V</mi><mi>C</mi></msub></mrow><mo stretchy="true">←</mo></mover><mo>&</mo><msub><mi>L</mi><mrow><mi>B</mi><mi>C</mi></mrow></msub></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow><mo>∠</mo><mi> </mi><mn>3</mn><mo>=</mo><mo>∠</mo><mi> </mi><mover accent="true"><mrow><msub><mi>V</mi><mi>A</mi></msub></mrow><mo stretchy="true">←</mo></mover><mo>&</mo><msub><mi>L</mi><mrow><mi>A</mi><mi>C</mi></mrow></msub><mo>+</mo><mo>∠</mo><mi> </mi><mover accent="true"><mrow><msub><mi>V</mi><mi>C</mi></msub></mrow><mo stretchy="true">←</mo></mover><mo>&</mo><msub><mi>L</mi><mrow><mi>A</mi><mi>C</mi></mrow></msub></mrow></mtd></mtr></mtable></mrow></math> where <math id="M13"><mrow><mo>∠</mo><mi> </mi><mover accent="true"><mrow><msub><mi>V</mi><mi>A</mi></msub></mrow><mo stretchy="true">←</mo></mover><mo>&</mo><msub><mi>L</mi><mrow><mi>A</mi><mi>B</mi></mrow></msub></mrow></math> is the acute angle of <math id="M14"><mrow><mover accent="true"><mrow><msub><mi>V</mi><mi>A</mi></msub></mrow><mo stretchy="true">←</mo></mover></mrow></math> and line L AB , as are the other angles. The minimum of three angles <math id="M15"><mrow><mo>∠</mo><mi> </mi><mn>1</mn></mrow></math> , <math id="M16"><mrow><mo>∠</mo><mn>2</mn></mrow></math> , and <math id="M17"><mrow><mo>∠</mo><mi> </mi><mn>3</mn></mrow></math> is selected, and the corresponding two clusters are set to one muon track, and the line direction is selected as the muon direction. A single-muon event may also have three clusters, which may confuse double-muon reconstruction. When a single-muon track is close to the edge of the LS sphere (the spatial distance between the muon track and LS center is from 13 to 16 m), the photon emissions on the muon track are all close to the PMTs and almost all PMTs are fired ( Fig. 10 Fig. 10 ). For this case, we still calculated the corresponding three angle values using Equation (3). The main feature is that all three angles are less than 5°. Therefore, we used this feature to distinguish between single and double muons. The midpoint ( Fig. 10 Fig. 10 ) is a false point, and the remaining two points are selected as incident and exit points. Reconstruction performance and discussion 1 To optimize and validate the reconstruction algorithm in this study, approximately 380,000 simulated muon samples were produced using the JUNO offline software [ 12 12 , 13 13 ]. The simulated muon samples included approximately 92% single-muon events, approximately 6% double-muon events, and less than 2% multiple-muon events. The reconstruction resolution and accuracy of the muon track form the basis for the muon veto strategy and background rejection. In addition, the reconstruction efficiency, robustness, and speed should be considered to evaluate the performance of the reconstruction algorithm, especially for experiments with high muon rate and large-volume detectors. Reconstruction performances of single and double muons 2 To obtain further exposure for the signal (neutrinos or others), the basic veto strategy is to apply a partial volume veto along the muon track for a period of time. The reconstructed and true track relative positions can be described using two parameters, Δ D and α . Here Δ D is defined as Δ D = D t - D r , where D r and D t are the distances between the reconstructed and true tracks from the center of the CD. We use α to denote the acute angle between two tracks ( Fig. 11 Fig. 11 ). The distributions of Δ D and α indicate the reconstruction performance. The single and double through-going-muon reconstruction performance at a distance of 16 m from the detector center is shown in Fig. 12 Fig. 12 . The distance to the center is divided into eight slices with a 2-m interval, which indicates a muon passing through different detector volumes. For a single muon, the average biases of Δ D in all slices (black squares in Fig. 12 Fig. 12 ) were less than 5 cm, and the average biases of α (the black squares in Fig. 12 Fig. 12 ) increase from 0.5° to 1.3° from the detector center to the edge. To reconstruct the double muons, as shown in Figs. 12c and d Figs. 12c and d , the average biases of Δ D increased from 2 to 15 cm, and the average biases of α increased from 0.6° to 1.5°. The reconstruction resolution for each slice is shown in Fig. 13 Fig. 13 . In the internal volume ( R ≤ 10 m), the resolution of Δ D and α is better than 20 cm and 0.5° for a single muon, and better than 30 cm and 0.6° for double muons. In the case of a muon passing through an external volume (10舁m ≤ R < 16 m), the resolution of Δ D and α increases to approximately 22 cm and 0.9° for a single muon and to approximately 35 cm and 1.2° for double muons, respectively. This is because the track length is short near the edge, and the corresponding clusters of muon incident and exit points may overlap; thus, there will be more deviations in the search for cluster seeds based on the charge-weight method. As mentioned in Sect. 2, a clipping muon (16舁m ≤ R < 17.7 m) creates a short track, and the two clusters easily overlap into one cluster; therefore, they are difficult to reconstruct into one track, so they cannot be evaluated by the parameters Δ D and α . Reconstruction efficiency and muon bundle veto strategy discussion 2 As mentioned in Sect. 3, to match two candidate seeds to a muon track, several optimization methods were applied, especially in the case of more than two cluster seeds. As a result, most muon tracks can be reconstructed well. However, there are still some cases in which the muon track cannot be reconstructed (reconstructed track unavailable) or the reconstructed track is far from the true track (Δ D or α is larger than the five standard deviations of their distributions). For a single and double through-going muons, the percentages of these two cases are 1.6% and 1.0% and 1.7% and 3%, respectively. After reconstruction, reconstruction efficiency can be defined as the ratio of successful reconstruction to the total number of events. Figure 14a Figure 14a shows the reconstruction efficiency of single and double through-going muons. A single through-going muon with high reconstruction efficiency has an average efficiency of 98.3% ± 0.1%. The average reconstruction efficiency of the double through-going muons was 95.5% ± 0.2%. In addition, as shown in Fig. 14b Fig. 14b , the reconstruction efficiencies of single and double through-going muons are stable when the PE cut varies from 70 PEs to 80 PEs. A muon veto strategy can be developed based on the reconstructed muon track. A simple approach is to veto a cylindrical volume along the muon trajectory in a given time period [ 3 3 ]. Considering the position dependence of the reconstruction bias and resolution (Sect. 4.1), the radius of the cylinder can be modified by increasing the reconstruction bias and resolution [ 27 27 ]. Furthermore, a distance-dependent veto strategy was introduced for the detection of 8 B solar neutrinos at JUNO [ 36 36 ], which investigates the time and distance distributions between the isotopes and muon tracks in detail and can significantly improve the signal-to-background-noise ratio. Considering the non-negligible contribution from the muon bundle and the higher estimated muon rate caused by the shallower overburden in JUNO, the veto strategy of JUNO must be optimized to obtain more exposure. Using the algorithm in this study, we can reconstruct different types of muons and develop their corresponding veto strategies in the future. For a muon event with an available reconstructed track, a distance-dependent cylindrical veto cut along the reconstructed track is applied. The rejected detector volume and veto time can be determined according to the space and time distributions from the muon track to 9 Li/ 8 He [ 36 36 ]. If a double-muon event is reconstructed with two nearby reconstructed tracks (track distance d < 3 m), a larger cylindrical radius is used in the above strategy. For a clipping or stopping muon, their tracks are difficult to reconstruct because of the heavy smearing of clusters; thus, a spherical veto using the cluster seed as the center are used for conservative exclusion. For multiple-muon events, the cluster features of the charge pattern are more complicated, and it is difficult to match them to the correct muon tracks. As a result, the tracking precision of multiple muons was much lower than that of single and double muons, so it is not shown in this paper. In reconstruction, if a muon event is tagged as multiple muons, a whole-volume veto is applied. In addition, for failed reconstruction events whose track or cluster seeds are unavailable, we also veto a whole volume. However, a new method (called neutron veto) is under development that can use the corresponding neutron to precisely determine the position of muon-induced isotopes. A neutron veto uses a spherical veto on neutrons due to isotope and neutron space correlations. This is expected to save more dead volume and allow a longer veto time window for long-lived isotope rejection. On the other hand, the neutron veto is also expected to have a good ability to reject 9 Li/ 8 He from the shower muon. Thus, in addition to the muon veto, the neutron veto can be used to extend the veto strategy in the future. Reconstruction speed performance 2 The time distribution of muon reconstruction is shown in Fig. 15 Fig. 15 , which was approximately 98 ms/event. The CPU (Intel Xeon E5-2660 v4) implementation of the reconstruction method in this study provides a significant processing speed improvement by a factor of 51 compared with the fastest light method (5000 ms/event) [ 28 28 , 30 30 ]. This speed makes it possible to apply the algorithm in online event classification, which requires fast online filtering for possible data rate reduction. Summary 1 This is the first reconstruction of muon bundles in a large-volume liquid scintillator detector. We developed a muon track reconstruction algorithm that mainly focuses on double muons, and it also shows good performance and potential for the reconstruction of a single muon. The algorithm was verified using the JUNO underground muon samples simulated by Geant4 in JUNO offline software. For the spherical geometry of JUNO's LS target, different muons have different characteristics. Muon events are classified into single-, double-, and multiple-muon categories, and each category has through-going, clipping, and stopping features. These types of muons were considered and reconstructed, which is similar to the real situation in the future when recording the data. The muon track incident and exit clusters on the PMT charge-hit pattern can be roughly searched, and their centers can be fitted using the ROOT tool “TSpectrum2.” However, the obtained clusters must be processed further and corrected. In addition to repeated charge weighting based on the surrounding PMTs, several optimization methods, such as rotation projection, charge smoothing, and geometry correction, have been developed. A complicated clustering strategy for double muons has been discussed in detail. The reconstructed results indicate that the algorithm can reconstruct muon tracks with a resolution of 20 cm from the detector center, 0.5° in angle, and a tracking efficiency of 98.3% ± 0.1% for single-muon events. For double-muon events, the resolutions of the distance to the detector center and angle were approximately 30 cm and 1.0°, respectively, and the tracking efficiency was 95.5% ± 0.2%. Finally, we have discussed an optimized veto strategy that includes a distance-dependent cylindrical veto and spherical veto for different muons, and it is expected to obtain more exposure. The reconstruction speed is faster than that of the reconstruction methods that use the fastest light information. This advantage can promote the deployment of the algorithm in online event classification, which requires fast calculation speed to classify events to perform fast online filtering for a possible data rate reduction. References F.P. An , J.Z. Bai , A.B. Balantekin et al. , Observation of electron-antineutrino disappearance at Daya Bay . Phys. Rev. Lett. 108 , 171803 ( 2012 ). doi: 10.1103/PhysRevLett.108.171803 http://doi.org/10.1103/PhysRevLett.108.171803 Observation of electron-antineutrino disappearance at Daya Bay G. Bellini , J. Benziger , D. Bick et al. , Cosmogenic Backgrounds in Borexino at 3800 m water-equivalent depth . JCAP 08 , 049 ( 2013 ). doi: 10.1088/1475-7516/2013/08/049 http://doi.org/10.1088/1475-7516/2013/08/049 Cosmogenic Backgrounds in Borexino at 3800 m water-equivalent depth F.P. An , G.P An , Q. An et al. , Neutrino Physics with JUNO . J. Phys. G 43 , no. 3 , 030401 ( 2016 ). doi: 10.1088/0954-3899/43/3/030401 http://doi.org/10.1088/0954-3899/43/3/030401 Neutrino Physics with JUNO S. Locke , A. Coffani , K. Abe et al. , New Methods and Simulations for Cosmogenic Induced Spallation Removal in Super-Kamiokande-IV . doi: 10.48550/arXiv.2112.00092 http://doi.org/10.48550/arXiv.2112.00092 New Methods and Simulations for Cosmogenic Induced Spallation Removal in Super-Kamiokande-IV D. Mei , A. Hime , Muon-induced background study for underground laboratories . Phys. Rev. D 73 , 053004 ( 2006 ). doi: 10.1103/PhysRevD.73.053004 http://doi.org/10.1103/PhysRevD.73.053004 Muon-induced background study for underground laboratories Y.C. Wu , X.Q. Hao , Q. Yue et al. , Measurement of Cosmic Ray Flux in China JinPing underground Laboratory . Chin. Phys. C 37 , no. 8 , 086001 ( 2013 ). doi: 10.1088/1674-1137/37/1/016001 http://doi.org/10.1088/1674-1137/37/1/016001 Measurement of Cosmic Ray Flux in China JinPing underground Laboratory P. Agnes , L. Agostino , I. F. M. Albuquerque et al. , The veto system of the DarkSide-50 experiment . JINST 11 , no. 03 , P03016 ( 2016 ). doi: 10.1088/1748-0221/11/03/p03016 http://doi.org/10.1088/1748-0221/11/03/p03016 The veto system of the DarkSide-50 experiment D.S. Akerib , R. Rosero , M. Yeh et al. , LUX-ZEPLIN (LZ) Conceptual Design Report . doi: 10.48550/arXiv.1509.02910 http://doi.org/10.48550/arXiv.1509.02910 LUX-ZEPLIN (LZ) Conceptual Design Report Y.L. Yan , W.X. Zhong , S.T. Lin et al. , Study on cosmogenic radioactive production in germanium as a background for future rare event search experiments . Nucl. Sci. Tech. 31 , 55 ( 2020 ). doi: 10.1007/s41365-020-00762-1 http://doi.org/10.1007/s41365-020-00762-1 Study on cosmogenic radioactive production in germanium as a background for future rare event search experiments A. Abusleme et al. , JUNO Physics and Detector . Prog. Part. Nucl. Phys. 123 , 103927 ( 2022 ). doi: 10.1016/j.ppnp.2021.103927 http://doi.org/10.1016/j.ppnp.2021.103927 JUNO Physics and Detector Y.P. Zhang , C.G. Yang . From Daya Bay experiment to Jiangmen Underground Neutrino Observatory experiment . Sci Sin-Phys Mech Astron , 2021 , 51 : 092002 . doi: 10.1360/SSPMA-2021-0088 http://doi.org/10.1360/SSPMA-2021-0088 (in Chinese) From Daya Bay experiment to Jiangmen Underground Neutrino Observatory experiment T. Li , X. Xia , X.T. Huang et al. , Design and Development of JUNO Event Data Model . Chin. Phys. C 41 , no. 6 , 066201 ( 2017 ). doi: 10.1088/1674-1137/41/6/066201 http://doi.org/10.1088/1674-1137/41/6/066201 Design and Development of JUNO Event Data Model K.J. Li , Z.Y. You , Y.M. Zhang et al. , GDML based geometry management system for offline software in JUNO . Nucl. Instrum. Meth. A 908 , 43 - 48 ( 2018 ). doi: 10.1016/j.nima.2018.08.008 http://doi.org/10.1016/j.nima.2018.08.008 GDML based geometry management system for offline software in JUNO T.K. Gaisser , T. Stanev , Muon bundles in underground detectors . Nucl. Instrum. Meth. A 235 , 183 - 192 ( 1985 ). doi: 10.1016/0168-9002(85)90260-8 http://doi.org/10.1016/0168-9002(85)90260-8 Muon bundles in underground detectors E. Scapparone , Study of muon bundles observed by the MACRO underground experiment . Nucl. Phys. B Proc. Suppl. 52 , 179 - 181 ( 1997 ). doi: 10.1016/S0920-5632(96)00875-4 http://doi.org/10.1016/S0920-5632(96)00875-4 Study of muon bundles observed by the MACRO underground experiment M. Grassi , J. Evslin , E. Ciuffoli et al. , Showering Cosmogenic Muons in A Large Liquid Scintillator . JHEP 09 , 049 ( 2014 ). doi: 10.1007/JHEP09(2014)049 http://doi.org/10.1007/JHEP09(2014)049 Showering Cosmogenic Muons in A Large Liquid Scintillator R.P. Kokoulin , N.S. Barbashina , A.G. Bogdanov et al. , Measuring the Cherenkov light yield from cosmic ray muon bundles in the water detector . Nucl. Instrum. Meth. A 952 , 161586 ( 2020 ). doi: 10.1016/j.nima.2018.11.045 http://doi.org/10.1016/j.nima.2018.11.045 Measuring the Cherenkov light yield from cosmic ray muon bundles in the water detector E.A. Yurina , N.S. Barbashina , A.G. Bogdanov et al. , Status of the NEVOD–DECOR Experiment on the Study of Muon Bundles Energy Deposit . Bull. Russ. Acad. Sci. Phys. 85 , no. 4 , 455 - 457 ( 2021 ). doi: 10.3103/S1062873821040390 http://doi.org/10.3103/S1062873821040390 Status of the NEVOD–DECOR Experiment on the Study of Muon Bundles Energy Deposit J.A. Aguilar , I. Al Samarai , A. Albert et al. , A fast algorithm for muon track reconstruction and its application to the ANTARES neutrino telescope . Astropart. Phys. 34 , 652 - 662 ( 2011 ). doi: 10.1016/j.astropartphys.2011.01.003 http://doi.org/10.1016/j.astropartphys.2011.01.003 A fast algorithm for muon track reconstruction and its application to the ANTARES neutrino telescope A. Stoessl , Muon track reconstruction and veto performance with D-Egg sensor for IceCube-Gen2 . PoS ICRC2017 , 1038 ( 2018 ). doi: 10.22323/1.301.1038 http://doi.org/10.22323/1.301.1038 Muon track reconstruction and veto performance with D-Egg sensor for IceCube-Gen2 V.A. Allakhverdyan , A. D. Avrorin , A. V. Avrorin et al. , Measuring muon tracks in Baikal-GVD using a fast reconstruction algorithm . Eur. Phys. J. C 81 , no. 11 , 1025 ( 2021 ). doi: 10.1140/epjc/s10052-021-09825-y http://doi.org/10.1140/epjc/s10052-021-09825-y Measuring muon tracks in Baikal-GVD using a fast reconstruction algorithm G. Bellini , J Benziger , D Bick et al. , Muon and Cosmogenic Neutron Detection in Borexino . JINST 6 , P05005 ( 2011 ). doi: 10.1088/1748-0221/6/05/P05005 http://doi.org/10.1088/1748-0221/6/05/P05005 Muon and Cosmogenic Neutron Detection in Borexino Y. Abe , J.C. dos Anjos , J.C. Barriere et al. , Precision Muon Reconstruction in Double Chooz . Nucl. Instrum. Meth. A 764 , 330 - 339 ( 2014 ). doi: 10.1016/j.nima.2014.07.058 http://doi.org/10.1016/j.nima.2014.07.058 Precision Muon Reconstruction in Double Chooz Zoa Conner , A study of solar neutrinos using the Super-Kamiokande detector , Dissertation ( UNIVERSITY OF MARYLAND COLLEGE PARK , 1997 ). A study of solar neutrinos using the Super-Kamiokande detector V.A. Allakhverdyan , A.D. Avrorin , A.V. Avrorin et al. , Technique for suppression of background cascades produced by atmospheric muon bundles in the Baikal-GVD . JINST 17 , no. 02 , C02013 ( 2022 ). doi: 10.1088/1748-0221/17/02/c02013 http://doi.org/10.1088/1748-0221/17/02/c02013 Technique for suppression of background cascades produced by atmospheric muon bundles in the Baikal-GVD S. Reck T. Eberl , U. Katz , Muon bundle reconstruction with KM3NeT/ORCA using graph convolutional networks . PoS ICRC2021 , 1048 ( 2021 ). doi: 10.22323/1.395.1048 http://doi.org/10.22323/1.395.1048 Muon bundle reconstruction with KM3NeT/ORCA using graph convolutional networks C. Genster , M. Schever , L. Ludhova et al. , Muon reconstruction with a geometrical model in JUNO . JINST 13 , no. 03 , T03003 ( 2018 ). doi: 10.1088/1748-0221/13/03/T03003 http://doi.org/10.1088/1748-0221/13/03/T03003 Muon reconstruction with a geometrical model in JUNO K. Zhang , M. He , W.D. Li et al. , Muon Tracking with the fastest light in the JUNO Central Detector. Radiation Detection Technology and Methods 2.1 ( 2018 ). doi: 10.1007/s41605-018-0040-8 http://doi.org/10.1007/s41605-018-0040-8 Muon Tracking with the fastest light in the JUNO Central Detector. Radiation Detection Technology and Methods 2.1 Z.Y. Li , Y.M. Zhang , G.F. Cao et al. , Event vertex and time reconstruction in large-volume liquid scintillator detectors . Nucl. Sci. Tech. 32 , 49 ( 2021 ). doi: 10.1007/s41365-021-00885-z http://doi.org/10.1007/s41365-021-00885-z Event vertex and time reconstruction in large-volume liquid scintillator detectors Y. Liu , W.D. Li , T. Lin et al. , Muon reconstruction with a convolutional neural network in the JUNO detector . Radiat Detect Technol Methods 5 , 364 - 372 ( 2021 ). doi: 10.1007/s41605-021-00259-4 http://doi.org/10.1007/s41605-021-00259-4 Muon reconstruction with a convolutional neural network in the JUNO detector M. Wurm , Studying neutrino properties in the future LENA experiment . Nucl. Phys. B Proc. Suppl. 237-238 , 314 - 316 ( 2013 ). doi: 10.1016/j.nuclphysbps.2013.04.114 http://doi.org/10.1016/j.nuclphysbps.2013.04.114 Studying neutrino properties in the future LENA experiment S. Zhang , J.S. Li , Y.J. Su et al. , A method for sharing dynamic geometry information in studies on liquid-based detectors . Nucl. Sci. Tech. 32 , 21 ( 2021 ). doi: 10.1007/s41365-021-00852-8 http://doi.org/10.1007/s41365-021-00852-8 A method for sharing dynamic geometry information in studies on liquid-based detectors CERN ROOT . An open-source data analysis framework , https://root.cern.chhttps://root.cern.ch https://root.cern.chhttps://root.cern.ch ; 2022 [accessed 15 March 2022 ]. An open-source data analysis framework CERN ROOT . An open-source data analysis framework , https://root.cern.ch/root/html526/TSpectrum2.html https://root.cern.ch/root/html526/TSpectrum2.html ; 2022 [accessed 15 March 2022 ]. An open-source data analysis framework M. Morháč , J. Kliman , V. Matoušek et al. , Identification of peaks in multidimensional coincidence γ -ray spectra . Nuclear Inst. and Methods in Physics Research, A , 2000 , 443 ( 1 ). doi: 10.1016/S0168-9002(99)01005-0 http://doi.org/10.1016/S0168-9002(99)01005-0 Identification of peaks in multidimensional coincidence γ-ray spectra A. Abusleme , T. Adam , S. Ahmad et al. , Feasibility and physics potential of detecting 8 B solar neutrinos at JUNO . Chin. Phys. C 45 , no. 2 , 023004 ( 2021 ). doi: 10.1088/1674-1137/abd92a http://doi.org/10.1088/1674-1137/abd92a Feasibility and physics potential of detecting 8B solar neutrinos at JUNO All authors contributed to the study conception and design. Material preparation, data collection and analysis were performed by Cheng-Feng Yang, Yong-Bo Huang and Ji-Lei Xu. The first draft of the manuscript was written by Cheng-Feng Yang, Yong-Bo Huang and Ji-Lei Xu, and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript. We thank the JUNO offline and reconstruction working group for many helpful discussions. 10.1007/s41365-022-01049-3 This work was supported by the National Natural Science Foundation of China (Nos. 12005044 and 11975258), the Strategic Priority Research Program of the Chinese Academy of Sciences (Nos. XDA10011200 and XDA10010900) . 2022-01-21 2022-04-19 2022-04-20 2022-06-15 2022-07-26 © China Science Publishing & Media Ltd. (Science Press), Shanghai Institute of Applied Physics, the Chinese Academy of Sciences, ChineseNuclear Society and Springer Nature Singapore Pte Ltd. 2022 This is a post-peer-review, pre-copyedit version of an article published in Nuclear Science and Techniques. The final authenticated version is available online at: http://dx.doi.org/10.1007/s41365-022-01049-3 http://dx.doi.org/10.1007/s41365-022-01049-3 . 2022 Nuclear Science and Techniques 05 33