Monte-Carlo based performance assessment of ASACUSA's antihydrogen detector

An antihydrogen detector consisting of a thin BGO disk and a surrounding plastic scintillator hodoscope has been developed. We have characterized the two-dimensional positions sensitivity of the thin BGO disk and energy deposition into the BGO was calibrated using cosmic rays by comparing experimental data with Monte-Carlo simulations. The particle tracks were defined by connecting BGO hit positions and hits on the surrounding hodoscope scintillator bars. The event rate was investigated as a function of the angles between the tracks and the energy deposition in the BGO for simulated antiproton events, and for measured and simulated cosmic ray events. Identification of the antihydrogen Monte Carlo events was performed using the energy deposited in the BGO and the particle tracks. The cosmic ray background was limited to 12 mHz with a detection efficiency of 81 %. The signal-to-noise ratio was improved from 0.22 s^{-1/2} obtained with the detector in 2012 to 0.26 s^{-1/2} in this work.


Introduction
Recently, antihydrogen (H) atoms have been produced [1] in a unique cusp trap [2,3,4] developed for the in-flight hyperfine spectroscopy of ground stateH atoms [5,6,7]. The most recent progress is reported in [8,9,10]. In 2012, the ASACUSA Cusp collaboration developed aH detector consisting of a BGO (Bi 4 Ge 3 O 12 ) scintillator disk in combination with a single anode photomultiplier (PMT) and 5 plastic scintillator plates. The detector was able to reject cosmic backgrounds with a high efficiency [11]. In order to further improve the background rejection efficiency, we have developed a newH detector. The single anode PMT has been replaced by 4 multi-anode PMTs (MAPMTs) for 2D photon readout of the BGO. The five plastic scintillators were replaced by a two-layer hodoscope with 32 plastic scintillator bars per layer to determine charged particle tracks with higher resolution [12].  Figure 1 shows a schematic diagram of the structure of the newH detector consisting of the thin BGO disk and the hodoscope. The BGO disk, has a diameter of 90 mm and a thickness of 5 mm and is housed on the vacuum side (10 −7 Pa) of a UHV viewport. The front surface of the BGO disk was coated with a carbon layer of thickness 0.7 µm to reduce multireflections of the light from scintillation on the surface. It was found that the carbon coating improved the position resolution by a factor of ∼ 2 in our previous device [13]. To achieve a position sensitive readout, 4 MAPMTs (Hamamatsu H8500C) each having 8 × 8 anodes with effective area of 49 mm × 49 mm were directly placed on the view port glass as shown in Fig. 2. The output of 8 × 8 anodes were amplified, digitized and stored by an amplifier unit (Clear Pulse 80190) which was a dedicated model for the H8500C and included 8 × 8 charge amplifiers and analogue-to-digital converters with 12 bit resolution.

2.H detector
The hodoscope consists of two layers of 32 plastic scintillator bars arranged in an octagonal configuration [12]. The scintillator bars are 300 mm × 20 mm × 5 mm for the inner layer and 450 mm × 35 mm × 5 mm for the outer. With face to face distances of 200 mm and 300 mm respectively for the inner and outer layers (see Fig. 1(a)). The solid angle covered by the scintillator bars in units of 4π seen from the center of the BGO is ω ∼ 80 %. Silicon photomultipliers (SiPM, KETEK PM3350TS) were connected to both ends of each bar. The output pulses from the SiPMs was amplified by dedicated front end modules described in detail elsewhere (see ref [14]) and recorded by 128 channels waveform digitizers (CAEN V1742). To obtain the relative sensitivity of each channel, 4 MAPMTs were assembled and irradiated by pulsed (200 ns width) light from a LED with the peak wavelength of 470 nm (see also Ref. [13]), the peak emission wavelength of the BGO scintillation light was 480 nm. The voltages applied to the MAPMTs were adjusted such that the total output charge of each MAPMT were equal. Figure 3 (a) shows an example of a 2D map of output charges from one of the 4 MAPMTs averaged over 10 4 pulses of LED light. The channel with the maximum output charge is arbitrarily set to 100 and all other channels are scaled accordingly. The relative gain of the channels of each MAPMT are evaluated using this mapping. This result was compared with the data sheet from the manufacturer where a tungsten filament lamp was used for calibration. By taking the ratio between the measured values and those from the data sheet for each anode, the distribution of the gain ratio as shown in Fig. 3 (b) was obtained. The standard deviation of this distribution is around 5%. Figures 4 (a) and (b) show the 2D charge distribution for two example cosmic rays events after offline gain matching of the individual MAPMT channels described above. In these examples, the BGO surface has been penetrated nearly perpendicular (a) or parallel (b). These figures demonstrate that position sensitive readout has been successfully implemented using 2D MAPMTs.

2D distribution of cosmic rays
In order to reconstruct particle tracks (to be discussed in detail later), we define the position of a hit on the BGO as the center of the anode giving the highest output charge. When the hit positions observed in the 2D distribution are outside of the BGO, such events are removed in the analysis because they are probablyČerenkov light generated in the glasses of the MAPMT or the viewport when high energy charged particles pass through them. The total charge Q is obtained by summing the charges from all channels (see Figs. 4 (a) and (b)).

Energy calibration of the BGO detector
For the energy calibration, cosmic rays were measured and the charge distribution f (Q) was compared with the energy deposition distribution g(E) calculated by a Monte-Carlo simulation using the GEANT4 toolkit 5 [15]. Blue solid circles in Fig. 5 show f (Q) for cosmic rays events when more than 2 inner hodoscope bars are hit in any coincident combination. A bar is considered to be 'hit' when there is a coincidence between the signals of the upstream and downstream SiPMs connected to the bar. 5 Geant4.9.6 Patch-02 was used. The simulation includes, the BGO detector together with the viewport and the vacuum duct. Cosmic rays are generated by the CRY package [16], however,Čerenkov light is not taken included. g(E) is convoluted with the energy resolution which is assumed to be proportional to where α is a fitting parameter. To compare f (Q) and G(E), the relation E = ηQ + E c is assumed, where η and E c are fitting parameters. E c corresponds to the contribution ofČerenkov light generated in the BGO and the glass. G(E) is fit to f (Q) using a least square method in the energy range from 4.5 to 50 MeV with 3 free parameters α, η and E c . The result is shown by the red line in Fig. 5. The goodness of fit is given by χ 2 /ndf = 1.65 for α=0.52 √ MeV, η=2.5 in units proportional to charge over energy and E c =0.50 MeV.
The cosmic ray flux distribution is known to follow cos 2 θ z , where θ z is the zenith angle. In this case, the average path length of cosmic rays in the BGO is evaluated to be approximately 8 mm. Considering that the energy deposited by a minimum ionizing particle (MIP) in the BGO is 9.0 MeV/cm [18], the peak at around 7 MeV in Fig. 5 is attributed to MIPs. As explained below (see subsection 5.3), only events for E > E th = 15 MeV are considered in the following sections. Figure 6 shows a simplified back view of theH detector where the signals from the 4 MAPMTs (shown in the square in the center) and the hodoscope bars (the outermost octagonal arrangement) for a typical cosmic event in the experiment are overlaid.

Identification ofH atoms and cosmic ray suppression
In this example, a spot like pattern is seen on the BGO. Green-colored hodoscope bars are simultaneously hit by the cosmic ray. The red shaded area in the upper left side shows the possible spatial range of the trajectory of the charged particle.
As a best guess, the bisector of the shaded area is taken to define a track which is shown by the thick black line. In the lower right side, 2 neighboring hodoscope bars are hit, in this case, the red shaded area is defined from the edges of the neighboring hodoscope bars and the BGO hit position. The particle track is also defined by a bisector of the shaded area. In this example, the number of tracks is k = 2.
Solid and dashed lines in Fig. 7 show the exprimental and simulated results of cosmic ray count rate n c k as a function of k. The highest n c k is observed for k = 2 and decreases by more than one order of magnitude for each additional unit of k. The track number analysis for simulations of antiprotons (p) irradiating the BGO disk uniformly with the typical count rate of theH atoms in the experiment Ip = 0.1 Hz results in the chain curve (np k ).
In this simulation, the CHIPS model was used in Geant4 for p annihilation. This model was previously tested with respect to the energy deposition analysis forH annihilation in the BGO [1,11]. The multiplicity of annihilation products fromp annihilation was studied using an emulsion detector and agreed with CHIPS results except for annihilation with heavy atoms 6 [19].
The distribution of the chain curve line shows a maximum again, but it decreases weakly as a function of k. When ap annihilates with a nucleus, approximately 3 charged and 2 neutral pions are produced on average [20]. Taking into account the charged pions and the solid angle ω covered by the hodoscope, np 2 is estimated to be 3ω 2 (1 − ω)Ip ∼ 0.04 Hz which can be compared to np 2 ∼ 0.03 Hz in Fig. 7.
As will be discussed in subsection 5.4, n c k can be decreased considerably from the dashed line and is around 5 × 10 −3 Hz as shown by the open circles in Fig. 7 for k = 1, 2 and 3.
Alternatively, np k does not decrease very much as seen from the open triangles. In these cases, n c k is well below np k . Further to this, n c 4 is more than one order of magnitude lower than the open circles and is negligibly small. Therefore we can assume that events for k ≥ 4 can be reasonably attributed top annihilation. In the following subsections, the events for k = 2, k = 3 and then k = 1 are considered.

Events for k = 2
To analyze 2-track events, the track direction is defined by the angle measured anticlockwise from a horizontal line on the x − y plane as shown in Fig. 8. The tracks are numbered in ascending ordered with increasing angle. The corresponding angles of the 1 st and the 2 nd tracks are named θ 1 and θ 2 (θ 1 < θ 2 ), respectively. Further we define θ 12 = θ 2 − θ 1 . Figures 9 (a) and (b) compare the 2D distribution of cosmic events as a function of θ 1 and θ 12 as obtained from experiment and simulation, respectively. It can be seen that the simulation result reproduces the experimental result well. A strong ridge is observed at θ 12 ∼ 180 degrees which corresponds to cosmic rays passing straight through the detector. On the other hand, in θ 1 direction, the distribution spread widely, centered at 90 degrees, which reflects the cosmic ray flux following a cos θ 2 z distribution. Figure 10 (a) shows the result of simulation of the 2D distribution of np 2 as a function of θ 1 and θ 12 . The distribution is very broad, as is expected fromp annihilations at low energy.
The n c 2 background can be decreased by removing events inside an ellipse defined by semi-minor axis ∆θ 2 and semi-major axis of 90 degrees in θ 12 and θ 1 , respectively (see Fig. 9). For example, the n c 2 background is reduced by one order of magnitude for ∆θ 2 = 10 degrees. Using the same cut, only about 10 % of theps are removed because of the different distributions of Fig. 9 and Fig. 10.
It is noted that in both Figs. 9 (a) and (b), we observe an additional small peak at θ 1 ∼ 270 degrees and θ 12 ∼ 40 degrees. Investigating the corresponding events in the simulation, it was found that energetic γ rays in the BGO produce electronpositron pairs and form this peak. The fraction of these events is 1% of the total events in Fig. 9 (a) and (b). In comparison, although the distribution of np 2 in Fig. 10 (a) is very broad, it has a peak at θ 12 ∼ 180 degrees as shown in Fig. 10 (b) which is the projection of (a) onto the θ 12 axis.
The preference at 180 degrees can be explained if we consider a specific type of event, one with three charged pions. When two charged pions hit the hodoscope, with the third pion escaping in a direction close to the beam axis, momentum conservation will favour θ 12 around 180 degrees.

Events for k = 3
The track directions for k = 3 are defined in the same manner in the case of k = 2. Because 3 tracks are involved, there are 3 ways to choose track pairs, which can then be described in an equivalent manner to 2-track events as shown in Figs. 11 (a)-(c). Figures 12 (a) and (b) show the experimental and simulated results of the 2D distributions of n c 3 , obtained by summing 3 distributions of θ 12 vs. θ 1 , θ 13 vs. θ 1 and θ 23 vs. θ 2 , i.e. every event is represented by 3 points corresponding to Figs. 11 (a)-(c) on the plot to conveniently summarize them. The simulation reproduces the experiment very well.
A peak at θ i j ∼ 180 degrees is observed, which is broader than the peak in Fig. 9 (a) and (b). By investigating the corresponding events in the simulation, it was found that a cosmic ray from above generates recoil electrons emitted downward in the BGO which forms the broad peak (see Fig. 11 (a) and (b)). Another peak is seen at θ i ∼ 270 degrees which is formed by the recoil electron together with the incident cosmic ray (see Fig. 11 (c)). Figure 13 shows Projection of (a) onto θ 12 axis. The ellipse is defined like in Fig. 9.
Equivalently to the case of k = 2, the cosmic ray events are expected to be reduced by removing events inside the ellipse (see Fig. 12) defined by semi-minor axis ∆θ 3 and semi-major axis of 90 degrees in θ i j and θ i , respectively. Whilst simultaneously not decreasing np 3 . It is noted that the peak at θ i ∼ 270 degrees in Figs. 12 (a) and (b) is not present when the events inside the ellipse are removed.

Signal-to-noise ratio (SNR) for k ≥ 2
The total rate of cosmic events after the data cut was obtained by taking the sum of n c k for k and is defined by N c i = Σ k≥i n c k . Figures 14 (a) and (b) show the experimental and simulated results of the 2D distributions of N c 2 , respectively, as a function of ∆θ 2 and ∆θ 3 .  The simulations reproduce the experimental data well. To evaluate the difference between the data and the simulation quantitatively, Figs 14 (c) and (d) show experimental and simulated results of N c 2 as a function of ∆θ 2 for ∆θ 3 = 0, 30, 60 and 90 degrees. The difference between the experimental and simulated results is less than 10%. In the later discussion, we analyze the simulation data. Figure 15 (a) shows the simulation result of the detection efficiency ofps 2 defined by i = Np i /Ip as a function of ∆θ 2 and ∆θ 3 , where Np i = Σ k≥i np k . It is shown that 2 decreases as ∆θ 2 and ∆θ 3 increase. We define the signal-to-noise ratio (SNR) by Figure 15 (b) shows the 2D distribution of x 2 as a function of ∆θ 2 and ∆θ 3 . The maximum x 2 is 0.24 s −1/2 at ∆θ 2 = 16 degrees and ∆θ 3 = 0 degrees with N c 2 = 6.9 mHz and 2 = 65 %. As is seen in Fig. 15 (b), to maximize x 2 , ∆θ 3 should be 0 for all values ∆θ 2 . This suggests that all events of np 3 can be identified asps. It is noted that this fact and the optimization of the SNR depend on Ip. Figure 16 shows x 2 as a function of E th , x 2 has a maximum at E th = 15 MeV, but varying only within 1 % in the range of 5 MeV < E th < 20 MeV. Therefore E th is not critical to the Energy threshold (MeV) Figure 16: x 2 as a function of E th , where x 2 is optimized for each E th using the procedure described in Section 5.3.
optimization of x 2 in this range. Figure 17 (a) shows the 2D distribution of n c 1 as a function of the deposition energy in the BGO E and θ 1 . We observe a ridge at E ∼ 7 MeV which corresponds to the MIP peak. In the θ 1 direction, the distribution spreads widely with a center at 90 and 270 degrees and its shape is attributed to the cosmic ray flux distribution of cos θ 2 z . The shape of the ridge appears to be an ellipse. However, the tail of the ridge seems to be more a triangular. Figure 17 (b) shows the 2D distribution of np 1 as a function of E and θ 1 . The distribution is very broad. n c 1 is expected to be reduced by removing events inside the triangle defined by ∆E and the base of 180 degrees as is seen in Fig. 17. Figures 18 (a) and (b) show the 2D distributions of N c 1 and 1 (k ≥ 1), respectively, as a function of ∆θ 2 and ∆E. It is seen that N c 1 and 1 decrease gradually as ∆θ 2 and ∆E increase. Figure 18 (c) shows x 1 as a function of ∆θ 2 and ∆E. The maximum of x 1 reached 0.26 s −1/2 at ∆θ 2 = 14 degrees and ∆E = 93 MeV with N c 1 = 12 mHz and 1 = 81 %. This is larger than the maximum value of x 2 , therefore, the analysis including the events for k = 1 improves the SNR. Comparing with theH detector developed in 2012 with x = 0.22 s −1/2 with the cosmic count rate of 4 mHz and the detection efficiency of 50 % for Ip = 0.1 Hz, the detector described in this work improves upon the SNR and the detection efficiency.

Conclusion
We have developed aH detector consisting of a thin BGO disk and a hodoscope. We have measured hit positions of cos- mic rays in the BGO disk and confirmed that the thin disk with a 2D readout by MAPMTs enables position sensitivity. The energy deposition in the BGO was calibrated by comparing cosmic ray data with Geant4 simulations. Charged particle tracks were determined by connecting the hit position on the BGO and hits on hodoscope bars. By removing the cosmic rays passing through the detector using the cut on ∆E, ∆θ 2 and E th , the background was reduced efficiently to N c 1 = 12 mHz with a detection efficiency of 1 = 81 %. The SNR was improved to x 1 = 0.26 s −1/2 which was compared to 0.22 s −1/2 for the detector used in 2012.