Large-angle production of charged pions by 3 GeV/c - 12 GeV/c protons on carbon, copper and tin targets

A measurement of the double-differential $\pi^{\pm}$ production cross-section in proton--carbon, proton--copper and proton--tin collisions in the range of pion momentum $100 \MeVc \leq p<800 \MeVc$ and angle $0.35 \rad \le \theta<2.15 \rad$ is presented. The data were taken with the HARP detector in the T9 beam line of the CERN PS. The pions were produced by proton beams in a momentum range from 3 \GeVc to 12 \GeVc hitting a target with a thickness of 5% of a nuclear interaction length. The tracking and identification of the produced particles was done using a small-radius cylindrical time projection chamber (TPC) placed in a solenoidal magnet. An elaborate system of detectors in the beam line ensured the identification of the incident particles. Results are shown for the double-differential cross-sections at four incident proton beam momenta (3 \GeVc, 5 \GeVc, 8 \GeVc and 12 \GeVc).


Introduction
The HARP experiment [1] makes use of a large-acceptance spectrometer for a systematic study of hadron production on a large range of target nuclei for beam momenta from 1.5 GeV/c to 15 GeV/c. The main motivations are the measurement of pion yields for a quantitative design of the proton driver of a future neutrino factory [2], a substantial improvement of the calculation of the atmospheric neutrino flux [3,4,5,6,7] and the measurement of particle yields as input for the flux calculation of accelerator neutrino experiments, such as K2K [8,9], MiniBooNE [10] and SciBooNE [11].
The measurement of the double-differential cross-section, d 2 σ π /dpdθ for π ± production by protons of 3 GeV/c, 5 GeV/c, 8 GeV/c and 12 GeV/c momentum impinging on a thin carbon, copper or tin target of 5% nuclear interaction length are presented.
Especially for carbon targets it is interesting to measure pion production cross-sections in the framework of the HARP measurement programme for neutrino flux calculations. Carbon targets are frequently used as hadron production targets in neutrino beam lines. In addition, measurements on carbon can be used to predict pion production off nitrogen and oxygen nuclei without a large extrapolation in the production models. The knowledge of the latter production cross-sections are needed to model atmospheric muon and neutrino fluxes. Owing to the relatively low incoming beam momenta the data are especially interesting for the calculation of hadron production in secondary interactions in extended production targets and in atmospheric flux calculations. The comparison of the measurements on copper and tin targets with the carbon data in this paper and with the tantalum data obtained with the same apparatus described in Ref. [12] can be used to check the dependence on the atomic number A in hadron production models. Copper and tin are interesting target materials as their atomic numbers are midway between light target materials, such as Be, Al and C (used in targets for conventional neutrino beams) and heavy targets such as Ta (relevant for the optimization of neutrino factory designes).
Data were taken in the T9 beam of the CERN PS. For this analysis, about 1,159,000 (1,066,000 and 1,284,000) incoming protons were selected which gave an interaction trigger in the Large Angle spectrometer, resulting in 235,000 (209,500 and 243,400) well-reconstructed secondary pion tracks for the carbon (copper and tin) target. The different settings have been taken within a short running period so that in their comparison detector variations are minimized.
The analysis proceeds by selecting tracks in the Time Projection Chamber (TPC) in events with incident beam protons. Momentum and polar angle measurements and particle identification are based on the measurements of track position and energy deposition in the TPC. An unfolding method is used to correct for experimental resolution, efficiency and acceptance and to obtain the double-differential pion production cross-sections, with a full error evaluation. A comparison with available data is presented. The analysis follows the same methods as the ones used for the determination of π ± production crosssections by protons on a tantalum target described in Ref. [12]. We refer to Ref. [12] for a detailed description of the analysis, only the main points and differences with respect to the latter are described here.

Experimental apparatus
The HARP detector is shown in Fig. 1. The forward spectrometer is built around a dipole magnet for momentum analysis, with large planar drift chambers (NDC) [13] for particle tracking and a time-of-flight wall (TOFW) [14], a threshold Cherenkov detector (CHE) and an electromagnetic calorimeter (ECAL) used for particle identification. The forward spectrometer covers an acceptance for tracks originating from the target with polar angles up to 250 mrad. This is well matched to the angular range of interest for the measurement of hadron production to calculate the properties of conventional accelerator neutrino beams [16,17]. In the large-angle region a cylindrical TPC with a radius of 408 mm is positioned in a solenoidal magnet with a field of 0.7 T. The target is inserted into the inner field cage of the TPC. The TPC is used for tracking, momentum determination and the measurement of the energy deposition dE/dx for particle identification [18]. A set of resistive plate chambers (RPC) form a barrel inside the solenoid around the TPC to measure the arrival time of the secondary particles [19]. Beam instrumentation provides identification of the incoming particle, the determination of the time when it hits the target, and the impact point and direction of the beam particle on the target. Several trigger detectors are installed to select events with an interaction and to define the normalization.
In addition to the data taken with the thin carbon, copper and tin targets of 5% nuclear interaction length (λ I ), runs were also taken with an empty target holder, a thin 2% λ I target and a thick 100% λ I target. Data taken with a liquid hydrogen target at 3 GeV/c, 5 GeV/c and 8 GeV/c incident beam momentum together with cosmic-ray data were used to provide an absolute calibration of the efficiency, momentum scale and resolution of the detector. Moreover, tracks produced in runs with Pb and Ta targets in the same period and with the same beam settings were used for the calibration of the detector, verification of the event reconstruction and analysis procedures (see Ref. [12] for further details).
The momentum of the T9 beam is known with a precision of the order of 1% [20]. The absolute normaliza-tion of the number of incident protons was performed using 300 000, 240 000 and 280 000 'incident-proton' triggers for the carbon, copper and tin data, respectively. These are triggers where the same selection on the beam particle was applied but no selection on the interaction was performed. The rate of this trigger was down-scaled by a factor 64. A cross-check of the absolute normalization was provided by counting tracks in the forward spectrometer.
A detailed description of the HARP apparatus is given in Ref. [15]. In this analysis the detector components of the large-angle spectrometer and the beam instrumentation are employed and are briefly summarized in the following.
A set of four multi-wire proportional chambers (MWPCs) measures the position and direction of the incoming beam particles with an accuracy of ≈1 mm in position and ≈0.2 mrad in angle per projection. At low momenta the precision of the prediction at the target is limited by multiple scattering. A beam time-of-flight system (BTOF) measures time difference of particles over a 21.4 m path-length. It is made of two identical scintillation hodoscopes, TOFA and TOFB (originally built for the NA52 experiment [21]), which, together with a small target-defining scintillator (TDS), also used for the trigger, provide particle identification at low energies. This provides separation of pions, kaons and protons up to 5 GeV/c and determines the initial time at the interaction vertex (t 0 ). The timing resolution of the combined BTOF system is about 70 ps. A system of two N 2 -filled Cherenkov detectors (BCA and BCB) is used to tag electrons at low energies and pions at higher energies. The electron and pion tagging efficiency is found to be close to 100%. The proton fraction in the incoming beam varies from 35% at 3 GeV/c to 92% at 12 GeV/c. The length of the accelerator spill is 400 ms with a typical intensity of 15 000 beam particles per spill. The average number of events recorded by the data acquisition ranges from 300 to 350 per spill for the four different beam momenta.
The target is placed inside the inner field cage (IFC) of the TPC such that, in addition to particles produced in the forward direction, backward-going tracks can be measured. All three targets have a nominal thickness of 5% λ I and a cylindrical shape with a nominal diameter of 30 mm. The 99.95% pure carbon target used for the measurement described here has a thickness of 18.94 mm with a variation of ±0.02 mm. Its density was measured to be 1.88 g/cm 3 . The copper target has a purity of 99.99% with a thickness of 7.52 mm with a variation of ±0.01 mm and a density of 8.92 g/cm 3 . The tin target has a purity of 99.99% with a thickness of 11.04 mm with a variation of ±0.04 mm and a density of 7.29 g/cm 3 . A set of trigger detectors completes the beam instrumentation: a thin scintillator slab covering the full aperture of the last quadrupole magnet in the beam line to start the trigger logic decision (BS); a small scintillator disk, TDS, positioned upstream of the target to ensure that only particles hitting the target cause a trigger; and 'halo' counters (scintillators with a hole to let the beam particles pass) to veto particles too far away from the beam axis. A cylindrical detector (inner trigger cylinder, ITC) made of six layers of 1 mm thick scintillating fibres is positioned inside the inner field cage of the TPC and surrounds the target. It provides full coverage of the acceptance of the TPC. The large-angle spectrometer consists of a TPC and a set of RPC detectors inside the solenoidal magnet. The TPC detector was designed to measure and identify tracks in the angular region from 0.25 rad to 2.5 rad from the beam axis. Charged particle identification (PID) can be achieved by measuring the ionization per unit length in the gas (dE/dx) as a function of the total momentum of the particle. Additional PID can be performed through a time-of-flight measurement with the RPCs.
In the present analysis, the TPC provides the measurement for the pattern recognition to find the particle tracks, and to measure their momentum through the curvature of their trajectory. It also provides PID using the measurement of energy deposition. The RPC system is used in this analysis to provide a calibration of the PID capabilities of the TPC.
In addition to the usual need for calibration of the detector, a number of hardware shortfalls, discovered mainly after the end of data-taking, had to be overcome to use the TPC data reliably in the analysis. The TPC is affected by a relatively large number of dead or noisy pads and static and dynamic distortions of the reconstructed trajectories. Static distortions are caused by the inhomogeneity of the electric field, due to an accidental mismatch between the inner and outer field cage (powered by two distinct HV supplies) and other sources. Dynamic distortions are caused instead by the build-up of ion-charge density in the drift volume during the 400 ms long beam spill. All these effects were fully studied and available corrections are described in detail in Ref. [12]. While methods to correct the dynamic distortions of the TPC tracks are being implemented, a pragmatic approach has been followed in the present analysis. Only the events corresponding to the early part of the spill, where the effects of the dynamic distortions are still small, are used 1 . The time interval between spills is large enough to drain all charges in the TPC related to the effect of the beam. The combined effect of the distortions on the kinematic quantities used in the analysis has been studied in detail and only that part of the data for which the systematic errors can be assessed with physical benchmarks was used, as explained in [12]. More than 40% of the recorded data can be used on average in the current analysis.
The absolute scale of the momentum determination is determined using elastic scattering data off a hydrogen target. The angle of the forward scattered particle (pion or proton) is used to give an absolute prediction for the momentum of the recoil proton. This prediction is compared with the measurement in the TPC. To study the stability of this measurement protons are selected in a narrow band with a relatively large dE/dx where the dE/dx depends strongly on momentum. The average momentum for the protons selected in this band remains stable within 3% as a function of time-in-spill over the part of the spill used for the analysis.

Data selection and analysis
The beam of positive particles used for this measurement contains mainly positrons, pions and protons, with small components of kaons, deuterons and heavier ions. Its composition depends on the selected beam momentum. The analysis proceeds by first selecting a beam proton hitting the target, not accompanied by other tracks. Then an event is required to be triggered by the ITC in order to be retained. After the event selection the sample of tracks to be used for analysis is defined. Tracks are only considered if they contain at least twelve space points out of a maximum of twenty. This cut is applied to ensure a good measurement of the track parameters and of the dE/dx. Furthermore, a quality requirement is applied on the fit to the helix. The latter requirement introduces a very small loss of efficiency. For tracks satisfying these conditions, a cut is made on d ′ 0 , the distance of closest approach to the extrapolated trajectory of the incoming beam particle in the plane perpendicular to the beam direction and z ′ 0 , the z-coordinate where the distance of the secondary track and the beam track is minimal. Finally, only tracks with momentum in the range between 100 MeV/c and 800 MeV/c are accepted. In addition, particles with transverse momentum below 55 MeV/c are removed. Table 1 shows the number of events and tracks at various stages of the selection. The total number of events taken by the data acquisition ("Total DAQ events") includes triggers of all types as well as calibration events; the number of "Protons on target" represents the count of the incoming beam trigger after off-line selection of accepted protons multiplied by the down-scale factor 64. The number of accepted events for this analysis ("Accepted protons with LAI (Large Angle Interaction)") is obtained using the same selection of incoming protons in coincidence with a trigger in the ITC. The large difference between the rows "Total DAQ events" and "Accepted protons with LAI" is due to the relatively large fraction of pions in the beam and to the larger number of triggers taken for the measurements with the forward spectrometer. These data will be the subject of other publications. The line "Maximum N evt " refers to the last number of events N evt in spill used to avoid dynamic distortion corrections, with the corresponding number of interaction triggers used in the analysis ("LAI in accepted spill part") and the fraction of the data used given under "Fraction of triggers used". The lines "Accepted momentum determination" and "In kinematic region and originating from target" give the number of tracks passing the momentum fit quality requirements and the selection of tracks originating in the target region. Finally, the rows "Negative particles", "Positive particles ", "π − selected with PID" and "π + selected with PID" show the number of accepted tracks with negative and positive charge and the ones passing in addition the pion PID criteria, respectively.
To give an impression of the complexity of the events, one can define an 'average multiplicity' as the ratio of the number of tracks with at least twelve hits in the TPC (regardless of their momentum, angle or spatial position) and the number of events accepted by the selection criteria with at least one such track. With this definition, the average multiplicity is 2.2, 2.6, 3.1 and 3.4 for the 3 GeV/c, 5 GeV/c, 8 GeV/c and 12 GeV/c beams in p-C data, respectively.
The double-differential cross-section for the production of a particle of type α can be expressed in the laboratory system as: where d 2 σα dpidθj is expressed in bins of true momentum (p i ), angle (θ j ) and particle type (α).
The factor A NAρt is the inverse of the number of target nuclei per unit area (A is the atomic mass, N A is the Avogadro number, ρ and t are the target density and thickness) 2 . The result is normalized to the number of incident protons on target N pot .
The 'raw yield' N α ′ i ′ j ′ is the number of particles of observed type α ′ in bins of reconstructed momentum (p i ′ ) and angle (θ j ′ ). These particles must satisfy the event, track and PID selection criteria. Although, owing to the stringent PID selection, the background from misidentified protons in the pion sample is small, the pion and proton raw yields (N α ′ i ′ j ′ , for α ′ = π − , π + , p) have been measured simultaneously. This makes it possible to correct for the small remaining proton background in the pion data without prior assumptions concerning the proton production cross-section.
The matrix M −1 ijαi ′ j ′ α ′ corrects for the efficiency and resolution of the detector. It unfolds the true variables ijα from the reconstructed variables i ′ j ′ α ′ with a Bayesian technique [22] and corrects the observed number of particles to take into account effects such as trigger efficiency, reconstruction efficiency, acceptance, absorption, pion decay, tertiary production, PID efficiency, PID misidentification and electron background. The method used to correct for the various effects is described in more detail in Ref. [12].
In order to predict the population of the migration matrix element M ijαi ′ j ′ α ′ , the resolution, efficiency and acceptance of the detector are obtained from the Monte Carlo. This is accurate provided the Monte Carlo simulation describes these quantities correctly. Where some deviations from the control samples measured from the data are found, the data are used to introduce (small) ad hoc corrections to the Monte Carlo.
Using the unfolding approach, possible known biases in the measurements are taken into account automatically as long as they are described by the Monte Carlo. For example the energy-loss of particles inside the target and material around the inner field cage is expressed as an average shift of the measured momentum distribution compared to the physical momentum. Known biases are therefore treated in the same way as resolution effects. In the experiment simulation, which is based on the GEANT4 toolkit [23], the materials in the beam-line and the detector are accurately described as well as the relevant features of the detector response and the digitization process. The Monte Carlo simulation compares well with data, as shown in Ref. [12].
The absolute normalization of the result is calculated in first instance relative to the number of incident beam particles accepted by the selection. After unfolding, the factor A NAρt is applied. The beam normalization using down-scaled incident-proton triggers has uncertainties smaller than 2% for all beam momentum settings.
The background due to interactions of the primary protons outside the target (called 'Empty target background') is measured using data taken without the target mounted in the target holder. Owing to the selection criteria which only accept events from the target region and the good definition of the interaction point this background is negligible 3 (< 10 −5 ).
The effects of these uncertainties on the final results are estimated by repeating the analysis with the relevant input modified within the estimated uncertainty intervals. In many cases this procedure requires the construction of a set of different migration matrices. The correlations of the variations between the cross-section bins are evaluated and expressed in the covariance matrix. Each systematic error source is represented by its own covariance matrix. The sum of these matrices describes the total systematic error.

Results
The measured double-differential cross-sections for the production of π + and π − in the laboratory system as a function of the momentum and the polar angle for each incident beam momentum are shown in Fig. 2 and 3, respectively. The error bars shown are the square-roots of the diagonal elements in the covariance matrix, where statistical and systematic uncertainties are combined in quadrature. Correlations cannot be shown in the figures. The correlation of the statistical errors (introduced by the unfolding procedure) are typically smaller than 20% for adjacent momentum bins and smaller for adjacent angular bins. The correlations of the systematic errors are larger, typically 80% for adjacent bins . Tables with the results of this analysis are also given in Appendix A. A discussion of the error evaluation is given below. The overall scale error (< 2%) is not shown. The measurements for the different beam momenta are overlaid in the same figure. For the 3 GeV/c data the point-to-point statistical error is larger than the systematic error, except for the lowest secondary momentum bin. Especially in the middle of the range (around 400 MeV/c), the systematic error is small. Thus the fluctuations between the points are expected to be of statistical nature. In the first angular bins the momentum resolution is relatively large compared to the bin size such that the unfolding procedure tends to display statistical fluctuations over two bins. Since the treatment of the data sets taken with different beam momenta is identical, structures visible in the spectra at 3 GeV/c and not visible in the other data sets are not likely to be artefacts of the efficiency corrections. Overall trends in the shapes, i.e. structures extending over more than two bins are, however, to be considered significant.
To better visualize the dependence on the incoming beam momentum, the same data averaged over the angular range (separately for the forward going and backward going tracks) covered by the analysis are shown separately for π + and π − in Fig. 4. The spectrum of pions produced in the backward direction is much steeper than that in the forward direction.
The increase of the pion yield per proton is visible in addition to a change of spectrum towards higher momentum of the secondaries produced by higher momentum beams in the forward direction.
The dependence of the integrated pion yields on the incident beam momentum is shown in Fig. 5 and compared with the p-Ta data taken with the same apparatus (Ref. [12]). The π + and π − yields integrated over the region 0.350 rad ≤ θ < 1.550 rad and 100 MeV/c ≤ p < 700 MeV/c are shown in the left panel and the data integrated over the region 0.350 rad ≤ θ < 0.950 rad and 250 MeV/c ≤ p < 500 MeV/c in the right panel. The beam energy dependence of the yields is clearly different in the p-C data compared to the p-Ta data. The dependence in the p-C data is much more flat with a saturation of the yield between 8 GeV/c and 12 GeV/c (in both integration regions). The π + and π − production yields exhibit a different behaviour. The p-Cu and p-Sn data are more similar to the p-Ta data than the p-C data indicating a smooth transition between light and heavy target nuclei.
The integrated π − /π + ratio in the forward direction is displayed in Fig. 6 as a function of secondary momentum. The previously published p-Ta data are reproduced in addition to the measurements on the three target nuclei presented in this paper. In the covered part of the momentum range more π + 's are produced than π − 's. The π − /π + ratio increases with increasing beam momentum and, depending on the beam momentum, a change of the sign of the slope of the ratio as a function of secondary momentum is visible in the p-C data. The latter feature is not present in the p-Cu, p-Sn and p-Ta [12] data. The ratio is closer to unity for the heavier target nuclei and a smaller variation with beam momentum is observed.
The dependence of the integrated pion yields on the atomic number A is shown in Fig. 7 combining the results with the p-Ta data (Ref. [12]) taken with the same apparatus and analysed using the same methods. The π + yields integrated over the region 0.350 rad ≤ θ < 1.550 rad and 100 MeV/c ≤ p < 700 MeV/c are shown in the left panel and the π − data integrated over the same region in the right panel for four different beam momenta. One observes a smooth behaviour of the integrated yields. The A-dependence is slightly different for π − and π + production, the latter saturating earlier, especially at lower beam momenta.

Systematic errors
The uncertainties are reported in some detail in Table 2 for the carbon target data and summarized for the copper and tin target data in Table 3. One observes that only for the 3 GeV/c beam is the statistical error similar in magnitude to the systematic error, while the statistical error is negligible for the 8 GeV/c and 12 GeV/c beams. The statistical error is calculated by error propagation as part of the unfolding procedure. It takes into account that the unfolding matrix is obtained from the data themselves 4 and hence contributes also to the statistical error. This procedure almost doubles the statistical error, but    avoids an important systematic error which would otherwise be introduced by assuming a cross-section model a priori to calculate the corrections.
The largest systematic error corresponds to the uncertainty in the absolute momentum scale, which was estimated to be around 3% using elastic scattering (see detailed discussion in [12]). It is difficult to better constrain the absolute momentum scale, since it depends on the knowledge of the beam momentum (known to 1%) and the measurement of the forward scattering angle in the elastic scattering interaction. At low momentum in the relatively small angle forward direction the uncertainty in the subtraction of the electron and positron background due to π 0 production is dominant. This uncertainty is split between the variation in the shape of the π 0 spectrum and the normalization using the recognized electrons. The assumption is made that the π 0 spectrum is similar to the spectrum of charged pions. Initial π − and π + spectra are obtained in an analysis without π 0 subtraction. The π − spectra are then used in the MC for the π 0 distributions. A full simulation of the production and decay into γ's with subsequent conversion in the detector materials is used to predict the background electron and positron tracks. In the region below 120 MeV/c a large fraction of the electrons can be unambiguously identified. These tracks are used as relative normalization between data and MC. The remaining background is then estimated from the distributions of the simulated electron and positron tracks which are accepted as pion tracks with the same criteria as used to select the data. These normalized distributions are subtracted from the data before the unfolding procedure is applied. Uncertainties in the assumption of the π 0 spectrum are taken into account by an alternative assumption that their spectrum follows the average of the π − and π + distribution. An additional systematic error of 10% is assigned to the normalization of the π 0 subtraction using the identified electrons and positrons.  The target region definition and the uncertainty in the PID efficiency and background from tertiaries are of similar size and are not negligible. Relatively small errors are introduced by the uncertainties in the absorption correction, absolute knowledge of the angular and the momentum resolution. The correction for tertiaries (particles produced in secondary interactions) is relatively large at low momenta and large angles. As expected, this region is most affected by this component.
As already mentioned above, the overall normalization has an uncertainty of 2%, and is not reported in the table. It is mainly due to the uncertainty in the efficiency that beam protons counted in the normalization actually hit the target, with smaller components from the target density and beam particle counting procedure.

Comparisons with earlier data
Very few pion production data sets are available in the literature for p-C, p-Cu and p-Sn interactions in this energy region. Our data can be compared with results from Ref. [24] and [25] where measurements of π − production are reported in 4.2 GeV/c and 10 GeV/c p-C interactions, respectively. The total number of π − observed in the above references is about 1300 (5650) in the 4.2(10) GeV/c data. In the papers cited above no tables of the double differential cross-sections were provided, the measurements being given in parametrized and graphical form only. The authors of Ref. [24] and [25] give the results as a simple exponential in the invariant cross-section: E A d 3 σ dp 3 , where E and p are the energy and momentum of the produced particle, respectively, and A the atomic number of the target nucleus 5 . Unfortunately, no absolute normalization is given numerically. To provide a comparison with these data, the parametrization  10.6 (9.9) 7.9 (6.5) 7.0 (5.7) 3.6 (3.5) 4.1 (4.6) 6.1 (6.8) 6.6 (6.5) 11.3 (10. was integrated over the angular bins used in our analysis and with an arbitrary overall normalization overlaid to our results. We compare the 4.2 GeV/c parametrization of Ref. [24] with our 5 GeV/c data and the Ref. [25] parametrization with our 12 GeV/c data. In the comparison with the 4.2 GeV/c parametrization the normalization c is a simple constant, while for the 10 GeV/c parametrization a smooth θ-dependence consistent with a graphical analysis of Ref. [25] was used. Thus only the comparison of the slopes with secondary momentum can be considered significant. Since the 8 GeV/c and 12 GeV/c p-C results are very similar, the lack of data with an exactly equal beam momentum does not play an important role. The results of this comparison are shown in Fig. 8. The shaded band gives the excursion of the parametrization due to the error in the slope parameters (±2σ) with an additional assumed 10% error on the absolute scale. The latter additional error takes into account the fact that the errors on the slopes fitted to the individual angular bins in the cited data are at least a factor of two larger than in the exponential slope obtained from their global parametrization. The agreement of our data with the simple parametrization is good. To judge the quality of the comparison, one should keep in mind that the statistics of Ref. [24] and [25] is much smaller (1300 π − and 5650 π − , respectively) than the statistics of the π − samples in our 5 GeV/c and 12 GeV/c data (18,000 and 43,000 π − , respectively). The errors on the slopes fitted to the individual angular bins in the cited data are at least a factor of two larger than in the exponential slope obtained from their global parametrization. The bands in the figure extend over the region where data from Ref. [24] and Ref. [25] are available.
Our p-C and p-Cu data can also be compared with π + and π − production measurements taken with 12 GeV/c incident protons from Ref. [26]. These data were taken with a magnetic spectrometer and only measurements at 90 degrees from the initial proton direction are available. The statistical point-to-point errors are quoted to be 3%, while the overall normalization has a 30% uncertainty due to the knowledge of the acceptance. In Fig. 9 their p-C data are shown together with the p-C data reported in this paper. The filled boxes show the data directly from Ref. [26], while the open boxes are scaled with a factor 0.72. This factor was defined by scaling the average of the π − and π + data from Ref. [26] at 179 MeV/c and 242 MeV/c to the HARP data averaged over the same region. The scale factor is within one standard deviation of the systematic normalization uncertainty of Ref. [26]. The latter data set compares well with the data described in this paper (filled circles) in the angular region 1.35 rad ≤ θ <1.55 rad at the same proton beam momentum.
In Fig. 10 the p-Cu data of Ref. [26] are shown together with the p-Cu results reported in this paper. The filled boxes show the data directly from Ref. [26], while the open boxes are scaled with a factor 0.91. This factor was defined in a similar procedure as described for the p-C data. The scale factor is very close to unity and well within one standard deviation of the systematic normalization uncertainty of Ref. [26]. Like for the p-C data, the full error bars including the 30% scale uncertainty of Ref. [26] are Figure 8: Comparison of the HARP p-C results with data from Ref. [24] and [25]. The left panel shows the comparison of the parametrization of the 4.2 GeV/c data of Ref. [24] with the 5 GeV/c data reported here; the right panel shows the comparison of the 10 GeV/c parametrization of [25] with the 12 GeV/c data. The absolute normalization of the parametrization was fixed to the data in both cases. The band shows the range allowed by varying the slope parameters given by [24] and [25] with two standard deviation and a 10% variation on the absolute scale. The angular ranges are shown in mrad in the panels.
drawn, for the scaled data only their quoted statistical error of 3% is shown. The agreement of the two data sets is excellent. The fact that the two scale factors are different may be due to the fact that the scale uncertainty in Ref. [26] holds separately for data sets taken with different target nuclei.
Available data at 12.3 GeV/c from the E910 experiment [27] are in reasonable agreement with our p-Cu results as shown in Fig. 11. In order to take into account the different angular binnings which prevent a direct comparison, a Sanford-Wang parametrization is fitted to our data. The fit is performed to the data redefined as d 2 σ π /dpdΩ . An area between two parametrizations is defined which contains our data points as shown in Fig. 11 (top panels). It is visible that the parametrization is not a perfect description to our data. Therefore, we define a band of ±15% around the best fit which contains almost all the HARP data points. The same parametrizations are then displayed in the binning of E910. While the shape of the distributions are similar for both π + and π − in HARP and E910 data sets, the absolute normalizations disagree by 5%-10%. For the individual data sets the systematic errors are between 5% and 10% depending on the range of secondary momentum. Since these errors are correlated between bins, the discrepancy in the π + and π − data separately are of the order of one standard deviation. However, the effects are opposite in π + and π − , giving a 15% difference in the π + /π − ratio between two experiments which is of the order of two standard deviations. This effect may point to an underestimation of systematic effects on the absolute normalization in one of the experiments or in the PID efficiency. Part of the difference is also due to an imperfect parametrization of our data sample. Owing to the symmetry of the HARP TPC, including its trigger counter, we do not expect a large systematic error in the HARP data between π + and π − production cross-sections. Figure 9: Comparison of the HARP results with π + and π − production data at 90 degrees from Ref. [26] taken with 12 GeV/c protons. The left panel shows the comparison of the π + production data of Ref. [26] with the data reported here; the right panel shows the comparison with the π − production data. The smaller filled boxes show the data directly from Ref. [26], while the open boxes are scaled as explained in the text. The latter data set compares well with the data described in this paper (filled circles) in the angular region 1.35 rad ≤ θ <1.55 rad.

Summary and Conclusions
An analysis of the production of pions at large angles with respect to the beam direction for protons of 3 GeV/c, 5 GeV/c, 8 GeV/c and 12 GeV/c impinging on thin (5% interaction length) carbon, copper and tin targets is described. The secondary pion yield is measured in a large angular and momentum range and double-differential cross-sections are obtained. A detailed error estimation has been discussed. Results on the dependence of pion production on the target atomic number A are also presented.
The use of a single detector for a range of beam momenta makes it possible to measure the dependence of the pion yield on the secondary particle momentum and emission angle θ with high precision. The A dependence of the cross-section can be studied using the combination of the present data with the data obtained with tantalum [12].
Very few pion production measurements in this energy range are reported in the literature. The only comparable results found in the literature agrees with the analysis described in this paper. Hadronic production models describing this energy range can now be compared with our new results and, if needed, improved. Data taken with different target materials and beam momenta will be presented in subsequent papers.

Acknowledgements
We gratefully acknowledge the help and support of the PS beam staff and of the numerous technical collaborators who contributed to the detector design, construction, commissioning and operation. In particular, we would like to thank G. Barichello Figure 10: Comparison of the HARP data with π + and π − production data at 90 degrees from Ref. [26] taken with 12 GeV/c protons. The left panel shows the comparison of the π + production data of Ref. [26] with the data reported here; the right panel shows the comparison with the π − production data. The smaller filled boxes show the data directly from Ref. [26], while the open boxes are scaled as explained in the text. The latter data set compares well with the data described in this paper (filled circles) in the angular region 1.35 rad ≤ θ <1.55 rad. . We gratefully acknowledge their support. This work was supported in part by the Swiss National Science Foundation and the Swiss Agency for Development and Cooperation in the framework of the programme SCOPES -Scientific co-operation between Eastern Europe and Switzerland. Figure 11: Comparison of the HARP data with π + and π − production data from Ref. [27] taken with 12.3 GeV/c protons. The top panels show a parametrization of the π + (left) and π − (right) production data described in this paper. The data have been normalized to represent d 2 σ π /dpdΩ . The shaded band represents the area between two parametrization which contain the data points. The bottom panels show the comparison of the same parametrization, now binned according to the E910 data. The bottom left (right) panel shows the π + (π − ) production data of Ref. [27]. The angular regions are indicated in mrad in the upper right-hand corner of each plot.
A Cross-section data Table 4: HARP results for the double-differential π + production cross-section in the laboratory system, d 2 σ π + /(dpdθ) for carbon. Each row refers to a different (p min ≤ p < p max , θ min ≤ θ < θ max ) bin, where p and θ are the pion momentum and polar angle, respectively. The central value as well as the square-root of the diagonal elements of the covariance matrix are given.  Table 5: HARP results for the double-differential π − production cross-section in the laboratory system, d 2 σ π − /(dpdθ) for carbon. Each row refers to a different (p min ≤ p < p max , θ min ≤ θ < θ max ) bin, where p and θ are the pion momentum and polar angle, respectively. The central value as well as the square-root of the diagonal elements of the covariance matrix are given.  Table 6: HARP results for the double-differential π + production cross-section in the laboratory system, d 2 σ π + /(dpdθ) for copper. Each row refers to a different (p min ≤ p < p max , θ min ≤ θ < θ max ) bin, where p and θ are the pion momentum and polar angle, respectively. The central value as well as the square-root of the diagonal elements of the covariance matrix are given.  Table 7: HARP results for the double-differential π − production cross-section in the laboratory system, d 2 σ π − /(dpdθ) for copper. Each row refers to a different (p min ≤ p < p max , θ min ≤ θ < θ max ) bin, where p and θ are the pion momentum and polar angle, respectively. The central value as well as the square-root of the diagonal elements of the covariance matrix are given.  Table 8: HARP results for the double-differential π + production cross-section in the laboratory system, d 2 σ π + /(dpdθ) for tin. Each row refers to a different (p min ≤ p < p max , θ min ≤ θ < θ max ) bin, where p and θ are the pion momentum and polar angle, respectively. The central value as well as the square-root of the diagonal elements of the covariance matrix are given.
θ min θ max p min p max d 2 σ π + /(dpdθ) (rad) (rad) (GeV/c) (GeV/c) (barn/(GeV/c rad)) 3 GeV/c   Table 9: HARP results for the double-differential π − production cross-section in the laboratory system, d 2 σ π − /(dpdθ) for tin. Each row refers to a different (p min ≤ p < p max , θ min ≤ θ < θ max ) bin, where p and θ are the pion momentum and polar angle, respectively. The central value as well as the square-root of the diagonal elements of the covariance matrix are given.