Tuning the GENIE Pion Production Model with MINERvA Data
P. Stowell, L. Pickering, C. Wilkinson, C.V.C. Wret, F. Akbar, D.A. Andrade, M. V. Ascencio, L. Bellantoni, A. Bercellie, M. Betancourt, A. Bodek, A. Bravar, H. Budd, G. Caceres, T. Cai, M.F. Carneiro, J. Chaves, H. da Motta, S.A. Dytman, G.A. Dıaz, J. Felix, L. Fields, A. Filkins, R. Fine, N. Fiza, R.Galindo, H. Gallagher, A. Ghosh, R. Gran, D.A. Harris, S. Henry, S. Jena, D. Jena, J. Kleykamp, M. Kordosky, D. Last, T. Le, X.-G. Lu, E. Maher, S. Manly, W.A. Mann, C. M. Marshall, K. S. McFarland, B. Messerly, J. Miller, J.G. Morfın, J. Mousseau, D. Naples, J.K. Nelson, C. Nguyen, A. Norrick, Nuruzzaman, V. Paolone, G.N. Perdue, M.A. Ramırez, R.D. Ransome, H. Ray, D. Rimal, P.A. Rodrigues, D. Ruterbories, H. Schellman, C.J. Solano Salinas, H. Su, M. Sultana, V. S. Syrotenko, E. Valencia, J. Wolcott, M. Wospakrik, B. Yaeggy, L. Zazueta, D. Zhang
TTuning the GENIE Pion Production Model with MINERvA Data
P. Stowell, L. Pickering,
2, 3
C. Wilkinson, C. Wret,
5, 3
F. Akbar, D. A. Andrade, M. V. Ascencio, L. Bellantoni, A. Bercellie, M. Betancourt, A. Bodek, A. Bravar, H. Budd, G. Caceres, T. Cai, M.F. Carneiro, J. Chaves, H. da Motta, S. A. Dytman, G.A. D´ıaz,
5, 8
J. Felix, L. Fields,
9, 15
A. Filkins, R. Fine, N. Fiza, R. Galindo, H. Gallagher, A. Ghosh,
18, 11
R. Gran, D. A. Harris, S. Henry, S. Jena, D. Jena, J. Kleykamp, M. Kordosky, D. Last, T. Le,
19, 21
X.-G. Lu, E. Maher, S. Manly, W. A. Mann, C. M. Marshall, ∗ K. S. McFarland,
5, 9
B. Messerly, J. Miller, J. G. Morf´ın, J. Mousseau, † D. Naples, J. K. Nelson, C. Nguyen, A. Norrick, Nuruzzaman,
21, 18
V. Paolone, G. N. Perdue,
9, 5
M. A. Ram´ırez, R. D. Ransome, H. Ray, D. Rimal, P. A. Rodrigues,
25, 5
D. Ruterbories, H. Schellman,
12, 15
C. J. Solano Salinas, H. Su, M. Sultana, V.S. Syrotenko, E. Valencia,
16, 7
J. Wolcott, ‡ M. Wospakrik, B. Yaeggy, L. Zazueta, and D. Zhang (The MINER ν A Collaboration) University of Sheffield, Dept. of Physics and Astronomy, Sheffield, United Kingdom Michigan State University, Dept. of Physics and Astronomy, East Lansing, Michigan 48824, USA Imperial College London, Dept. of Physics, London, United Kingdom University of Bern, Albert Einstein Center for Fundamental Physics, LHEP, Bern, Switzerland University of Rochester, Rochester, New York 14627 USA AMU Campus, Aligarh, Uttar Pradesh 202001, India Campus Le´on y Campus Guanajuato, Universidad de Guanajuato, Lascurainde Retana No. 5, Colonia Centro, Guanajuato 36000, Guanajuato M´exico. Secci´on F´ısica, Departamento de Ciencias, Pontificia Universidad Cat´olica del Per´u, Apartado 1761, Lima, Per´u Fermi National Accelerator Laboratory, Batavia, Illinois 60510, USA University of Geneva, 1211 Geneva 4, Switzerland Centro Brasileiro de Pesquisas F´ısicas, Rua Dr. Xavier Sigaud 150, Urca, Rio de Janeiro, Rio de Janeiro, 22290-180, Brazil Dept. of Physics, Oregon State University, Corvallis, Oregon 97331, USA Dept. of Physics and Astronomy, University of Pennsylvania, Philadelphia, PA 19104 Dept. of Physics and Astronomy, University of Pittsburgh, Pittsburgh, Pennsylvania 15260, USA Northwestern University, Evanston, Illinois 60208 Dept. of Physics, College of William & Mary, Williamsburg, Virginia 23187, USA IISER, Mohali, Knowledge city, Sector 81, Manauli PO 140306 Departamento de F´ısica, Universidad T´ecnica Federico Santa Mar´ıa, Avenida Espa˜na 1680 Casilla 110-V, Valpara´ıso, Chile Physics Dept., Tufts University, Medford, Massachusetts 02155, USA Dept. of Physics, University of Minnesota – Duluth, Duluth, Minnesota 55812, USA Rutgers, The State University of New Jersey, Piscataway, New Jersey 08854, USA Oxford University, Dept. of Physics, Oxford, United Kingdom Massachusetts College of Liberal Arts, 375 Church Street, North Adams, MA 01247 University of Florida, Dept. of Physics, Gainesville, FL 32611 University of Mississippi, Oxford, Mississippi 38677, USA Universidad Nacional de Ingenier´ıa, Apartado 31139, Lima, Per´u (Dated: October 3, 2019)Faced with unresolved tensions between neutrino interaction measurements at few-GeV neutrinoenergies, current experiments are forced to accept large systematic uncertainties to cover discrep-ancies between their data and model predictions. In this paper, the widely used pion productionmodel in GENIE is compared to four MINERvA charged current pion production measurementsusing NUISANCE. Tunings, i.e. , adjustments of model parameters, to help match GENIE to MIN-ERvA and older bubble chamber data are presented here. We find that scattering off nuclear targetsas measured in MINERvA is not in good agreement with expectations based upon scattering offnucleon (hydrogen or deuterium) targets in existing bubble chamber data. An additional ad hoccorrection for the low- Q region, where collective nuclear effects are expected to be large, is pre-sented. While these tunings and corrections improve the agreement of GENIE with the data, themodeling is imperfect. The development of these tunings within the NUISANCE framework allowsfor straightforward extensions to other neutrino event generators and models, and allows omittingand including new data sets as they become available. ∗ now at Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA a r X i v : . [ h e p - e x ] O c t † now at University of Michigan, Ann Arbor, MI 48109, USA ‡ now at Tufts University, Medford, MA 02155, USA I. INTRODUCTION
In recent years, experimental groups have started topublish neutrino interaction cross-section measurementson nuclear targets in terms of measurable final state par-ticle content, instead of inferred initial interaction chan-nels. This avoids the problem of correcting for complexnuclear effects to make a measurement in terms of the ini-tial interaction channels. For example, events with onlya single pion can be produced by the decay of hadronicresonances formed at the primary neutrino interaction,followed by loss of a nucleon from the resonance’s de-cay as a result of final state interactions (FSI) withinthe nuclear medium. Such events can also be producedby other sequences of interactions, such as a deep inelas-tic collision where only a single pion is produced afterFSI. A measurement of charged current events with oneidentified pion in the final state is a benchmark for mod-els, independent of the details of how each model assessesany particular interaction channel’s contributions to thatfinal state. The limitation of giving results in terms of fi-nal state particle content is that FSI are important, andresult in the contribution of many different interactionchannels into a specific final state.There are tensions between published data from theT2K, MiniBooNE, and MINERvA experiments [1–5].These tensions exist in the charged current production ofboth zero and one pion final states, and a model has yetto emerge that can reliably simulate all experiments atonce. This is troubling, as current and future neutrino os-cillation experiments require a cross section model whichis predictive across the range of energies covered by theseexperiments and for a variety of targets.The differences in neutrino fluxes, scattering targets,available phase space and signal definitions between ex-periments make it difficult to diagnose the exact causes ofdisagreement within the global data set. In particular, asresults must be averaged over the neutrino flux distribu-tion of each experiment, it is difficult to disentangle theenergy dependence of an observed deficiency in a particu-lar model, and decide how uncertainties should be prop-agated in neutrino energy. Tensions between measure-ments from a single experiment can uncover fundamentalproblems with a model which should be addressed, be-fore considering the more difficult issue of developing, orempirically tuning a model which fits data from multipleexperiments. In this work, we employ published MIN-ERvA pion production data. The cross-section measure-ments utilized in this effort have not been reanalyzed ormodified in any way.NUISANCE [6] was developed to provide the neutrinoscattering community with a flexible framework in whichvarious neutrino interaction generators can be validatedand empirically tuned to data. Its structure allows forgenerator tunings to be easily adapted to account forchanges in the underlying model or data. In this work,the default pion production model in the GENIE [7, 8]neutrino interaction generator is tuned to MINERvA data. Although more sophisticated pion production mod-els exist (e.g. [9–12]), GENIE is widely used by theneutrino scattering community, and its model uncertain-ties have a central importance to the field. Althoughthe work is only directly applicable to one generator,the methods developed in this paper are easily adapt-able to different generators. All the data and methodsare publicly available and integrated into the open sourceNUISANCE framework, facilitating similar studies usingother generators and models.In Section II, the data are reviewed and the goodness-of-fit test statistic is defined for the tuning process. Sec-tion III describes the default GENIE pion productionmodel, and reviews comparisons of this model to data.In Section IV, the parameter reweighting package inGENIE is discussed along with the specific parameterstuned therewith. We also discuss other corrections to theGENIE model made to improve agreement with bubblechamber data [13, 14]. In Section V, we tune additionalsystematic parameters in GENIE to improve agreementwith the MINERvA data in combination with the bubblechamber data. In Section VI, additional low- Q ad hoccorrections are added to the model to resolve observedtensions, motivated by the need for similar correctionsobserved at both MINOS [15] and MiniBooNE [16]. Fi-nally, in Section VII we present our conclusions. II. DATA INCLUDED IN THE FITS
We tune to four of MINERvA’s published chargedcurrent pion production measurements taken ona polystyrene scintillator target: ν µ CC1 π ± [17], ν µ CC N π ± [18], ν µ CC1 π [19] and ¯ ν µ CC1 π [18], sum-marized in Table I . The MINERvA detector [20] doesnot determine the polarity of charged pions. The fractionof π − in ν µ CC1 π ± sample is small ( ∼ ν µ CC1 π ± and ν µ CC N π ± signal definition allows forany number of neutral pions. Approximately 3% of theMINERvA ν µ CC1 π ± signal events have at least one neu-tral pion in the final state. All four analyses include sig-nal definition cuts on the true “reconstructed” mass ofthe hadronic system assuming the struck nucleon is atrest, W rec , and the true neutrino energy E ν .The kinematic variable distributions used in this workare the momentum and angle of the outgoing muon withrespect to the incoming neutrino beam, p µ and θ µ , andthe kinetic energy and angle of the outgoing pion withrespect to the incoming neutrino beam, T π and θ π . Inthe ν µ CC N π ± channel, where there is at least one π ± in the final state, there is one entry in the distribu-tions of θ π and T π for each π ± in an event. The dataare reported as efficiency corrected results unfolded to In “ ν µ CC Nπ ± ”, the N indicates one or more identified pionsand does not refer to a nucleon. Channel ν µ CC1 π ± [17] ν µ CC Nπ ± [18] ν µ CC1 π [19] ¯ ν µ CC1 π [18]N bins p µ bins θ µ bins T π bins θ π
14 14 11 11N bins total 38 39 35 36Signal definition 1 π ± , ≥ π > π ± , ≥ π π , 0 π ± π , 0 π ± µ − µ − µ − µ + W rec < . W rec < . W rec < . W rec < . θ µ < ◦ —TABLE I. Summary of the measurements used in this analysis. W rec is the true reconstructed hadronic mass assuming thestruck nucleon is at rest. None of the measurements veto on activity other than the µ and π in their signal definition, and allselections require 1 . < E ν <
20 GeV. true kinematic variables, which may introduce modeldependence. This is notably problematic in regions oflow efficiency—present in the charged pion channels at θ π ∼ ◦ , T π <
50 MeV and T π >
350 MeV, wherethe signal efficiency is zero [17]. The pion selection cuts,not present in the signal definition, remove about 50%of signal events, with little dependence upon the muonvariables, but a clear impact on the shape of the pionkinematic variables.The published cross-sections are one dimensional withcorrelations provided between the bins within each dis-tribution. No correlations are provided between mea-surements of different final states, or between differentone-dimensional projections of the same measurement.These correlations are expected to be large, coming pre-dominantly from flux and detector uncertainties. Addi-tionally, the ν µ CC1 π ± event sample is a subset ( ∼ ν µ CC N π ± event sample, and including both channelsintroduces a statistical correlation. Not assessing corre-lations between the distributions, while common practicein this field, is a limitation when tuning models to mul-tiple data sets. It introduces a bias in the χ statisticthat is difficult to quantify, and requires imposing ad hocuncertainties [4] as the test-statistic is not expected tofollow a χ distribution for the given degrees of freedom.The covariance matrices contain a flux-dominated nor-malization component which we expect to be fully cor-related across all distributions. To account for the cor-related uncertainty, we use the full covariance matrix, M ij , for the p µ distribution and shape-only covariancematrices, S ij , for the other three distributions in eachof the topologies. Whilst any distribution could set thenormalization constraint, the shape of the p µ distribu-tion for each channel was chosen since it was found to berelatively insensitive to model variations and had goodshape agreement with the data. The joint χ is thereforedefined as the sum of the full p µ χ and shape-only θ µ , T π and θ π χ ’s: χ = N pµ (cid:88) ij ∆ i ( M − ) ij ∆ j + N k (cid:88) kij ∆ Sk,i ( S − ) ij ∆ Sk,j (1)where i and j are bin indices,∆ i = d p µ ,i − m p µ ,i (2)∆ Sk,i = d k,i − (cid:32) m k,i × (cid:80) j d k,j (cid:80) j m k,j (cid:33) , (3)and d k,i and m k,i are the data and MC values respectivelyfor the i th bin in the k th distribution. The shape-only co-variance matrices are provided in the public data releasefor the ν µ CC1 π ± and ν µ CC1 π measurements, and themethod of Ref. [21] (section 10.6.3) was used to extractthem for the ν µ CC N π ± and ¯ ν µ CC1 π channels. III. PION PRODUCTION IN GENIE
This analysis begins with version 2.12.6 of GENIE,which is close to what is used by MINERvA, T2K, NOvAand MicroBooNE. We use the Smith-Moniz relativis-tic Fermi gas (RFG) model [22] with an added highmomentum tail as per Bodek and Ritchie [23]. TheValencia random phase approximation screening [24] isapplied as a weight to quasielastic events. The two-particle two-hole process is simulated using the Valen-cia model [11, 25]. MINERvA currently uses a modifi-cation of v2.8.4 [19, 26, 27] with an increased rate forthe Valencia two-particle two-hole process; that modifi-cation is not used here. An important difference in sin-gle pion production between v2.8.x and v2.12.x is theangular distributions of single pion events in the Rein-Sehgal model, discussed below. A sample of 2.5 millionevents were generated using the MINERvA flux predic-tions [28], a polystyrene target and the official GENIE2.12.6 splines [29].To simulate pion production, GENIE uses the Rein-Sehgal (RS) model [30] with a hadronic invariant masscut of W ≤ . N (1990) were not includeddue to their unclear experimental status at the time ofimplementation. Resonance-resonance and resonance-nonresonance interference terms are not included. Lep-ton mass terms are only included in calculating phasespace limits and are neglected when calculating thecross sections. A discussion of the limitations of thissimplification can be found in Ref. [31]. In earlierversions—including v2.8.4—the pion-nucleon distribu-tion was isotropic in the resonance rest frame, but waschanged in 2.12.x. Here we use the non-isotropic modelas our default and reweight to the isotropic distribu-tion, explained later. The RS nonresonant backgroundis not used by GENIE; rather, a deep inelastic scattering(DIS) model is extended to cover that invariant mass re-gion. The DIS model uses the Bodek-Yang parametriza-tion [32], and the AGKY model to describe hadroniza-tion [33]. In the AGKY model, the KNO model [34] isused for W ≤ . W ≥ . ∼ T π <
500 MeV to agree with data [27]. This correc-tion also moves the shape of the T π spectrum closer tothe predictions of the Berger-Sehgal coherent model [38].The ν µ CC1 π ± channel has a small contribution from co-herent production in the lowest Q bins but the inclusionof this suppression has only a small effect on the MC pre-dictions. To maintain a model similar to that currentlybeing used by MINERvA, this suppression is included inthe analysis presented in Sections IV and onwards.The “hA Intranuke” effective cascade model [39] isused to model pion and nucleon FSI. In this model, theeffect of intranuclear scattering is parameterized as asingle cascade step applied to each particle emanatingfrom the primary interaction. This model steps hadronsthrough a nucleus of radius r ∼ A / and a nuclear den-sity function derived from electron scattering data. Thehadron’s mean free path is determined from tabulatedhadron-proton and hadron-neutron cross sections [40].The probability to interact with the nucleus is high; itis e.g. ∼
73% for a pion from an E ν = 3 GeV quasielasticevent in carbon. When a FSI occurs, the possible inter-actions (absorption, pion production, knockout, chargeexchange, elastic scatter) are chosen according to their proportions for iron.Default GENIE predictions separated by nucleon levelinteraction channels for the MINERvA data are shown inFig. 1. The shape of the p µ distributions agree well withthe data for all four measurements. However, the modeloverestimates the cross section for π ± production and asa result the χ for the ν µ CC1 π ± and ν µ CC N π ± , givenin the fourth column (“Default”) of Table II, are large.The model overestimates θ µ below < ◦ in the π chan-nels, although it does correctly predict the shape of the θ µ distribution in the π ± channels. The model under-estimates the production rate at large θ µ in ν µ CC1 π .The shape of the T π distribution is in larger disagree-ment for ν µ CC N π ± data than for ν µ CC1 π ± . Since the ν µ CC N π ± distributions summed over all identified π ± ,redistributing kinetic energy between π ± in events withmore than one π ± could resolve some of this tension.The π channels are under-predicted at low T π . Fi-nally, GENIE predictions are too high in magnitude at θ π ≈ ◦ in both the ν µ CC1 π ± and ν µ CC N π ± channels,and the prediction has the wrong shape in the ν µ CC1 π ± channel. Comparisons using the transport theory basedGiBUU model [41] show similar shape disagreements de-spite GiBUU’s use of an advanced semiclassical cascademodel to simulate FSI [1].Each of the measurements are shown as MC/data ra-tio distributions in Fig. 2. Similar comparisons betweenthe MiniBooNE and MINERvA experiments are foundin Ref. [5]. The shape-only data sets ( θ µ , θ π , T π ) werenormalized to match the data before the ratio was takenand the error bars in Fig. 2 reflect the extracted shape-only uncertainties on the data, so that the distributionsreflect their contributions to the total χ . IV. TUNABLE PARAMETERS IN THE GENIEMODEL
The GENIE event generator allows assessment of sys-tematic uncertainties through the GENIE reweightingpackage. A large number of event weighting “dials” areincluded to allow model uncertainties to be evaluated.The dials adjusted in this note are summarized in Ta-ble III and are chosen because of their connection to thekinematic variables and interaction modes studied herein.Experiments often use variations in the charged-current resonant axial mass, M resA , as a systematic uncer-tainty which varies both the normalization and Q shapeof resonant interactions along with variations in a to-tal resonant cross-section normalization dial, NormRes.Variations in NormRes approximates the behaviour ofvarying F A (0) in the axial form factor in the Rein-Sehgalmodel. Since low θ µ correlates with low Q , variations in M resA have the largest effect on the shape of the muon an-gular distributions as shown in Fig. 3, and have a smalleffect on the θ π spectrum.Dials are available to vary the normalization of thenonresonant 1 π production channels in GENIE ( e.g. Distribution Channel N bins
Default ANL/BNL FrAbs Tune FrInel Tune p µ (Rate) ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ν µ CC1 π θ µ (Shape) ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ν µ CC1 π T π (Shape) ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ν µ CC1 π θ π (Shape) ν µ CC1 π ±
14 25.4 26.5 13.0 12.6 ν µ CC Nπ ±
14 11.7 11.1 6.9 6.2 ν µ CC1 π
11 13.5 15.0 8.3 8.9¯ ν µ CC1 π
11 5.7 5.9 3.4 3.5Total χ
148 275.6 312.7 242.3 240.7TABLE II. Channel by channel contributions to the χ at different stages of the tuning process.Parameter Default Value GENIE parameter nameCC Resonant Axial Mass ( M resA ) 1 . ± .
22 GeV
MaCCRES
CC Resonant Normalization (NormRes) 100 ±
20 %
NormCCRES
CC1 π Nonresonant Normalization (NonRes1 π ) 100 ±
50 %
NonRESBGvnCC1piNonRESBGvpCC1piNonRESBGvbarnCC1piNonRESBGvbarpCC1pi
CC2 π Nonresonant Normalization (NonRes2 π ) 100 ±
50 %
NonRESBGvnCC2piNonRESBGvpCC1piNonRESBGvbarnCC1piNonRESBGvbarpCC1pi
Pion Angular Emission ( π -iso) 0 (RS) Theta Delta2Npi
Pion Absorption FSI Fraction (FrAbs) 100 ±
30 %
FrAbs pi
Pion Inelastic FSI Fraction (FrInel) 100 ±
40 %
FrInel pi
TABLE III. Summary of the GENIE dials optimized in this note, their default values and the uncertainties recommended bythe GENIE collaboration. We do not use the defaults for M resA , NormRes and NonRes1 π and instead impose central valuesand uncertainties from tunings to ANL and BNL data as described in the text. NonRESBGvnCC1pi , NonRESBGvpCC1pi ) but each dial in-troduces similar modifications to the predictions. Toreduce the number of free parameters in the fit de-scribed in Section V, these dials were grouped into a sin-gle background scaling for nonresonant 1 π production,NonRes1 π , following the approach in Ref. [13, 14]. Asimilar treatment was also applied to nonresonant 2 π production, NonRes2 π , with the neutrino and antineu- trino related parameters assumed to be 100% correlatedin both cases. The effects of varying the nonresonantcontributions are shown in Fig. 4. Variations in theNonRes2 π dial introduce a large change in normaliza-tion for the ν µ CC N π ± channel and has a minor effect inthe other single pion channels as the fraction of multi- π events is small.Reanalysis of data from ANL and BNL bubble cham- /Nbin) c DATA (MC Shape p CC1 on p + p CC1 on n + p CC1 p CCN Other1p1h/2p2h (GeV) m P / G e V ) ( c m m / dp s d - · (19.1/8) + p CC1 m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (7.1/9) + p CC1 m n (MeV) p T
100 200 300 / M e V ) ( c m p / d T s d - · (2.9/7) + p CC1 m n (degrees) p q / d e g r ee s ) ( c m pq / d s d - · (25.4/14) + p CC1 m n (GeV) m P / G e V ) ( c m m / dp s d - · (35.4/9) + p CCN m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (4.5/9) + p CCN m n (MeV) p T
100 200 300 / M e V ) ( c m p / d T p ) d N F ( / T - · (39.8/7) + p CCN m n (degrees) p q / d e g r ee s ) ( c m pq / d p ) d N F ( / T - · (11.7/14) + p CCN m n (GeV) m P / G e V ) ( c m m / dp s d - · (11.1/8) p CC1 m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (35.1/9) p CC1 m n (MeV) p T / M e V ) ( c m p / d T s d - · (28.3/7) p CC1 m n (degrees) p q / d e g r ee s ) ( c m pq / d s d - · (13.5/11) p CC1 m n (GeV) m P / G e V ) ( c m m / dp s d - · (7.4/9) p CC1 m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (9.3/9) p CC1 m n (MeV) p T / M e V ) ( c m p / d T s d - · (19.3/7) p CC1 m n (degrees) p q / d e g r ee s ) ( c m pq / d s d - · (5.7/11) p CC1 m n FIG. 1. GENIE 2.12.6 Default model predictions compared to MINERvA data. Colors correspond to particle content at thenucleon interaction. “Other” is dominated by coherent pion production. “MC Shape” shows the total MC prediction afterit has been normalized to match the total data normalization. In the case of the shape-only distributions ( θ µ , T π , θ π ) theshape-only χ /N bins values are shown. All cross sections are per nucleon. bers has provided a tuning of GENIE’s single pion pro-duction model on free nucleons. The work showed thata small shift in M resA was required to model the low- Q region and a large suppression of the nonresonant π pro-duction ( − π + and π production. The reanalysis usedthe measured ratios of the rates of single π productionto CCQE measurements to cancel errors in the flux. We note that by using CCQE data multiple times, they in-troduce hidden correlations which may have a small ef-fect on the postfit uncertainties. However, as the singlepion statistical uncertainties at ANL [42] and BNL [43]were magnitudes higher than the CCQE statistical un-certainty [44, 45], the effect was neglected in that work,and is also neglected here. The resulting parameter tunesshown in Table IV and Fig. 5 have been partially adopted (GeV) m p M C / D a t a R a ti o ( R A TE ) + p CC1 m n + p CCN m n p CC1 m n p CC1 m n (deg.) m q M C / D a t a R a ti o ( S HA P E ) (MeV) p T M C / D a t a R a ti o ( S HA P E ) (deg.) p q M C / D a t a R a ti o ( S HA P E ) FIG. 2. MC/data ratios for the default GENIE predictions. The p µ distribution provides a rate comparison in the χ calculation;the other distributions are treated as shape-only, i.e. , the MC is normalized to match the data and the uncertainties are fromthe shape-only covariance matrix. by MINERvA and NOvA which both apply the nonres-onant rescaling of 43% but leave the other parametersunchanged. Parameter GENIE default ANL/BNL tune M resA [GeV] 1 . ± .
22 0 . ± . ±
20 115 ± π [%] 100 ±
50 43 ± Fig. 6 shows MINERvA data and the predictions ofGENIE when its output has been reweighted to re-flect the parameter changes of Table IV. The channel-by-channel contributions to the χ are given in thefifth column (“ANL/BNL”) of Table II. Incorporatingthe parameter changes improves the total normaliza-tion agreement in the p µ distributions for ν µ CC1 π ± and ν µ CC N π ± . The χ for the p µ distribution is alsoimproved in the ¯ ν µ CC1 π channel, even though theANL/BNL data is from neutrino interactions only. The χ for the p µ distribution in the ν µ CC1 π channel issomewhat worse as the parameter tunes reduce the pre-dicted nucleon ν µ CC1 π cross section. The modifica-tion of M resA shifts the θ µ predictions to lower values,increasing the χ contributions. The T π and θ π distribu-tions change mostly by normalization, having a smallereffect on the χ . The overall agreement of GENIE withMINERvA data is not improved by incorporating theANL/BNL information. Indeed, the total χ increases,largely because of the χ contributions from the θ µ dis-tributions.GENIE provides a dial that influences the resonances’decay into the pion-nucleon system in the resonance restframe, π -iso, and allows events to be reweighted continu-ously between the default anisotropic distribution ( π -iso= 0) and the isotropic distribution ( π -iso = 1). The µ θ / nu c l e on ) ( c m µ θ / d σ d − × Data + π CC1 µ ν Nominal Variations σ MaCCRES 1 Shape Var. σ MaCCRES 1
FIG. 3. The effect varying the M resA dial on the default GE-NIE prediction for θ µ . The red bands show the variation tothe total rate and shape. The green bands are obtained bynormalizing the reweighted curves to the default predictedrate to highlight the smaller effect the dial has on the shapeof the distributions. Adler angle is highly sensitive to the π -iso parameterand has been measured by neutrino induced pion produc-tion experiments on single nucleons, such as ANL [42],BNL [43], BEBC [47, 48] and FNAL [49]. Nucleondata strongly prefers an anisotropic process, as shownin Fig. 7. Nonetheless, π -iso has some impact, albeit onethat does depend on how FSI are modelled, on the shapeof MINERvA θ π and T π distributions, seen in the bottomof Fig. 7, and was therefore included in this work.The GENIE hA model for FSI has uncertainties fromthe π − A cross-section data to which the model wastuned. The total π − A cross section has a strongerconstraint than each of the individual interaction crosssections, so GENIE provides dials to vary the fractionalcontribution of each component. The available fractionaldials are pion absorption (FrAbs), pion inelastic scatter-ing (FrInel), pion elastic scattering, pion charge exchangeand pion production. V. TUNING THE GENIE MODEL
Fig. 6 and Table II show the unsatisfactory agreementof the GENIE prediction against MINERvA data. Thedisagreement worsens after incorporating the prior con-straint from ANL and BNL data; this correction, basedon nucleon data, is inadequate. This section describesfits that improve the agreement with MINERvA data. The Adler angle is the angle between the pion and the threemomentum transfer in the resonance rest frame [46]. (GeV) p T / nu c l e on / G e V ) ( c m p / d T s d - · Data p CC1 m n Nominal Variations s p NonRES1 Shape Var. s p NonRES1 (degrees) π θ / d e g r ee s / nu c l e on ) ( c m π θ / d π ) d N Φ ( / T − × Data + π CCN µ ν Nominal Variations σ π NonRES2 Shape Var. σ π NonRES2
FIG. 4. The effect of varying the NonRes1 π dial on the defaultGENIE prediction for T π (top) and of varying the NonRes2 π dial on the prediction for θ π (bottom). -1-0.8-0.6-0.4-0.200.20.40.60.81 RESA M DIS RES
RESA M DISRES 1.0 -0.3 -0.9-0.3 1.0 0.2-0.9 0.2 1.0
FIG. 5. Correlation matrix from tuning GENIE to reproducethe ANL/BNL pion production measurements included in our χ penalty term. . /Nbin) c DATA (MC ShapeNominal MC Shape p CC1 on p + p CC1 on n + p CC1 p CCN Other1p1h/2p2h (GeV) m P / G e V ) ( c m m / dp s d - · (13.8/8) + p CC1 m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (12.4/9) + p CC1 m n (MeV) p T
100 200 300 / M e V ) ( c m p / d T s d - · (2.6/7) + p CC1 m n (degrees) p q / d e g r ee s ) ( c m pq / d s d - · (26.5/14) + p CC1 m n (GeV) m P / G e V ) ( c m m / dp s d - · (19.5/9) + p CCN m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (10.4/9) + p CCN m n (MeV) p T
100 200 300 / M e V ) ( c m p / d T p ) d N F ( / T - · (34.7/7) + p CCN m n (degrees) p q / d e g r ee s ) ( c m pq / d p ) d N F ( / T - · (11.1/14) + p CCN m n (GeV) m P / G e V ) ( c m m / dp s d - · (19.6/8) p CC1 m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (71.5/9) p CC1 m n (MeV) p T / M e V ) ( c m p / d T s d - · (31.4/7) p CC1 m n (degrees) p q / d e g r ee s ) ( c m pq / d s d - · (15.0/11) p CC1 m n (GeV) m P / G e V ) ( c m m / dp s d - · (6.4/9) p CC1 m n (degrees) m q / d e g r ee s ) ( c m mq / d s d - · (14.0/9) p CC1 m n (MeV) p T / M e V ) ( c m p / d T s d - · (17.9/7) p CC1 m n (degrees) p q / d e g r ee s ) ( c m pq / d s d - · (5.9/11) p CC1 m n FIG. 6. GENIE ANL/BNL single pion tuning model predictions compared to MINERvA data. The distributions have beenweighted to the ANL/BNL tuning parameter set, and have had the coherent pion correction applied. Colors correspond toparticle content at the nucleon interaction. “Other” is dominated by coherent pion production. “MC Shape” shows the totalMC prediction after it has been normalized to match the total data normalization. In the case of the shape-only distributions( θ µ , T π , θ π ) the shape-only χ /N bins values are shown. All cross sections are per nucleon. The parameters M resA , NormRes and NonRes1 π are in-cluded in the fits with a penalty term added to the χ from the ANL and BNL data. The penalty term uses thecovariance, M , shown in Fig. 5: χ = N =3 (cid:88) i,j ( x i − f i ) (cid:0) M − (cid:1) ij ( x j − f j ) , (4)where x i are the parameter values i at each iterationof the fit, and f i are the parameter values from the fit1 Adler q cos - - N u m b e r o f e v e n t s
1p Data + p CC1 m n BNL RS EjectionIsotropic Ejection p q / nu c l e on ) ( c m pq / d s d - · Data + p CC1 m n RS EjectionIsotropic Ejection
FIG. 7. Effect of varying the π -iso dial on pion angular dis-tributions for BNL ν µ CC1 π +
1p data (top) and MINERvA ν µ CC1 π + (bottom). to ANL and BNL data. The GENIE default model isstrongly disfavoured with χ = 299.3, but changingNonRes1 π to 43% while leaving all the other parametersat their default values reduces the χ to 21.8, showingthat the largest tension is due to the NonRes1 π parame-ter.The π -iso dial is allowed to vary in the range 0–1 inthe fit, corresponding to a continuous variation betweenan RS angular distribution and an isotropic distribu-tion for ∆(1232) decay. To avoid the normalization ofthe ν µ CC N π ± measurement pulling parameters in the ν µ CC1 π ± model, the NonRes2 π dial was allowed to varybetween 0-300% of the nominal value.When varying one of the five hA pion FSI dials, GE-NIE automatically adjusts the remaining parameters topreserve the total pion cross-section and maintain agree-ment with pion-nucleus scattering data. This “cushion”technique introduced instabilities in the χ surface, so itwas not possible to include multiple pion FSI parametersin a simultaneous fit. Instead we performed fits with only one of the FSI parameters floating. No χ penalty termswere added for the FSI dial in either tuning: the param-eters were driven solely by MINERvA data. The chargeexchange and pion production dials had small contribu-tions to the overall χ for the selected data, forcing theparameters to be inflated beyond +3 σ of GENIE’s rec-ommendation, with large post-fit uncertainties. Further-more, the pion elastic scattering parameter is stronglyconstrained by external data, so its 1 σ variation has asmall impact on the MINERvA distributions. The non-FSI fit parameters’ (e.g. M resA ) central values and uncer-tainties all agreed for the five fits. Here we present theresults from the FrAbs and FrInel fits.The NUISANCE interface to MINUIT2 [50] was usedto perform the fits. At each iteration, the GENIE-ReWeight package was used on an event-by-event basisto update the MC predictions before the total χ was cal-culated. The uncertainties in the fitted parameters weredetermined using the HESSE routine in MINUIT2. Thebest fit results from the joint tuning are shown in Ta-ble V. Fig. 8 shows the ratios of the best fit prediction tothe data for all four kinematic variables of interest whenthe pion absorption FSI parameter (FrAbs) is floated inthe fit; Fig. 9 is the same, but when the inelastic scatter-ing FSI parameter (FrInel) is floated. Notably, the twoFSI fits are very similar in both minimum χ and best-fitparameter values.Comparing to the results of the ANL and BNL reanaly-sis, larger values of M resA and smaller values of NormReswere found by the fit, pulling the parameters closer toGENIE nominal. The NonRes1 π parameter is stronglybound by the bubble chamber data and the MINERvAdata did little to improve on this constraint. The penaltyterm contributed to the χ by 9.3 for the FrAbs fit and11.1 for the FrInel fit. This is a significant improvementover the default, but indicates that there is mild tensionbetween the nucleon and nuclear data. The post-fit cor-relation matrices are provided in Fig 10. The ANL/BNLinput correlations are largely maintained in our fit.Tables VI and VII show the results when individualMINERvA data were tuned in separate fits. Since threeof the four channels were removed in these fits the con-straint from data is weakened and the total χ is steeredby the bubble chamber χ penalty. The individual chan-nel fits also found values at the 300% limit for NonRes2 π dial, except in the ν µ CC N π ± channel, where the resultwas unchanged by the fit. Only the ν µ CC N π ± chan-nel has a significant contribution from nonresonant 2 π ± production. In the other fits, the parameter is largelyunconstrained and has little impact on the fitted distri-butions. The χ per degree of freedom is indicative ofa poor fit in the ν µ CC N π ± and ν µ CC1 π channels, butnot in the ν µ CC1 π ± or ¯ ν µ CC1 π channels. Furthermore,the ν µ CC1 π shows the strongest χ penalty, indicatingtension with the ANL/BNL prior. Given the differentkinematic regions covered by the channels (see Table I)and the different physics ( e.g. fraction of coherent pionproduction) it is difficult to infer what combination of2 Parameter Default Value ANL/BNL Value FrAbs Fit Result FrInel Result M resA (GeV) 1 . ± .
22 0 . ± .
05 1 . ± .
04 1 . ± . ±
20 115 ± ± ± π (%) 100 ±
50 43 ± ± ± π (%) 100 ±
50 - 166 ±
32 161 ± π -iso 0 = RS - 1 = Iso (limit) 1 = Iso (limit)FrAbs (%) 100 ±
30 - 109 ±
16 -FrInel (%) 100 ±
40 - - 109 ± χ χ χ DoF
148 148 145 145TABLE V. Fit results from tuning GENIE parameters in NUISANCE. The “ANL/BNL Value” column shows the contributionswhen parameters are fixed at values of Table IV. effects are at work. Isotropic emission was preferred inall fits, driven by the θ π distributions. Disagreements inthe θ π spectrum are clearly seen in the data/MC ratiosof Fig. 8 and 9, and the large χ values observed for the ν µ CC N π ± and ν µ CC1 π channels.The individual χ contributions in the joint tuning bestfit, shown in sixth and seventh columns (“FrAbs Tune”and “FrInel Tune”) of Table II, show that not all distri-butions in all channels benefit from the model variations,as the default GENIE fits have a better χ for some dis-tributions. In particular, the ν µ CC1 π channel distribu-tions have worse agreement after the tuning, with onlythe θ π distribution improving in χ , whereas all chan-nels benefit from the shift to isotropic emission. Whilethere is an overall improvement over the ANL/BNL tunewhen comparing the combined χ results, Figs. 8 and 9show that there are still unresolved shape disagreementsin both the T π and θ µ kinematics.The tension between MINERvA’s nuclear data and theconstraints from ANL and BNL nucleon data is diffi-cult to confidently pinpoint; the lack of lepton mass ef-fects [51], modification to the resonance propagator inthe nucleus [52, 53], missing diagrams describing thenon-resonant background contributions [9], dynamicalcoupled channels [54], interactions on correlated initialstates, and the pion FSI model [1] are all part of an in-complete list of possible culprits. VI. AD HOC Q SUPRESSION
Further modifications beyond the standard GENIE di-als are required to resolve the observed tensions. Fig. 11(not used for any tuning) shows the Q distributions ob-served at MINERvA in for our tunes. The data is belowthe predictions of the tunes of Section V at low values of Q . There are also also differences at low θ µ , as shown inFig. 8 and 9. Measurements of ν µ CC1 π ± and ν µ CC1 π interactions on mineral oil at MiniBooNE have shown a data/MC shape discrepancy for the RS implementationin the NUANCE model [55, 56] in both Q and cos θ µ distributions [16, 57]. In the MINOS quasielastic analy-sis [15] on iron, which used NeuGen [58], a similar dis-agreement was observed when studying pion productiondominated sidebands. Indeed, concerns about low- Q modeling date back almost a decade [31]. The data fromMINOS and MiniBooNE experiments and the MINERvAdata on CH studied herein suggest that the RS imple-mentation common to each of the generators needs to besuppressed at low Q . Collective effects, which are usu-ally modeled in the random phase approximation, areknown to affect the Q distribution of neutrino-nucleusreactions at low Q . Motivated by these considerations,we attempted to improve the θ µ modeling by introducinga Q -dependent correction to the model.The MINOS collaboration suppression was expressedas R = A { − (cid:112) Q /Q } , (5)where the free parameters A = 1 .
010 and Q =0 .
156 GeV were empirically extracted from bin-by-bin fitsin Q to the data, and a hard cut-off at Q < . was imposed.We chose an empirical function so that the shape ofthe suppression preferred by each of the MINERvA chan-nels could be extracted. The empirical correction func-tion is applied to events with a resonance decay insidethe nucleus giving rise to a pion. Our suppression termis defined by choosing 3 points ( x i , R i ) i =1 , , between0 . < x < . . < R < .
0, where x ≡ Q . Moti-vated the ANL/BNL curves in Fig. 11, the correction isassumed to approach unity as Q approaches 0 . ,providing the constraint ( x , R ) = (0 . , . R by assuming a simple interpolation between the points3 Parameter ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ¯ ν µ CC1 π M resA (GeV) 0.97 ± . ± .
05 1 . ± .
05 0 . ± . ± ± ± ± π (%) 43 ± ± ± ± π (%) 300 (limit) 99 ±
30 300 (limit) 300 (limit) π -iso 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit)FrAbs (%) 156 ±
53 128 ±
34 126 ±
17 82 ± χ χ χ DoF
35 36 32 33TABLE VI. Individual channel tuning results when the FrAbs dial is treated as the free FSI parameter.Parameter ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ¯ ν µ CC1 π M resA (GeV) 0.97 ± . ± .
05 1 . ± .
05 0 . ± . ± ± ± ± π (%) 42 ± ± ± ± π (%) 300 (limit) 110 ±
30 300 (limit) 300 (limit) π -iso 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit)FrInel (%) 117 ±
54 127 ±
33 0 (limit) 80 ± χ χ χ DoF
35 36 32 33TABLE VII. Individual channel tuning results when the FrInel dial is treated as the free FSI parameter. ( x , . x , R ), and (0 . , . R ( Q < x ) = R ( Q − x )( Q − x )( x − x )( x − x )+ ( Q − x )( Q − x )( x − x )( x − x ) . (6)This interpolation function is then used to calculate thecorrection for each event as w ( Q ) = 1 − (1 − R )(1 − R ( Q )) . (7)where R defines the magnitude of the correction func-tion at the intercept, x = 0 . x is chosen to be Q = 0 .
35 GeV so that R describes the curvature atthe centre point of the correction. Expressing the weightswith Equations 6 and 7 ensures that the magnitude at x always lies between R and 1.0, avoiding parameter setswith large unphysical peaks in the correction function.Additionally, the squared term in Equation 7 ensuresthat w ( Q ) → . x → x , avoiding discontinuoussteps in the weighting function at x . The fitted pa-rameters R and R were limited to 0 . < R < . . < R < . e.g. doublepeaks. The fit results are shown in Table VIII. The correctionfrom the fit with FrAbs taken as a free parameter arecompared to the MINOS low- Q correction in Fig. 12.Our fits obtain a suppression factor that is similar to theMINOS one, with almost identical suppression at Q = 0,albeit with less curvature, particularly in the ν µ CC1 π ± and ν µ CC N π ± channels. The correction factors from thefit with FrInel or FrAbs as free parameters give similarresults.The correlation matrices for the fits including a Q de-pendent suppression are provided in Fig 13. Again, theANL/BNL input prior covariance is maintained. The pa-rameters largely correlate in the same way for the FrAbsand FrInel fit, and for the FrInel fit the R and R pa-rameters are negatively correlated.Figure 14 (Fig. 15) shows the ratio of the resulting fitsto the MINERvA data when FrAbs (FrInel) is taken asa free parameter. As anticipated, the predictions nowhave better agreement with the data in regards to the θ µ distribution, and the χ values are improved by theintroduction of our ad hoc low- Q correction. Other fitparameters are for the most part unchanged by the in-troduction of the low- Q correction. Furthermore, M resA and NormRes are closer to their values when fitting ANLand BNL data, indicating the Q correction alleviates the4 (GeV) m p M C / D a t a R a ti o ( R A TE ) + p CC1 m n + p CCN m n p CC1 m n p CC1 m n (deg.) m q M C / D a t a R a ti o ( S HA P E ) (MeV) p T M C / D a t a R a ti o ( S HA P E ) (deg.) p q M C / D a t a R a ti o ( S HA P E ) FIG. 8. MC/data ratios at the best fit point for the FrAbs joint tuning. tension between nucleon and nuclear modeling. Fig. 11shows the comparison of all our models directly againstMINERvA data in Q . Although the tuning sees im-provement in the χ for the ν µ CC1 π and ¯ ν µ CC1 π dis-tributions, the ν µ CC1 π ± and ν µ CC N π ± distributionsget worse, hinting at tensions in the charged and neu-tral pion production channels.Tables IX and X show the results of the fits to in-dividual channels, and Table XI shows the breakdownof contributions to the χ from the individual channels.The best fit χ value was significantly improved for eachchannel tuning when using a low- Q suppression and theextracted parameters were consistent with the ANL/BNLtunings. Pion kinematic distributions are not improved,and in some cases are slightly worse, as a result of includ-ing the low- Q suppression. It is clear from Table VIII(or by comparing Tables VI and IX) that the low- Q suppression has a similar effect in the fit to the FrAbsparameter. When the low- Q suppression is introduced,FrAbs tends to consistently lower values. It is also clearthat the ¯ ν µ CC1 π channel favors stronger low- Q sup-pression than the other channels. VII. CONCLUDING REMARKS
We have adjusted the parameters of the GENIE modelthat are important for pion production to match MIN-ERvA data in the ν µ CC1 π ± , ν µ CC N π ± , ν µ CC1 π and¯ ν µ CC1 π channels, using the NUISANCE framework.We incorporate existing results which informs the GE-NIE model using ANL and BNL bubble chamber datafrom scattering off protons and deuterons. Fits of se-lected GENIE model parameters were done using thekinematic distributions p µ , θ µ , T π and θ π . Parameterfits were performed with either the fraction of pions ab-sorbed or the fraction of pions inelastically scattered inFSI as a floating parameter, with broadly similar conclu-sions for the two cases.The results of the fit (see Table V) show that the tuningimproves the GENIE pion production model significantly,but tensions remain. The pull on the ANL/BNL priordemonstrates a tension between MINERvA nuclear tar-get data and the light-target bubble-chamber data setsused to make the prior, indicating a deficiency in theGENIE nuclear model which cannot be fixed by modify-ing the available reweighting dials. Additionally, fittingto individual MINERvA pion production channels pro-5 (GeV) m p M C / D a t a R a ti o ( R A TE ) + p CC1 m n + p CCN m n p CC1 m n p CC1 m n (deg.) m q M C / D a t a R a ti o ( S HA P E ) (MeV) p T M C / D a t a R a ti o ( S HA P E ) (deg.) p q M C / D a t a R a ti o ( S HA P E ) FIG. 9. MC/data ratios at the best fit point for the FrInel joint tuning. C o rr e l a t i on - - - - - r e s A M N o r m . R e s . p N on R e s . - i s o p F r A b s p N on R e s . resA MNorm. Res. p Non Res. 1-iso p FrAbs p Non Res. 2 C o rr e l a t i on - - - - - r e s A M N o r m . R e s . p N on R e s . - i s o p F r I ne l p N on R e s . resA MNorm. Res. p Non Res. 1-iso p FrInel p Non Res. 2
FIG. 10. Correlation matrix from from tuning GENIE parameters in NUISANCE with FrAbs included as a fit parameter (left)and with FrInel included as a fit parameter (right). ) (GeV Q / nu c l eon ) / G e V ( c m / d Q s d - · – p CC1 m n MINERvA =10.6 c Default, =13.1 c ANL/BNL, =5.5 c FrAbs, =5.8 c FrInel, =17.4 c , FrAbs+Q =13.8 c , FrInel+Q ) (GeV Q / nu c l eon ) / G e V ( c m / d Q s d - · – p CCN m n MINERvA =14.8 c Default, =16.7 c ANL/BNL, =9.9 c FrAbs, =10.5 c FrInel, =30.5 c , FrAbs+Q =28.2 c , FrInel+Q ) (GeV Q / nu c l eon ) / G e V ( c m / d Q s d - · p CC1 m n MINERvA =35.1 c Default, =62.0 c ANL/BNL, =46.4 c FrAbs, =47.2 c FrInel, =33.7 c , FrAbs+Q =34.9 c , FrInel+Q ) (GeV Q / nu c l eon ) / G e V ( c m / d Q s d - · p CC1 m n MINERvA =10.5 c Default, =20.8 c ANL/BNL, =11.7 c FrAbs, =11.8 c FrInel, =6.5 c , FrAbs+Q =6.4 c , FrInel+Q
FIG. 11. Comparisons of the nominal and tuned models to MINERvA ν µ CC1 π ± (left top), ν µ CC Nπ ± (right top), ν µ CC1 π (left bottom) and ¯ ν µ CC1 π (right bottom) distributions in Q . The χ is computed using the full covariance matrices. Thedistributions were not explicitly used in the tuning procedure. Parameter FrAbs Tune FrAbs + low- Q Tune FrInel Tune FrInel + low- Q Tune M resA ( GeV ) 1 . ± .
04 0 . ± .
02 1 . ± .
04 0 . ± . ± ± ± ± π (%) 43 ± ± ± ± π (%) 166 ±
32 99 ±
31 161 ±
33 120 ± π -iso 1 . ±
16 48 ±
21 - -FrInel (%) - - 109 ±
24 132 ± R - 0 . ± .
06 - 0 . ± . R - 0 . . ± . χ χ χ DoF
145 143 145 143TABLE VIII. Ad hoc low- Q suppression model tuning results compared to the tuning results without the low- Q suppression. ) (GeV Q C o rr ec ti on F ac t o r Joint FrAbs Fit FrAbs Fit + p CC1 m n FrAbs Fit + p CCN m n MINOS Parametrisation ) (GeV Q C o rr ec ti on F ac t o r Joint FrAbs Fit FrAbs Fit p CC1 m n FrAbs Fit p CC1 m n MINOS Parametrisation
FIG. 12. Extracted low- Q suppression factors from the FrAbs + low- Q tuning to each channel. The left and right plotscompare the results for the charged and neutral pion production channels respectively. Shown in red is the uncertainty bandextracted from the joint fit to all 4 channels simultaneously. duces different best-fit parameters, demonstrating thatGENIE cannot describe the different exclusive channelsin a consistent manner with the available dials (shownin Tables VI and VII). Because the four channels coverdifferent kinematic regions (see Table I) and contain dif-ferent physics ( e.g. different coherent pion productioncontributions or nonresonant processes), it is difficult topinpoint the origin of the discrepancy between the modeland the different MINERvA data sets.Following experimental hints of discrepancies at low- Q for a variety of cross-section measurements on nucleartargets, an additional empirical low- Q suppression wasintroduced and the fits were repeated. Although the datashowed a preference for a strong suppression at low- Q and the agreement improved for θ µ and Q distributions,tensions remain. In particular, fits to individual MIN-ERvA channels still produced different results, and favordifferent parameter values for the low- Q suppression.The main conclusion of this work is that currentneutrino experiments operating in the few–GeV regionshould think critically about single pion production mod-els and uncertainties, as the Monte Carlo models whichare currently widely used in the field are unable to ex-plain multiple data sets, even when they are from a singleexperiment.A key strength of this analysis is its developmentwithin the NUISANCE framework, allowing it to be eas-ily repeated with alternate model assumptions, neutrino8 C o rr e l a t i on - - - - - r e s A M N o r m . R e s . p N on R e s . - i s o p F r A b s p N on R e s . Lag R Lag R resA MNorm. Res. p Non Res. 1-iso p FrAbs p Non Res. 2 Lag R Lag R C o rr e l a t i on - - - - - r e s A M N o r m . R e s . p N on R e s . - i s o p F r I ne l p N on R e s . Lag R Lag R resA MNorm. Res. p Non Res. 1-iso p FrInel p Non Res. 2 Lag R Lag R
FIG. 13. Correlation matrix from tuning GENIE parameters with an ad hoc low- Q supression with FrAbs included as a fitparameter (left) and with FrInel included as a fit parameter (right). (GeV) m p M C / D a t a R a ti o ( R A TE ) + p CC1 m n + p CCN m n p CC1 m n p CC1 m n (deg.) m q M C / D a t a R a ti o ( S HA P E ) (MeV) p T M C / D a t a R a ti o ( S HA P E ) (deg.) p q M C / D a t a R a ti o ( S HA P E ) FIG. 14. MC/data ratios at the best fit points from the FrAbs tuning with low- Q suppression included. (GeV) m p M C / D a t a R a ti o ( R A TE ) + p CC1 m n + p CCN m n p CC1 m n p CC1 m n (deg.) m q M C / D a t a R a ti o ( S HA P E ) (MeV) p T M C / D a t a R a ti o ( S HA P E ) (deg.) p q M C / D a t a R a ti o ( S HA P E ) FIG. 15. MC/data ratios at the best fit points from the FrInel tuning with low- Q suppression included.Parameter ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ¯ ν µ CC1 π M resA (GeV) 0 . ± .
02 0 . ± .
02 0 . ± .
05 0 . ± . ± ± ± ± π (%) 43 ± ± ± ± π (%) 300 (limit) 70 ±
28 300 (limit) 300 (limit) π -iso 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit)FrAbs (%) 92 ±
65 79 ±
40 74 ±
22 34 ± R . ± .
16 0 . ± .
13 0 . ± .
14 0 . ± . R .
50 (limit) 0 .
50 (limit) 0 . ± .
31 1 .
00 (limit)MINERvA χ χ χ DoF
33 34 30 31TABLE IX. Individual channel FrAbs + low- Q tuning results. Parameter ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ¯ ν µ CC1 π M resA (GeV) 0 . ± .
02 0 . ± .
02 0 . ± .
05 0 . ± . ± ± ± ± π (%) 43 ± ± ± ± π (%) 300 (limit) 78 ±
28 300 (limit) 300 (limit) π -iso 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit) 1 = Iso (limit)FrInel (%) 179 ±
63 173 ±
37 8 ±
125 103 ± R . ± .
14 0 . ± .
13 0 . ± .
17 0 . ± . R .
50 (limit) 0 .
50 (limit) 0 . ± .
37 1 .
00 (limit)MINERvA χ χ χ DoF
33 34 30 31TABLE X. Individual channel FrInel + low- Q tuning results.Distribution Channel N bins FrAbs Tune FrAbs + low- Q Tune FrInel Tune FrInel + low- Q Tune p µ (Rate) ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ν µ CC1 π θ µ (Shape) ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ν µ CC1 π T π (Shape) ν µ CC1 π ± ν µ CC Nπ ± ν µ CC1 π ν µ CC1 π θ π (Shape) ν µ CC1 π ±
14 13.0 13.4 12.6 12.6 ν µ CC Nπ ±
14 6.9 7.0 6.2 6.3 ν µ CC1 π
11 8.3 12.2 8.9 9.4¯ ν µ CC1 π
11 3.4 4.4 3.5 3.7Total χ
148 242.3 212.2 240.7 215.7TABLE XI. Channel by channel contributions to the χ at for the GENIE tunings with and without the low- Q correctionincluded. interaction generators, and different data. The develop-ments presented here will be used in future iterations ofthis work, as the MINERvA collaboration works towardsa GENIE model that provides a good description of alltheir available data, and can be easily applied to othermeasurements and experiments. ACKNOWLEDGEMENTS
P.S., L.P., and C.V.C.W. would like to thank the UKScience and Technology Facilities Council (STFC) forPh.D. funding support. P.S. acknowledges the FermilabNeutrino Physics Center for the scholarship that funded this work and thanks Fermilab and the MINERvA col-laboration for their hospitality during this work. C.W.acknowledges the support of the Swiss National ScienceFoundation and SERI.This document was prepared by members of the MIN-ERvA Collaboration using the resources of the FermiNational Accelerator Laboratory (Fermilab), a U.S. De-partment of Energy, Office of Science, HEP User Facil-ity. Fermilab is managed by Fermi Research Alliance,LLC (FRA), acting under Contract No. DE-AC02-07CH11359. These resources included support for theMINERvA construction project, and support for con-struction also was granted by the United States NationalScience Foundation under Award No. PHY-06197271and by the University of Rochester. Support for par-ticipating scientists was provided by NSF and DOE(USA); by CAPES and CNPq (Brazil); by CoNaCyT(Mexico); by Proyecto Basal FB 0821, CONICYT PIAACT1413, Fondecyt 3170845 and 11130133 (Chile); byCONCYTEC, DGI-PUCP, and IDI/IGI-UNI (Peru); and by the Latin American Center for Physics (CLAF); NCNOpus Grant No. 2016/21/B/ST2/01092 (Poland). Wethank the MINOS Collaboration for use of its near de-tector data. Finally, we thank the staff of Fermilab forsupport of the beam line, the detector, and computinginfrastructure. [1] U. Mosel and K. Gallmeister, Phys. Rev. C (2017).[2] L. Alvarez-Ruso et al. , Prog. Part. Nucl. Phys. (2018).[3] K. N. R. Gonz´alez-Jim´enez and N. Jachowicz, Phys. Rev.D (2018).[4] C. Wilkinson et al. , Phys. Rev. D (2016).[5] K. Mahn, C. Marshall, and C. Wilkinson, Ann. Rev.Nucl. Part. Sci. (2018).[6] P. Stowell et al. , J. Instrum. (2017).[7] C. Andreopoulos et al. , Nucl. Instrum. Methods Phys.Res. (2009).[8] C. Andreopoulos et al. , arXiv:1510.05494 (2015).[9] M. Kabirnezhad, Phys. Rev. D (2018).[10] T. Leitner, O. Buss, L. Alvarez-Ruso, and U. Mosel,Phys. Rev. C (2009).[11] J. Nieves, I. R. Simo, and M. J. V. Vacas, Phys. Rev. C (2011).[12] E. Hernadez and J. Nieves, Phys. Rev. D (2017).[13] P. Rodrigues, C. Wilkinson, and K. McFarland, Eur.Phys. J. C (2016).[14] C. Wilkinson, P. Rodrigues, S. Cartwright, L. Thompson,and K. McFarland, Phys. Rev. D (2014).[15] P. Adamson et al. (MINOS), Phys. Rev. D (2015).[16] A. A. Aguilar-Arevalo et al. (MiniBooNE), Phys. Rev. D (2011).[17] B. Eberly et al. (MINERvA), Phys. Rev. D (2015).[18] C. L. McGivern et al. (MINERvA), Phys. Rev. D (2016).[19] O. Altinok et al. (MINERvA), Phys. Rev. D (2017),with errata of summer 2019.[20] L. Aliaga et al. (MINERvA), Nucl. Instrum. MethodsPhys. Res., Sect. A (2014).[21] T. Katori, Ph. D. thesis, University of Indiana (2008).[22] R. A. Smith and E. J. Moniz, Nucl. Phys. B (1972).[23] A. Bodek and J. Ritchie, Phys. Rev. D (1981).[24] J. Nieves, J. E. Amaro, and M. Valverde, Phys. Rev. C (2004).[25] R. Gran, J. Nieves, F. Sanchez, and M. J. V. Vacas,Phys. Rev. D (2013).[26] C. E. Patrick et al. (MINERvA), Phys. Rev. D (2018).[27] A. Mislivec et al. (MINERvA), Phys. Rev. D (2018).[28] L. Aliaga et al. (MINERvA), Phys. Rev. D (2016).[29] C. Andreopoulos et al. (1981).[31] K. M. Graczyk and J. T. Sobczyk, Phys. Rev. D (2008).[32] A. Bodek and U. K. Yang, J. Phys. G (2003).[33] T. Yang, C. Andreopoulos, H. Gallagher, K. Hoffmann,and P. Kehayias, Eur. Phys. J. C (2009).[34] Z. Koba, H. B. Nielsen, and P. Olesen, Nucl. Phys. B (1972).[35] S. M. T. Sjostrand and P. Z. Skands, Comput. Phys.Commun. (2008).[36] D. Rein and L. M. Sehgal, Nucl. Phys. B (1983).[37] D. Rein and L. M. Sehgal, Phys. Lett. B (2007).[38] C. Berger and L. M. Sehgal, Phys. Rev. D (2009).[39] S. Dytman, Acta Phys. Polon. B (2009).[40] R. A. Arndt, W. J. Briscoe, I. I. Strakovsky, and R. L.Workman, Phys. Rev. C (2006).[41] O. Buss et al. , Phys. Rept. (2012).[42] G. M. Radecky and others (ANL), Phys. Rev. D (1982).[43] T. Kitagaki and others (BNL), Phys. Rev. D (1986).[44] S. J. Barish and others (ANL), Phys. Rev. D (1979).[45] T. Kitagaki and others (BNL), Phys. Rev. D (1990).[46] S. L. Adler, Annals of Physics (1968).[47] P. Allen et al. (Aachen-Bonn-CERN-Munchen-Oxford),Nuclear Physics B (1980).[48] D. Allasia et al. (Amsterdam-Bologna-Padova-Pisa-Saclay-Torino), Z. Phys. C (1983).[49] J. Bell. et al. , Phys. Rev. Lett. (1978).[50] F. James and M. Roos, Comput. Phys. Commun. (1975).[51] C. Berger and L. M. Sehgal, Phys. Rev. D (2007).[52] S. K. Singh et al. , Physics Letters B (1998).[53] E. Oset and L. L. Salcedo, Nuclear Physics A (1987).[54] S. X. Nakamura et al. , Reports on Progress in Physics (2017).[55] D. Casper, Nucl. Phys. Proc. Suppl. (2002).[56] P. Przewlocki, Acta Phys. Polon. B (2009).[57] A. A. Aquilar-Arevalo et al. (MiniBooNE), Phys. Rev. D (2011).[58] H. Gallagher, Nucl. Phys. Proc. Suppl.112