Measurements of total production cross sections for $\pi^{+}$+C, $\pi^{+}$+Al, $K^{+}$+C, and $K^{+}$+Al at 60 GeV/c and $\pi^{+}$+C and $\pi^{+}$+Al at 31 GeV/c

This paper presents several measurements of total production cross sections and total inelastic cross sections for the following reactions: $\pi^{+}$+C, $\pi^{+}$+Al, $K^{+}$+C, $K^{+}$+Al at 60 GeV/c, $\pi^{+}$+C and $\pi^{+}$+Al at 31 GeV/c . The measurements were made using the NA61/SHINE spectrometer at the CERN SPS. Comparisons with previous measurements are given and good agreement is seen. These interaction cross sections measurements are a key ingredient for neutrino flux prediction from the reinteractions of secondary hadrons in current and future accelerator-based long-baseline neutrino experiments.


I. INTRODUCTION
The NA61 or SPS Heavy Ion and Neutrino Experiment (SHINE) [1] at the CERN Super Proton Synchrotron (SPS) has a broad physics program that includes heavy ion physics, cosmic ray physics, and neutrino physics. Long-baseline neutrino beams are typically initiated by high-energy protons that strike a long target, yielding hadrons that can decay to neutrinos or can reinteract in the target or in the aluminum focusing horns, potentially producing additional neutrinoyielding hadrons. NA61/SHINE has already been very successful at measuring the yields of secondary hadrons generated by 31 GeV=c protons on carbon targets [2,3] for the Tokai-to-Kamioka (T2K) long-baseline neutrino oscillation experiment [4]. Data at higher energies are now being collected to benefit other neutrino experiments, particularly MINERνA [5], NOνA [6] that use the current NuMI neutrino beam line at Fermilab, and the proposed DUNE experiment [7] which will use the planned LBNF beam line. The NuMI beam line is initiated by 120 GeV=c protons on a carbon target, while LBNF will use 60-120 GeV=c protons on a carbon or beryllium target.
In addition to the interactions of the primary protons in the neutrino beam targets, a significant fraction of the neutrinos result from hadrons coming from the reinteractions of 10-60 GeV=c energy pions, protons, and kaons in the carbon target, aluminum horns, and other beam line materials. For example, in the on-axis low-energy beam configuration in NuMI, there are on average ∼1.4 hadronic interactions per ν μ [8]. For the current medium-energy NuMI beam configuration and for DUNE, this increases slightly. Uncertainties on the total cross sections for these ancestor hadrons results in uncertainties on the total hadron production rate and production location. Especially for kaon interactions, the existing data is somewhat limited and not very well reproduced by Monte Carlo (MC) [8].
During the fall of 2015, NA61/SHINE recorded interactions of positively charged protons, pions, and kaons on thin carbon and aluminum targets. In the case of pions, interactions were recorded at beam momenta of 31 GeV=c and 60 GeV=c. Kaons were recorded with a beam momentum of 60 GeV=c only, and protons at 31 GeV=c only. The NA61/SHINE vertex magnets were not operational during this period. Therefore, final state particles could not be identified and spectral measurements could not be extracted from this data run. As a result of this setup, data taking was optimized for making measurements of the total production and total inelastic cross sections for each interaction. In the future, NA61 will extract hadron production spectra from data collected more recently with magnetic fields.
The total cross section of hadron-nucleus interactions σ tot can be defined in terms of the inelastic σ inel and coherent elastic σ el cross sections: The inelastic cross section σ inel is defined as the sum of all processes due to strong interactions except coherent nuclear elastic scattering. The production processes are defined as those in which new hadrons are produced. The inelastic processes additionally include interactions which only result in the disintegration of the target nucleus (quasielastic interactions). Taking into account quasielastic scattering as a subset of the inelastic scattering process, one can define the production cross section σ prod in terms of the quasielastic cross section σ qe as This paper is organized as follows: Section II describes the experimental apparatus. Section III presents the event selection to ensure the quality of the measurements. Section IV presents the procedure for measuring σ inel and σ prod cross sections. Section V describes the corrections to the raw trigger probability. Section VI discusses systematic uncertainties. The final results and discussion are presented in Secs. VII and VIII.

II. EXPERIMENTAL SETUP, BEAMS,
AND DATA COLLECTED NA61/SHINE receives a secondary hadron beam from the 400 GeV=c SPS proton beam. The primary proton beam strikes a beryllium target 535 m upstream generating the secondary beam. A magnet system is then used to select the desired beam momentum. Unwanted positrons and electrons are absorbed by two 4 mm lead absorbers.
The NA61/SHINE detector [1] is shown in Fig. 1. In standard operation, it comprises four large time projection chambers (TPCs) and a time-of-flight (ToF) system allowing NA61/SHINE to make spectral measurements of produced hadrons. Two of the TPCs, Vertex TPC 1 (VTPC-1) and Vertex TPC 2 (VTPC-2), are contained within superconducting magnets, capable of generating a combined maximum bending power of 9 T · m. However these magnets were not operational during the 2015 run presented here. Downstream of the VTPCs are the Main TPC Left (MTPC-L) and Main TPC Right (MTPC-R). Additionally, a smaller TPC, the Gap TPC (GTPC), is positioned along the beam axis between the two VTPCs. Two side time-of-flight walls, ToF-left and ToFright, were present. The Projectile Spectator Detector (PSD), a forward hadron calorimeter, sits downstream of the ToF system.
The most critical systems for the analyses of the 2015 data presented here are the trigger system and the beam position detectors (BPDs). The NA61/SHINE trigger system uses two scintillator counters (S1 and S2) to trigger on beam particles. The S1 counter provides the start time for all counters. Three veto scintillation counters (V0, V1 and V1 p ) each with a 1 cm diameter hole are used to remove divergent beam particles upstream of the target. The S4 scintillator with a 1 cm radius sits downstream of the target and is used to determine whether or not an interaction has occurred. A Cherenkov differential counter with achromatic ring focus (CEDAR) [9,10] and a threshold Cherenkov counter (THC) select beam particles of the desired species. The CEDAR focuses the Cherenkov ring from a beam particle onto a ring of 8 PMTs. The pressure is set to a fixed value so that only particles of the desired species will trigger the PMTs, and typically a coincidence of at least 6 PMTs is required to tag a particle for the trigger. Pressure scans taken of the CEDARs are shown in Fig. 2. For these 2015 data at 31 GeV=c the beam was composed of approximately 87% pions, 11% protons, and 2% kaons, and the CEDAR pressure was set to 1.7 bar for the pion beam data, and 3.32 bar for the proton beam data. At 60 GeV=c the beam was composed of approximately 74% pions, 23% protons, and 3% kaons. The CEDAR pressure was set to 1.78 bar for the kaon beam data and 1.68 bar for the pion beam data.
The beam particles are selected by defining the beam trigger (T beam ) as the coincidence of S1 ∧ S2 ∧ V0 ∧ V1 ∧ V1 p ∧ CEDAR ∧ THC. The interaction trigger (T int ) is defined by the coincidence of T beam ∧ S4 to select beam particles which have interacted with the target. A correction factor will be discussed in detail in Sec. VA to correct for interactions that hit the S4. Three BPDs, which are proportional wire chambers, are located 30.39, 9.09, and 0.89 m upstream of the target and determine the location of the incident beam particle to an accuracy of ∼100 μm.
For these 2015 data, the interactions of p, π þ , and K þ beams were measured on thin carbon and aluminum targets. The carbon target was composed of graphite of density ρ ¼ 1.84 g=cm 3 with dimensions of 25 mmðWÞ × 25 mmðHÞ× 20 mmðLÞ, corresponding to roughly 4% of a protonnuclear interaction length. The aluminum target has a density of ρ ¼ 2.70 g=cm 3 with dimensions of 25 mm ðWÞ× 25 mm ðHÞ × 14.8 mm ðLÞ, corresponding to roughly 3.6% of a proton-nuclear interaction length.

A. Event selection
Several cuts were applied to events to ensure the purity of the measurement and to control the systematic effects caused by beam divergence. First, the so-called WFA (wave form analyzer) cut was used to remove events in which multiple beam particles pass through the beam line in a small time frame. The WFA determines the timing of beam particles that pass through the S1 scintillator, with a resolution of 100 nsec. If another beam particle passes through the beam line close in time to the triggered beam particle, it could cause a false trigger in the S4 scintillator. In order to mitigate this effect, a conservative cut of AE2 μs was applied to the time window to ensure that only one particle is allowed to pass through the S1 in a 4 μs time window around the selected beam particle.
The trajectories of the incoming beam particles are measured by three BPDs, located along the beam line upstream of the target as shown in Fig. 1. The measurements from the BPDs are especially important for estimating the effects of beam divergence on the cross section measurements. To understand these effects, tracks are fitted to the reconstructed BPD clusters, and the tracks are extrapolated to the S4 plane. The so-called "good BPD" cut requires that the event includes a cluster in the mostdownstream BPD and that a track was successfully fit to the BPDs. Figures 3 and 4 show the resulting BPD  extrapolation to the S4 plane for the interactions studied. The left plots show the extrapolated positions for all beam particles that pass the beam trigger, and the right plots the extrapolated positions for beam particles that pass the interaction trigger, which requires an anticoincidence with the S4 scintillator. It can be seen from these figures that the 31 GeV=c beams were much wider than the 60 GeV=c beams. From these figures, it is also evident that the V1 veto counter (which is close to the most-downstream BPD) and S4 were not well aligned, especially for the 31 GeV=c beam. The beam was wide enough that a significant fraction of the beam particles have trajectories missing the S4, as can be see in the right plot of Fig. 3 where a halo of beam particles misses the edge of the S4, mimicking the interaction trigger. This leads to an apparent interaction rate higher than the actual interaction rate. To reduce this effect, a radial cut was applied to the BPD tracks extrapolated to the S4, and this is indicated by the red circles on Figs. 3 and 4.
The number of events after the described selection cuts for the interactions: 60 GeV=c K þ and π þ and 31 GeV=c π þ with C and Al targets (target inserted) and with the targets removed (target removed) are shown in Tables I-III.

IV. INTERACTION TRIGGER CROSS SECTIONS
In general, the probability of a beam particle interaction inside of a thin target is proportional to the thickness L of the target and the number density of the target nuclei n. Thus, the interaction probability P can be defined by taking into account the thin target approximation and by defining the interaction cross section σ as   Target  Inserted  Removed  Inserted  Removed   Total  593k  195k  535k  234k  WFA  591k  195k  532k  233k  Good BPD  547k  180k  491k  215k  Radial cut  437k  143k  367k  159k   TABLE III. Event selection table for K þ þ C and K þ þ Al at 60 GeV=c. Target  Inserted  Removed  Inserted  Removed   Total  505k  239k  339k  156k  WFA  503k  238k  337k  155k  Good BPD  466k  221k  312k  144k  Radial cut  463k  219k  310k  143k   TABLE II. Event selection table for π þ þ C and π þ þ Al at 60 GeV=c.
Interaction π þ þ C π þ þ Al Target  Inserted  Removed  Inserted  Removed   Total  528k  247k  459k  286k  WFA  513k  240k  448k  279k  Good BPD  479k  225k  417k  260k  Radial cut  463k  217k  405k  252k The density of nuclei n can be calculated in terms of N A , ρ, and A, which are Avogadro's number, the material density, and the atomic number, respectively. The counts of beam and interaction triggers as described in Sec. II can be used to estimate the trigger probability as follows: where NðT beam Þ is the number of beam events passing the event selection cuts and NðT beam ∧ T int Þ is the number of selected beam events which also have an interaction trigger. In order to correct for events in which the beam particle interacts outside of the target, data were also recorded with the target removed from the beam (target removed) by rotating the target holder out of the path of the beam. Figure 5 shows an example of the trigger interaction probabilities for each run for the π þ þ C at 60 GeV=c data set. Table IV gives the total trigger interaction probabilities for the data sets used in this paper for both the target inserted and target removed data. The kaon target removed interaction probabilities are larger than those for pions due to the fact that ∼1% of the beam kaons will decay between BPD3 and S4. Taking into account the trigger probabilities with the target inserted (I) and the target removed (R), P I Tint and P R Tint , the interaction probability P int can be obtained: Equation (3) leads to the definition of the trigger cross section σ trig , calculated with P int and the effective target length L eff , which accounts for the exponential beam attenuation: The effective target length can be calculated with the absorption length, where By combining Eqs. (6)-(8), one can simplify the equation for σ trig as V. CORRECTION FACTORS

A. S4 trigger correction factors
The trigger cross section contains the interactions where the resulting particles miss the S4 scintillator counter that is   downstream of the target. But even when there has been an interaction in the target, there is a possibility that a forwardgoing particle will strike the S4 counter. Moreover, not all elastically scattered beam particles strike the S4. Corrections must be applied to account for these effects. Combining Eqs. (1) and (2), the trigger cross section can be related to the production cross section through MC correction factors as follows: where f prod , f qe , and f el are the fractions of production, quasielastic, and elastic events that miss the S4 counter. σ qe and σ el are also estimated from Monte Carlo. Equation (10) can be rewritten to obtain σ prod and σ inel as and A GEANT4 detector simulation [11][12][13] was used to estimate the MC correction factors discussed above. The FTFP_BERT physics list with GEANT4 version of 10.2.p03 was used to estimate correction factors as presented in Table V.

B. Beam composition correction factors
In the case of π þ beams, a correction must also be applied to account for contamination from μ þ and e þ . The CEDAR and threshold Cherenkov detectors do not have the power to completely discriminate positrons and muons from pions at 31 GeV=c and 60 GeV=c as shown in [9,10]. Fortunately, it was possible to estimate the amount of positron contamination with the TPC system and the PSD. During the neutrino data taking in 2016, a special maximum field data run was taken during which the 60 GeV=c h þ beam was bent into the MTPC-L.
Positrons deposit most of their energy in the first 2 out of 10 longitudinal sections of the PSD, while pions penetrate deeper. A pure pion sample is obtained by selecting beam particles that deposit less than 20% of their total energy in the first section of the PSD. To determine the e þ and π þ compositions of the beam, a sum of two Gaussians is fit to the dE=dx distribution. From the fit, the positron contamination was determined to be 2% AE 2% for the 60 GeV=c beam. Figure 6 shows the resulting fit to the maximum field data. The GEANT4 MC simulation is used to determine the effect of the positrons on the trigger cross section. The resulting corrections applied to σ prod (σ inel ) are þ2.2% (þ2.1%) for π þ þ C at 60 GeV=c and þ1.8% (1.7%) for π þ þ Al at 60 GeV=c.
In the case of 31 GeV=c, the potential for positron contamination was reduced by requiring that the CEDAR had a more stringent sevenfold coincidence signal. No special data run was undertaken with the 31 GeV=c beam to measure the positron contamination, so no correction is applied. But this contamination will be taken into account later as an asymmetric systematic uncertainty.
For the pion beams at both 31 GeV=c and 60 GeV=c, a small number of muons are also present in the beam due to the decays of pions upstream of the target. Many of these muons diverge from the beam and will strike the veto counters, but beam simulations at both momenta show the muon fraction that will pass the veto counters and trigger our beam counters is about 1.5 AE 0.5% of the pion beam. A correction for the muon component of the beam is applied to the 31 GeV=c and 60 GeV=c pion beam interactions.
For the kaon beam, any kaons that decay upstream of the CEDAR will not satisfy the beam selection and will not be selected as good beam particles. Only kaon decays downstream of the CEDAR where the decay products head  6. The binned data shows the dE=dx distribution of the maximum field data set for the 60 GeV=c π þ beam. Overlaid is the sum of Gaussians fit to the histogram as well as the individual π þ and e þ components. From this fit, the positron contribution was estimated to be 2%.
toward the S4 will pass the beam selection and "good BPD" cut. It was estimated that only 0.1% of the CEDAR-tagged kaons will decay with decay products that pass these cuts. Therefore, no correction is applied for kaon decays in the beam line.

VI. SYSTEMATIC UNCERTAINTIES
A. Target density uncertainty The uncertainty on the target density affects the calculation of the trigger cross section as shown in Eq. (9). The density uncertainty for each target is estimated by calculating the standard deviation of the target densities determined from measurements of the mass and dimensions of the machined target samples. A 0.65% uncertainty on the density of carbon and a 0.29% uncertainty on the density of aluminum were used.

B. Out-of-target interactions
As shown in Eq. (5), the measured interaction rates are corrected for interactions occurring outside of the target by measuring the trigger rates with the target both inserted and removed. To look for possible additional systematic effects, two special runs were undertaken with the target holder in the "I" position and with the target holder in the "R" position, but with no target attached. The data were taken with 31 GeV=c and 60 GeV=c π þ .
In the case of the 31 GeV=c target holder data, there was no significant difference between the trigger probability of the empty target holder data and the target removed data. However, in the case of the 60 GeV=c data, a high trigger probability in the target holder "I" run was observed. These additional out-of-target interactions may be related to the beam conditions during those runs. An asymmetric uncertainty was assigned for the 60 GeV=c interactions.

C. S4 size uncertainty and efficiency
Another systematic uncertainty comes from the uncertainty in the size of the S4 scintillator. The diameter of the S4 has previously been found to have an uncertainty of AE0.40 mm. In order to propagate this uncertainty to σ inel and σ prod , two additional MC simulation samples with the S4 diameter modified were generated.
Previous NA61/SHINE analyses have found that S4 inefficiency is negligibly small [14] and this analysis also found no S4 inefficiency by looking at GTPC tracks in target removed data. The S4 inefficiency is concluded to be less than 0.1% and neither an uncertainty nor a correction relating to the S4 scintillator efficiency is applied to the results.

D. Beam composition uncertainty
As was mentioned in Sec. V B, for interactions with the 60 GeV=c π þ beam, a correction was applied to reflect the small amount of positrons in the beam. To be conservative,  100% of this correction is assumed as a systematic uncertainty. For π þ at 31 GeV=c, no correction is applied, but an uncertainty is reported accounting for a 1% positron contamination.
As was also mentioned in Sec. V B, the muon fraction in the pion beam is estimated to be 1.5% for both the 31 GeV=c and 60 GeV=c π þ beams and a correction was applied. An uncertainty of 0.5% is applied to this correction.
The CEDAR counter has a high purity of identifying kaons using a sixfold coincidence for 60 GeV=c beams. The lower limit on the purity of the kaon beam was calculated to be 99.4% according to the CEDAR gas pressure scan data. The estimated systematic error from this source is applied to the total systematic uncertainty.

E. Model uncertainties
The S4 correction factors f prod , f inel , f el and f qe as well as the cross sections σ qe and σ el were estimated with GEANT4 MC simulations using the FTFP_BERT physics list. In order to estimate the model uncertainties associated with these correction factors, the correction factors were recalculated with three additional physics lists: QBBC, QGSP_BERT and FTF_BIC. Using these additional physics lists, the model dependency on the total cross section measurements was studied.
All of systematic uncertainties on the production and inelastic cross section measurements are summarized in Tables VI and VII respectively.

VII. RESULTS
Several production cross sections have been measured in this analysis: π þ þ C (π þ þ Al) at 31 GeV=c is found to be 158.3 mb (310.4 mb), π þ þ C (π þ þ Al) at 60 GeV=c is found to be 171.6 mb (321.0 mb), and K þ þ C (K þ þ Al) at 60 GeV=c is found to be 144.5 mb (284.0 mb), respectively. Statistical, systematic, and physics model uncertainties are estimated separately and are summarized in Table VIII. π þ and K þ at 60 GeV=c measurements are compared with the results of Carrol et al. [15] as shown in Fig. 7. These NA61 results are consistent within our errors with the previous measurements, and our error bands are smaller, especially for the kaons.
Several inelastic cross sections have also been determined in this analysis: π þ þ C (π þ þ Al) at 31 GeV=c is found to be 177.0 mb (340.0 mb), π þ þ C (π þ þ Al) at 60 GeV=c is found to be 188.2 mb (347.0 mb), and K þ þ C (K þ þ Al) at 60 GeV=c is found to be 159.0 mb (307.5 mb), respectively. Statistical, systematic, and physics model uncertainties are estimated separately and are summarized in Table IX. These measurements are compared with the results of Denisov et al. [16] as shown in Fig. 8. These NA61 results are consistent within errors with the existing measurements at 30 GeV=c.
Additionally, a short data run of interactions of 31 GeV=c protons with carbon was analyzed as a crosscheck with the previous higher statistics NA61/SHINE total cross section results from the 2009 T2K data run [2].   The total production (total inelastic) cross section was found to be 229.8

VIII. SUMMARY
In summary, the production and inelastic cross sections of π þ and K þ on carbon and aluminum targets have been measured with the NA61/SHINE experiment. The production cross section with π þ beams at 31 GeV=c was measured for the first time with a precision of about 2%. At 60 GeV=c the measured production cross sections are comparable to previous results for π þ and K þ and the precision was improved to about 3% and 2%, respectively. Inelastic cross section measurements with π þ and K þ beams at 60 GeV=c were measured for first time with precisions of about 3% and 2%, respectively. For the inelastic production cross section for π þ at 31 GeV=c reasonable agreement with a previous measurement was found. Especially for π þ beams, the measurements here are limited by positron contamination in the beam and steps will be taken in future data taking to better limit this uncertainty.
The current uncertainties on the neutrino fluxes in the NuMI neutrino beam at Fermilab from the MINERνA Collaboration [8] rely on measurements of the inelastic cross section (which is termed the "absorption" cross section in the MINERνA paper). For π þ þ C and π þ þ Al they assumed an uncertainty of 5%, while for the K þ þ C and K þ þ Al cross sections they assumed a 10%-30% uncertainty, which is significantly larger than the systematic uncertainties determined in this paper. Thus, these data will greatly reduce the uncertainty on the neutrino flux prediction in NuMI due to kaon interactions.