Measurement of the muon flux from 400 GeV/c protons interacting in a thick molybdenum/tungsten target

The SHiP experiment is proposed to search for very weakly interacting particles beyond the Standard Model which are produced in a 400 GeV/c proton beam dump at the CERN SPS. About 1011\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10^{11}$$\end{document} muons per spill will be produced in the dump. To design the experiment such that the muon-induced background is minimized, a precise knowledge of the muon spectrum is required. To validate the muon flux generated by our Pythia and GEANT4 based Monte Carlo simulation (FairShip), we have measured the muon flux emanating from a SHiP-like target at the SPS. This target, consisting of 13 interaction lengths of slabs of molybdenum and tungsten, followed by a 2.4 m iron hadron absorber was placed in the H4 400 GeV/c proton beam line. To identify muons and to measure the momentum spectrum, a spectrometer instrumented with drift tubes and a muon tagger were used. During a 3-week period a dataset for analysis corresponding to (3.27±0.07)×1011\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(3.27\pm 0.07)~\times ~10^{11}$$\end{document} protons on target was recorded. This amounts to approximatively 1% of a SHiP spill.


Introduction
The aim of the SHiP experiment [1] is to search for very weakly interacting particles beyond the Standard Model which are produced by the interaction of 400 GeV/c protons from the CERN SPS with a beam dump. The SPS will deliver 4 × 10 13 protons on target (POT) per spill, with the aim of accumulating 2 × 10 20 POT during five years of operation. The target is composed of a mixture of TZM (Titanium-Zirconium doped Molybdenum, 3.6λ 1 ), W (9.2λ) and Ta (0.5λ) to increase the charm cross-section relative to the total cross-section and to reduce the probability that long-lived hadrons decay.
An essential task for the experiment is to keep the Standard Model background level to less than 0.1 event after 2 × 10 20 POT. About 10 11 muons per spill will be produced in the dump, mainly from the decay of π, K , ρ, ω and charmed mesons. These muons would give rise to a serious background for many hidden particle searches, and hence their flux has to be reduced as much as possible. To achieve this, SHiP will employ a novel magnetic shielding concept [2] that will suppress the background by five orders of magnitude. The design of this shield relies on the precise knowledge of the kinematics of the produced muons, in particular the muons with a large momentum (>100 GeV/c) and a large transverse momentum (>3 GeV/c) as they can escape the shield and end up in the detector acceptance.
To validate the muon spectrum as predicted by our simulation, and hence the design of the shield, the SHiP Collaboration measured the muon flux in the experiment in the 400 GeV/c proton beam at the H4 beam line of the SPS at CERN in July 2018 [3].

Spectrometer
The experimental setup, as implemented in FairShip (the SHiP software framework), is shown in Fig. 1. A cylindrical SHiP-like 2 target (10 cm diameter and 154.3 cm length) was followed by a hadron absorber made of iron blocks (240 × 240 × 240 cm 3 ) and surrounded by iron and concrete shielding blocks. The dimensions of the hadron absorber were optimised to stop pions and kaons while keeping a good p T acceptance of traversing muons. The SPS beam counters (XSCI.022.480/481, S0 in Fig. 1) and beam counter S1 were used to count the number of POT seen by the experiment.
A spectrometer was placed downstream of the hadron absorber. It consisted of four drift-tube stations (T1-T4, modified from the OPERA experiment [4]) with two stations   [5]. The drift-tubes were arranged in modules of 48 tubes, staggered in four layers of twelve tubes with a total width of approximately 50 cm. The four modules of height 110 cm making up stations T1 and T2 were arranged in a stereo setup (x − u views for T1 and v − x views for T2), with a stereo angle of 60 • . T3 and T4 had only x views and were made of four modules of 160 cm height.
The drift-tube trigger (S2) consisted of two scintillator planes, placed before (S2a) and behind (S2b) the first two tracking stations.
A muon tagger was placed behind the two downstream drift-tube stations. It consisted of five planes of single-gap resistive plate chambers (RPCs), operated in avalanche mode, interleaved with 1 × 80 cm and 3 × 40 cm thick iron slabs. In addition to this, a 80 cm thick iron slab was positioned immediately upstream of the first chamber. The active area of the RPCs was 190 cm × 120 cm and each chamber was read out by two panels of x/y strips with a 1 cm pitch.
The two upstream tracking stations were centered on the beam line, whereas the two downstream stations and the RPCs were centered on the Goliath magnet 3 opening to maximize the acceptance.
The data acquisition was triggered by the coincidence of S1 and S2. For more details on the DAQ framework, see [6], and for a description of the trigger and the DAQ conditions during data taking, see [7].
The protons were delivered in 4.8 s duration spills (slow extraction). There were either one or two spills per SPS supercycle, with intensities ∼ 3 × 10 6 protons per second. The 1sigma width of the beam spot was 2 mm. For physics analysis, 20128 useful spills were recorded with the full magnetic field of 1.5 T, with 2.81×10 11 raw S1 counts. After normalization 3 The centre of the Goliath magnet is 17.86 cm above the beam line.
(see Sect. 3.1) this corresponds to (3.25 ± 0.07) × 10 11 POT. Additional data were taken with the magnetic field switched off for detector alignment and tracking efficiency measurement.

Normalization
The calculation of the number of POT delivered to the experiment must take the different signal widths and dead times of the various scintillators into account. Moreover, some protons from the so-called halo, might fall outside the acceptance of S1 and will only be registered by S0.
In low-intensity runs these effects are small. We select some spills of these runs and split them into 50 slices of 0.1 s. We then determine the number of POT per slice and count the number of reconstructed muons in each slice, which should be independent of the intensity. By leaving the dead times as free parameters in a straight line fit, we find [8] that the number of POT required to have an event with at least one reconstructed muon is 710 ± 15. The systematic error of 15 POT accounts for the variation between the runs used for the normalization. The statistical error is negligible.
The efficiency of the trigger relies on the efficiency of detecting a muon signal in two scintillator planes S2a and S2b (see Fig. 1 and [8]). Each plane is equipped with two photo-multipliers (PMs), and the signal of each of the PMs is recorded for each event. The calculated trigger inefficiency is less than 1‰ and is hence neglected. Multiplying the number of reconstructed muons found in the 20128 spills by 710 we calculated that this data set corresponds to (3.25 ± 0.07) × 10 11 POT.

Tracking
For the drift-tubes, the relation between the measured drifttime and the distance of the track to the wire (the "r -t" relation) is obtained from the Time to Digital Converter (TDC) distribution by assuming a uniformly illuminated tube. When reconstructing the data, the r -t relations are established first by looking the TDC distributions of simple events (i.e. events with at least 2 and a maximum of 6 hits per tracking station). In the simulation, the true drift radius is smeared with the expected resolution. The pattern recognition subsequently selects hits and clusters to form track candidates and provides the starting values for the track fit. The RPC pattern recognition proceeds similarly. Drift-tube tracks are then extrapolated to RPC tracks and tagged as muons if they have hits in at least three RPC stations. Figure 2 shows a two-muon event in the event display.

Momentum resolution
The expected drift-tube hit resolution based on the OPERA results is 270 µm [4]. However, due to residual misalignment and imperfect r -t relations, the measured hit resolution was slightly worse, 373 µm, as shown in Fig. 3. To study the impact of degraded spatial drift-tube resolution the momentum distribution from the simulation was folded with additional smearing as shown in Fig. 4. The tails towards large momentum p and p T are caused mainly by tracks fitted with wrong drift times due to background hits.
From Fig. 4 we conclude that the momentum resolution is not strongly affected by the degraded resolution of the drift-tubes that is observed. The effect of the degraded drifttube resolution is therefore negligible for our studies of the momentum spectrum. To account for residual effects in the track reconstruction, the resolution in the simulation was set to 350 µm.

Tracking efficiencies
The tracking efficiency in the simulation depends on the station occupancy, and in data and simulation the occupancies are different (apparently caused by different amounts of delta rays). By taking this into account, the efficiency in the simulation is reduced from 96.6 to 94.8%.
To determine the tracking efficiency in data, we use the RPCs to identify muon tracks in the data with the magnetic field turned off. We then take the difference between the tracking efficiency in the simulation with magnetic field off (96.9%) and the measured efficiency (93.6%) as the systematic error: 3.3%. For more details on the analysis and reconstruction, see [9].

Comparison with the simulation
A large sample of muons was generated (with Pythia6, Pythia8 [10] and GEANT4 [11] in FairShip) for the background studies of SHiP, corresponding to the number of POT as shown in Table 1. The energy cuts (E min ) of 1 GeV and 10 GeV were imposed to save computing time. The primary proton nucleon interactions are simulated by Pythia8 (using the default tune). The emerging particles are transported by GEANT4 through the target and hadron absorber producing a dataset of also referred to as "mbias" events. A special setting of GEANT4 was used to switch on muon interactions to produce rare dimuon decays of low-  mass resonances. Since GEANT4 does not have production of heavy flavour in particle interactions, an extra procedure was devised to simulate heavy-flavour production not only in the primary pN collision but also in collisions of secondary particles with the target nucleons. For performance reasons, this was done with Pythia6. The mbias and charm/beauty datasets were combined by removing the heavy-flavour contribution from the mbias and inserting the cascade data with appropriate weights. The details of the full heavy-flavour production for both the primary and cascade interactions are described in [12].

Results
The main objective of this study is to validate our simulations for the muon background estimation for the SHiP experiment. For this purpose, we compare the reconstructed momentum distributions ( p and p T ) from data and simulation.
As discussed in the previous section (see also Fig. 4), the events outside the limits ( p > 350 GeV/c or p T > 5 GeV/c) are dominated by wrongly reconstructed trajectories due to background hits and the limited precision of the tracking detector. In SHiP, where the hadron absorber is 5 m long, only muons with momentum p > 5 GeV/c have sufficient energy to traverse the entire absorber. We therefore restrict our comparison to 5 GeV/c < p < 300 GeV/c and p T < 4 GeV/c. For momenta below 10 GeV/c, we only rely on the reconstruction with the tracking detector, since they do not reach the RPC stations. Above 10 GeV/c we require the matching between drift-tube and RPC tracks. Figures 5 and 6 show the p and p T distributions of muon tracks. The distributions are normalized to the number of POT for data (see Sect. 3.1) and simulation respectively. For the simulated sample, muons from some individual sources are also shown in addition to their sum.
In Fig. 7, we show the p T distributions in slices of p. Table 2 shows a numerical comparison of the number of tracks in the different momentum bins.  Figure 8 shows the muon p − p T distribution in data. Figure 9 gives a view of the differences between data and simulation in the p − p T plane. Plotted is the difference between number of data and simulated tracks divided by the sum of the tracks in data and simulation in bins of p and p T .
For momenta above 150 GeV/c, the simulation underestimates tracks with larger p T , while the total number of tracks predicted are in agreement within 20%. The difference between data and simulation is probably caused by a different amount of muons from pion and kaon decays. It was seen that by increasing the contribution of muons from pion and kaon decays in the simulation the difference between data and simulation was reduced.
The FLUKA [13,14] generator is used to determine the radiation levels in the SHiP environment. To benchmark FLUKA with typical settings used for radiological estimates related to muons in the SHiP environment, the muon flux setup was implemented in FLUKA and the simulation with this setup was compared to that made with Pythia/GEANT4. The results of this comparison are given in Annex 1. This independent prediction provides additional support for the validity of the SHiP background simulation.

Conclusions
We have measured the muon flux from 400 GeV/c protons impinging on a heavy tungsten/molybdenum target. The physics processes underlying this are a combination of the production of muons through decays of non-interacting pions and kaons, the production and decays of charm particles and low-mass resonances, and the transportation of the muons through 2.4 m iron. Some 20-30% differences in the absolute rates are observed. The simulation underestimates contributions to larger transverse momentum for higher muon momenta. Given the complexity of the underlying processes, the agreement between the prediction by the simulation and the measured rate is remarkable. Systematic errors for the track reconstruction (3%) and POT normalization (15 POT)/μ-event have been studied and estimated.
A further understanding of the simulation and the data will be obtained with an analysis of di-muon events, the results of which will be the subject of a future publication. We thank M. Al-Turany, F. Uhlig. S. Neubert and A. Gheata their assistance with FairRoot. We acknowledge G. Eulisse and P.A. Munkes for help with Alibuild. The measurements reported in this paper would not have been possible without a significant financial contribution from CERN. In addition, several member institutes made large financial and in-kind contributions to the construction of the target and the spectrometer sub detectors, as well as providing expert manpower for commissioning, data taking and analysis. This help is gratefully acknowledged.

Data Availability Statement
This manuscript has no associated data or the data will not be deposited. [Authors' comment: The raw datasets analysed during the current study are available from the corresponding author on reasonable request.] Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/. Funded by SCOAP 3 .

Appendix A: FLUKA-GEANT4 comparison
Appendix A.1: Simulation samples The geometry of the muon flux spectrometer was reproduced in FLUKA with a few approximations [15]. A large sample of muons was generated with FLUKA for simulating primary proton nucleon interactions as well as the transport of the emerging particles. This sample was used for the comparison with GEANT4. For performance reasons three samples were made with different momentum thresholds (set for all particles). This increased the statistics in the corresponding momentum bins. The number of POT for the three samples is shown in Table 3.
To be consistent with the GEANT4 simulations done for SHiP, the comparison is limited to 5 GeV/c < p < 300 GeV/c and p T < 4 GeV/c . The primary proton-nuclei interactions are simulated and transported through the target and hadron absorber by FLUKA. Special settings of FLUKA were used to include: • full simulation of muon nuclear interactions and production of secondary hadrons; • delta ray production from muons (>10 MeV); • pair production and bremsstrahlung by high-energy muons; • full transport and decay of charmed hadrons and tau leptons; • decays of pions, kaons and muons described with maximum accuracy and polarisation.
The physics settings utilised in the FLUKA simulations were chosen such as to activate all relevant processes like charm decays and most accurate pion and kaon decay descriptions, and to be as close as possible to the physics lists employed in the GEANT4 simulations.

Appendix A.2: Results
In this section, we compare the reconstructed momentum distributions, p and p T , between FLUKA and GEANT4.
Tracks are considered to be muons if they have hits in the T1, T2, T3 and T4 stations. The distributions are taken at the T1 station and normalized to the number of POT.
As shown in Fig. 5, FLUKA predicts a lower rate compared to GEANT4. In the momentum range 5 GeV/c < p < 200 GeV/c, the agreement between the two simulations is at the level of ∼ 20%, above 200 GeV/c there is a discrepancy of a factor ∼ 3.
As shown in Fig. 6, FLUKA predicts a lower rate compared to GEANT4. In the transverse momentum range 0 < p T < 1 GeV/c the agreement between the two simulations is at the level of ∼ 20%, while above 1 GeV/c, there is a discrepancy of a factor ∼ 3 (Figs. 10, 11).
It should be noted that FLUKA does not allow users to change the underlying physics models or cross sections themselves. The uncertainties shown are therefore purely statistical. Given the complexity of the processes underlying the production of muons and the approximations included in the geometry implementations, the agreement between the FLUKA and GEANT4 simulations is reasonable. The differences between FLUKA and GEANT4 over the full muon momentum and transverse momentum spectra are within a factor 3. The large discrepancies of up to a factor 2-3 are mostly in the tails documenting the systematic differences between the FLUKA and GEANT4 models in these regions. Therefore a safety factor of 3 is recommended for future radiological estimates related to muons in the SHiP environment.