The Ghent Hybrid model in NuWro: a new neutrino single-pion production model in the GeV regime

JHEP12(2024)141 Published for SISSA by Springer Received: May 9, 2024 Revised: November 13, 2024 Accepted: November 18, 2024 Published: December 17, 2024 The Ghent Hybrid model in NUWRO: a new neutrino single-pion production model in the GeV regime Qiyu Yan ,a,b Kajetan Niewczas ,cAlexis Nikolakopoulos ,d,c Raúl González-Jiménez ,eNatalie Jachowicz ,cXianguo Lu ,bJan Sobczyk f and Yangheng Zheng a aSchool of Physical Sciences, University of Chinese Academy of Sciences, Beijing 100049, China bDepartment of Physics, University of Warwick, Coventry, CV4 7AL U.K. cDepartment of Physics and Astronomy, Ghent University, Proeftuinstraat 86, B-9000 Gent, Belgium d Theoretical Physics Department, Fermilab, Batavia IL, U.S.A. eDepartamento de Física Atómica, Molecular y Nuclear, Universidad de Sevilla, 41080 Sevilla, Spain fInstitute of Theoretical Physics, University of Wrocław, pl. M. Borna 9, 50-204, Wrocław, Poland E-mail: [email protected],[email protected], [email protected],[email protected],[email protected], [email protected],[email protected],[email protected] Abstract: Neutrino-induced single-pion production constitutes an essential interaction channel in modern neutrino oscillation experiments, with its products building up a significant fraction of the observable hadronic final states. Frameworks of oscillation analyses strongly rely on Monte Carlo neutrino event generators, which provide theoretical predictions of neutrino interactions on nuclear targets. Thus, it is crucial to integrate state-of-the-art single-pion production models with Monte Carlo simulations to prepare for the upcoming systematics-dominated landscape of neutrino measurements. In this work, we present the implementation of the Ghent Hybrid model for neutrino-induced single-pion production in the NuWro Monte Carlo event generator. The interaction dynamics includes coherently-added contributions from nucleon resonances and a non-resonant background, merged into the pythia branching predictions in the deep-inelastic regime, as instrumented by NuWro. This neutrino-nucleon interaction model is fully incorporated into the nuclear framework of the generator, allowing it to account for the influence of both initial- and final-state nuclear medium effects. We compare the predictions of this integrated implementation with recent pion production data from accelerator-based neutrino experiments. The results Open Access,©The Authors. Article funded by SCOAP3.https://doi.org/10.1007/JHEP12(2024)141 JHEP12(2024)141 of the novel model show improved agreement of the generator predictions with the data and point to the significance of the refined treatment of the description of pion-production processes beyond the ∆region. Keywords: Neutrino Interactions, Neutrino Mixing ArXiv ePrint: 2405.05212 JHEP12(2024)141 Contents 1 Introduction 1 2 The Ghent Hybrid model for single-pion production 3 3 NUWRO event generator 4 3.1 Description of the ∆(1232) resonance 4 3.2 ∆-PYTHIA SPP model 5 3.3 NUWRO algorithm 5 3.4 Hybrid model in NUWRO 6 3.5 Nuclear modeling 8 4 Comparison with data 9 4.1 Pion kinematics 9 4.2 Transverse kinematics imbalance 11 5 Conclusions 16 A Comparison with T2K π+TKI data 18 B Comparison with MINERvA TKI data using ESF 19 C Single pion production dσ dWwith 3 GeV neutrino 19 1 Introduction In the energy range of a few GeV, as explored by accelerator-based neutrino experiments such as DUNE [ 1 ], Hyper-Kamiokande [ 2 ], NOvA [ 3 ], and T2K [ 4 ], the understanding of the neutrinonucleus scattering cross section is limited to a precision of, at best, 10% [ 1 ], corresponding to a 3–7% contribution to the overall experimental systematic error [ 5 ]. Alongside the normalization of the neutrino flux, this limitation represents a significant source of uncertainty in oscillation analyses. Achieving the measurement of the CP violation phase, which can shed light on mechanisms explaining the prevalence of matter over antimatter in the current universe, requires reducing systematic errors to as low as 1% [ 1 ]. Neutrino-induced single-pion production (SPP) is an important reaction channel in these experiments. Its significance extends to atmospheric neutrino programs in Super-Kamiokande [ 6 ], Hyper-Kamiokande, and JUNO [ 7 ], where GeV neutrinos play an important role in determining the neutrino mass ordering. We will focus on discussing SPP in the NuWro Monte Carlo (MC) generator. NuWro is one of the major MC generators extensively used in studies by experimental groups, and the accuracy of its predictions is vital for the precision of neutrino oscillation measurements. Pion production in neutrino interactions occurs via two primary mechanisms: resonance production and non-resonance production. The current NuWro implementation utilizes a dedicated ∆resonance production and decay model for the resonant contribution. The final state – 1 – JHEP12(2024)141 Figure 1. Diagrammatic definitions of the ∆-pythia (∆-P) and Hybrid-pythia (H-P) models, with the former combining the ∆model and pythia in NuWro, while the latter replaces the ∆model with the Ghent SPP Hybrid model. All models are arranged by increasing W , as indicated by the long arrow. angular distribution is taken from the experimental results from ANL or BNL [ 8 ]. This model accurately describes the inclusive cross section for pion production through ∆resonance formation [ 8 ]. However, the reliability of the NuWro approach decreases in the higher invariant mass ( W ) region due to the absence of explicit contributions from other resonances. To address this limitation, NuWro incorporates pythia [ 9 ] to model the final state originating from other resonances within the framework of non-resonant production, using inclusive cross sections from the Bodek-Yang approach of quark parton distribution functions [ 10 ]. Additionally, pythia is employed in the ∆region to account for pion production through non-resonant processes. However, using Bodek-Yang approach and pythia hadronization models at low W and low squared-four momentum transfer ( Q2 ) raises concern. In this kinematic regime, interactions do not primarily occur on quarks, and the underlying assumptions of the Lund model, the theoretical foundation of pythia, are barely met [ 11 ]. To address these drawbacks, in ref. [ 12 ], different approaches to implement in MC generators the model of refs. [ 13 , 14 ] for SPP off the nucleon were discussed. The idea is to model neutrino-induced SPP across the entire kinematic range, extending up to potentially unlimited W values. Consequently, there is no reliance on the hadronization model for the SPP channel, enabling modeling of the differential cross section across the complete kinematic spectrum. In this study, we incorporate this methodology into the NuWro event generator framework, ensuring its seamless integration with the nuclear model and the intranuclear cascade. To evaluate its efficacy, this new model implementation is compared against existing datasets. Specifically, we will utilize the MINERvA π+ sample [ 15 ] and MicroBooNE π0 data [ 16 ], which are mostly sensitive to the ∆resonance region, and the MINERvA π0 transverse kinematic imbalance (TKI) [ 17 ] data, which has significant contributions from heavier resonances. The paper is structured as follows. In section 2we summarize the main features of the Ghent Hybrid model. Section 3details the default description of SPP processes in NuWro and the implementation of the transition to the regime where multi-pion and – 2 – JHEP12(2024)141 deep inelastic scattering (DIS) events become important, both in standard NuWro and in the new implementation incorporating the Hybrid model (figure 1). The resulting theory predictions are compared to the experimental measurements in section 4, followed by the conclusions in section 5. 2 The Ghent Hybrid model for single-pion production The Hybrid model [ 13 ] for electroweak SPP developed by the Ghent group consists of two main components, each targeting distinct kinematic regions. The low-energy part of the model (LEM) includes contributions from resonances and background based on treelevel diagrams described in ref. [ 18 ]. It has been systematically benchmarked against data and the MAID07 [ 19 ] and dynamical coupled channels models [ 20 , 21 ] with satisfactory agreement [ 13 , 14 , 22 ]. For the ∆resonance, the vector form factors of ref. [ 23 ] are used. The axial form factors were determined through analysis of bubble chamber data detailed in ref. [ 24 ]. The vector form factors for the higher mass resonances ( P11 (1440) , D13 (1520) , S11 (1535)) are described in ref. [ 14 ], drawing upon findings from refs. [ 19 , 23 , 25 ]. The axial form factors for higher-mass resonances are reported in ref. [ 13 ]. At Q2 = 0, they are determined based on partially conserved axial current (PCAC) considerations, disregarding any undetermined axial couplings. The Q2 -dependence of the axial form factors is taken from ref. [ 23 ]. At high energies ( W above 1 . 5-2GeV), the LEM model starts exhibiting anomalous behavior because it relies solely on tree-level diagrams. To overcome this limitation, the Hybrid model employs a description of the non-resonant background based on Regge phenomenology (ReChi). In this approach, the tree-level propagator of t -channel meson exchanges in the lowenergy background is replaced by a Regge propagator [ 26 ]. This approach offers an efficient description of the forward-scattering process, which largely dominates in the high-energy regime [ 26 – 28 ]. Additionally, the high-energy behavior of the resonances is regulated by the inclusion of cut-off hadronic form factors [ 29 ]. The transition between the low- and high-energy background is implemented, at the amplitude level, in a phenomenological way as a function of W : Jµ Hybrid =Jµ RES + cos2ϕ(W)Jµ LEM + sin2ϕ(W)Jµ ReChi,(2.1) where Jµ Hybrid is the full hadronic current, while Jµ RES , Jµ LEM , and Jµ ReChi denote the contributions from resonances, low-energy background, and high-energy background, respectively. The transition function ϕ ( W )depends on W and is defined as follows: ϕ(W) = π 2  1−1 1 + exp W−W0 L  .(2.2) Here, W0 and L represent the center and width of the transition, respectively, fixed at 1 . 5 GeV and 0 . 1 GeV . This means that below 1 . 4 GeV , the prescription is essentially provided by LEM, whereas above 1 . 6 GeV , the strength stems predominantly from the ReChi model. The Hybrid model has been utilized in several studies to predict lepton-induced SPP on the nucleon and nucleus. The latter studies employed the relativistic plane-wave impulse approximation (RPWIA) and the relativistic distorted-wave impulse approximation (RDWIA) – 3 – JHEP12(2024)141 as the nuclear framework [ 14 , 30 – 33 ]. In both approaches, the initial (bound) nucleons are described as relativistic mean-field wave functions, i.e., solutions of the Dirac equation with relativistic potentials. In RPWIA, the final (knocked-out) nucleon is treated as a plane wave, while in RDWIA, it is represented by a distorted wave, accounting for elastic final-state interactions (FSIs). In both cases, the pion was treated as a plane wave. The implementation of the Hybrid Model in NuWro for neutrino-nucleon pion production was presented in ref. [ 12 ]. In this work, we incorporate the model into NuWro’s nuclear framework, which employs a local Fermi gas (or the factorized spectral function approach) along with an intranuclear cascade model for the propagation of hadrons. 3 NUWRO event generator The NuWro MC event generator covers neutrino energies ranging from approximately 100 MeV to 100 GeV [ 34 ]. Neutrino-nucleus interactions are typically modeled within the impulse approximation scheme, wherein the initial interaction occurs on a bound nucleon followed by re-interactions of resulting hadrons inside the nucleus. Inelastic initial interactions are categorized as either resonant production (RES) for invariant hadronic masses W≤ 1 . 6 GeV or deep inelastic scattering (DIS) for W > 1 . 6 GeV [ 35 ]. NuWro is optimized for sub-GeV neutrino beams, with RES dominated by the ∆(1232) excitation and its subsequent decay leading to SPP: νN′→∆ℓ, ∆→πN,(3.1) where N (′) and ℓ represent the final (initial) nucleon and the scattered lepton, respectively. Motivated by quark-hadron duality, NuWro has implemented an explicit model for ∆ excitation, while the strength provided by heavier resonances is included in an effective manner [ 36 ]. Additional inelastic channels, such as two-pion production, are modelled with pythia. By “pythia”, we mean that the lepton inclusive cross section is evaluated using the Bodek-Yang approach [ 10 ], and the final hadronic state is obtained using the hadronization routines of pythia. In the case of SPP, a transition region from RES to DIS is modeled to ensure a smooth transition between both regimes. To achieve a quantitative comparison between the predictions of the SPP model and experimental data, accurate models for both the initial state and FSIs must be incorporated. The effects of FSIs alter the observed spectra compared to the original model predictions. For instance, processes like pion absorption or production can induce changes in the observed final-state topology. 3.1 Description of the ∆(1232) resonance The NuWro ∆excitation model is formulated in terms of form factors obtained through a simultaneous fit to both ANL and BNL SPP data [ 8 ]. This fitting procedure was conducted for the νµp→µ−pπ+ reaction, assuming negligible contribution from the non-resonance background. Consequently, in the NuWro neutrino νµp→µ−pπ+ channel, there is no non-resonant background. However, in the νµn→µ−nπ+ and νµn→µ−pπ0 channels (and in analogous channels for ¯νµ scattering, related by isospin symmetry), a non-resonant – 4 – JHEP12(2024)141 Channel α0 νp →ℓ−pπ+0.0 ¯νn →ℓ+nπ− νn →ℓ−nπ+0.2 ¯νp →ℓ+pπ− νn →ℓ−pπ00.3 ¯νp →ℓ+nπ0 Table 1. NuWro default values of the α0parameter (see eq. (3.4)). background is added incoherently as a fraction of the pythia contribution, its size guided by the experimental data. In the paragraphs below, we describe the merging of both models to account for inelastic events in NuWro. 3.2 ∆-PYTHIA SPP model When modeling SPP channels, NuWro implements the following cross section formula [ 36 ]: σSPP =β(W)σ∆+α(W)σpythia:SPP,(3.2) where σ∆ stands for the ∆contribution and σpythia:SPP represents the SPP component of the DIS cross section as modeled by pythia. σpythia:SPP is defined with single pion production functions fSPP ( W )extracted from pythia as probabilities to get SPP final state as an outcome of the hadronization: σpythia:SPP =fSPP(W)·σDIS,(3.3) with σDIS obtained from the Bodek-Yang approach. The blending of the two is controlled by the functions α ( W )and β ( W ). In NuWro, the default option is that β ( W ) = 1 −α ( W ) (more general assumption is available as well), where α ( W )is a linear function of W : α(W) =          W−Wth Wmin−Wth α0W < Wmin, W−Wmin+α0(Wmax−W) Wmax−Wmin Wmin ≤W≤Wmax, 1Wmax < W. (3.4) Wth = M + mπ is a threshold for pion production (with M the nucleon mass), Wmin and Wmax define the transition region, and their default values are Wmin = 1 . 3GeV and Wmax = 1 . 6GeV. The linear function α ( W )has the properties: α ( Wth ) = 0, α ( Wmin ) = α0 , and α ( Wmax )=1. The parameter α0 defines the size of the non-resonant background. Its value is separately selected for each SPP channel, as shown in table 1. 3.3 NUWRO algorithm In the RES region W < Wmax NuWro combines SPP model with other inelastic channels, such as two-pion production, extracted from pythia. The generation of inelastic events is achieved through the following steps: – 5 – JHEP12(2024)141 1. Random phase space sampling: the algorithm begins by randomly selecting a point within the two-dimensional available kinematic phase space ( W, Q2 )for a given neutrino energy, with the constraint W < Wmax . This sampled point defines the W value for the resulting hadronic system. 2. Hadronization and event weight: at the selected phase space point, pythia is invoked for hadronization, where the hadronic energy and momentum are converted into final-state particles. (a) Non-SPP Event: if the final state is not SPP, the event is accepted and assigned with a weight equal to d2σDIS dQ2dW multiplied by the available phase space in ( W, Q2 ). (b) SPP Event: if the final state generated by pythia is SPP, the event is assigned with a weight 1 fSPP(W) d2σSPP dQ2dWmultiplied by the available phase space. The blending defined in eq. (3.2) is realized by probabilistically choosing an either ∆or pythia origin of the event according to β ( W ) : α ( W )ratio. If the former is chosen, the kinematics of the final hadronic state is replaced by the outcome of ∆ decay, incorporating information regarding angular correlations. If the latter is selected, the event is stored as it is. 3.4 Hybrid model in NUWRO Considering the Hybrid model’s ability to describe SPP across a wide kinematic range, an expansion of NuWro’s current transition region to higher values of W appears reasonable. A motivation to do this follows from what is seen in figure 2, where we compare the predicted SPP cross sections ( dσ dW ) using different models (another set of predictions with lower neutrino energy 3 GeV is shown in figure 16 in appendix C). Around W = 1 . 5 GeV , the Hybrid model exhibits a second peak, particularly noticeable in the p π0 and n π+ channels (figures 2b and 2c), attributed to contributions from the second family of nucleon resonances, D13 (1520), and S11 (1535). The ∆model does not incorporate these contributions. In the results of the ∆-P model, a structure can be observed that might mistakenly appear to arise from the second resonance region. However, this structure is merely an artifact of the transition region and holds no physical significance, as discussed in ref. [ 36 ]. Moreover, the ∆model requires an earlier transition to the pythia model to include the cross section strength in the second resonance region and beyond, which is not accounted for in the pure ∆model. To minimize or eliminate Pythia contributions to the SPP channel, the transition region would ideally be placed at very high W values. However, in the current NuWro algorithm governing single pion production and DIS events, such a high W transition makes the code significantly more time consuming. A remedy would be to replace the current RES/DIS transition defined according to values of W , by SPP/DIS transition defined by topology of final states, but this requires a change of the basic structure of the code. In this study of the Hybrid model, we set Wmin = 2 . 8GeV and Wmax = 3 . 2GeV. Additionally, as the Hybrid model already contains contributions from the non-resonant background, in its implementation in NuWro, we set α0 = 0 (see eq. (3.4)). Henceforth, we will refer to the NuWro ∆and Hybrid models incorporating the transition to pythia as ∆-pythia and Hybrid-pythia, respectively. – 6 – JHEP12(2024)141 1 2 3 4 (GeV)W 0 20 40 60 39− 10× /GeV) 2 (cmW/d σ d + π p − µ → p µ ν(a) 8 GeV, H-P Hybrid -P∆ ∆ PYTHIA 1 2 3 4 (GeV)W 0 5 10 15 39− 10× /GeV) 2 (cmW/d σ d 0 π p − µ → n µ ν(b) 8 GeV, H-P Hybrid -P∆ ∆ PYTHIA 1 2 3 4 (GeV)W 0 5 10 15 39− 10× /GeV) 2 (cmW/d σ d + π n − µ → n µ ν(c) 8 GeV, H-P Hybrid -P∆ ∆ PYTHIA Figure 2. Predictions for SPP by an 8 GeV neutrino off the nucleon for the models discussed in this paper for three possible channels: (a) νµp→µ−pπ+ , (b) νµn→µ−pπ0 , and (c) νµn→µ−nπ+ . The transition regions in W are 1.3–1.6 GeV and 2.8–3.2 GeV for the ∆-P and H-P models, respectively. The wiggle in ∆-P and H-P curves is due to numerical fluctuations in evaluating fSPP ( W ). Note the cross sections at the ∆peaks predicted by the ∆model: 9:2:1, as expected from isospin arguments [37–39]. – 7 – JHEP12(2024)141 •For CC-0π[54]: 1. νµ + A→µ− + p + X , requiring one µ and at least one proton in the final state, with no pions. 2. 1.5 GeV/c< pµ<10 GeV/c, θµ<20◦. 3. 0.45 GeV/c< pp<1.2 GeV/c, θp<70◦. For the 0 π measurement, the contribution from pion production is relatively minor: arising from pion absorption, a process addressed by the FSI mechanism included in the NuWro cascade model. Therefore, the difference between the predictions of the ∆and Hybrid model can be anticipated to be minimal. On the other hand, in the π0 measurement, there is no restriction based on Wrest . The predicted W distributions are shown in figure 7b. The presence of a resonance peak near W = 1 . 5 GeV and the ability to extend to higher W values in the Hybrid model significantly alter the distributions. The contribution from the ∆resonance in the ∆-P model ceases at W = 1 . 6 GeV . As explained before (section 3.4), an apparent presence of a contribution from the second resonance region in the ∆-P model results is an artifact arising from the RES-DIS transition region. In contrast, the Hybrid model exhibits a pronounced peak in the second resonance region, which is not present in the ∆model. As a consequence, pythia is extrapolated down to a low W -value in ∆-P, resulting in a less realistic picture, with significantly more strength in the second resonance region compared to the H-P model. Figure 8illustrates the comparison between the predictions of the ∆and Hybrid models and the MINERvA 0 π and π0 measurements. As mentioned before, in the 0 π measurement, pion production contributes to the event sample via pion absorption during FSI (upper panels), therefore, the change of pion production model has only limited impact on the δαT and pN distributions. Conversely, when we compare to the π0 measurement (lower panels in figure 8), the model predictions become noticeably different. Firstly, a change in the normalisation is observed in δαT . Since the shape of δαT is dictated by FSI, which is independent of the pion production modeling, we do not expect to see shape differences for the ∆and Hybrid models predictions. Note that, as was pointed out in the original paper [ 17 ], it is interesting that the δαT shape from the two measurements happens to be similar (figure 9a). Now, this similarity is also captured by the models. Secondly, in the π0 measurement depicted in figure 8, the Hybrid model’s prediction for the pN distribution exhibits a notable reduction in the Fermi motion peak compared to the ∆model, leading to better agreement with the data. It is noteworthy that the pN shape is similar between the 0 π and π0 measurements (figure 9b). However, despite the improvement over the ∆model, the Hybrid model still faces challenges in fully capturing the shape of the Fermi motion peak. Further insight is gained by replacing the initial state of LFG with ESF for the π0 measurement (see appendix Bfor full comparison plots with ESF). As shown in figure 10, the Hybrid model also shows improvement with ESF. However, while ESF results in a lower peak compared to LFG, it yields a higher χ2 , indicating that ESF also struggles to accurately capture the shape of the Fermi motion peak. Part of the improvement observed in the Hybrid model indeed arises from its refined depiction of interactions in higher W regions. This can be demonstrated by varying the – 14 – JHEP12(2024)141 0 50 100 150 (degree) T αδ 0 5 10 3− 10× (1/degree)σ/ T αδ /d σ d (a) MINERvA, area-normalised, LFG -P∆ π0 H-Pπ0 Dataπ0 -P∆ 0 π H-P 0 π Data 0 π 0 0.2 0.4 0.6 0.8 1 )c (GeV/ N p 0 1 2 3 4 )c (1/GeV/σ/ N p/d σ d (b) MINERvA, area-normalised, LFG -P∆ π0 H-Pπ0 Dataπ0 -P∆ 0 π H-P 0 π Data 0 π 0 0.2 0.4 0.6 0.8 1 )c (GeV/ N p 0 1 2 3 4 )c (1/GeV/σ/ N p/d σ d (b) MINERvA, area-normalised, LFG -P∆ π0 H-Pπ0 Dataπ0 -P∆ 0 π H-P 0 π Data 0 π Figure 9. Shape-only comparison summarising figure 8. Recall that LFG is only applied to non-QE channels; the initial state for QE is modelled by SF. 0 0.2 0.4 0.6 0.8 )c (GeV/ N p 0 1 2 3 4 39− 10× )c/GeV/ 2 (cm N p/d σ d 0 πMINERvA, /ndf=164/12 2 χ-P LFG ∆ /ndf=92/12 2 χH-P LFG /ndf=202/12 2 χ-P ESF ∆ /ndf=134/12 2 χH-P ESF Data Figure 10. Additional comparison to pNwith the initial state ESF. – 15 – JHEP12(2024)141 starting point ( Wmin ) and the stopping point ( Wmax ) of the H-P transition region (figure 11). Adjusting the transition window, ( Wmin, Wmax ), in the H-P model to the higher end gives more phase space to the Hybrid model and less to pythia. Clear trends are observed in figure 11: as both Wmin and Wmax increase, the χ2 decreases, indicating that the Hybrid model provides an improved description of the data compared to pythia, in the kinematic region of interest for the MINERvA π0 TKI analysis. When the transition begins at a sufficiently high W ( > 2 . 4 GeV ), the model consistently yields an optimized chi-squared value. This is consistent with the dσ dW results in figure 7, which indicate that the transition primarily occurs in the high- W region, thereby minimizing pythia’s influence on the final prediction. An intuitive illustration of this effect is shown in figure 12, where dσ dpN is calculated for a fixed transition width of 0 . 4 GeV at various W values. As W increases, the cross section decreases. In this particular case, it tends to converge to the prediction by the pure Hybrid mode. This reduced Fermi peak contributes to a more accurate description of the data. 5 Conclusions Neutrino pion production is an important interaction mechanism at GeV neutrino facilities. In this work, we describe the implementation of the Ghent single-pion production (SPP) model [ 13 ] in the neutrino event generator NuWro. Unlike the ∆-pythia (∆-P) implementation, until now used in NuWro, the resulting Hybrid-pythia (H-P) model essentially removes contributions from pythia to the SPP channel. The description of SPP is improved through extra resonance contributions in the second resonance region, interference with the non-resonant background of ref. [ 18 ], and a Regge description at high W . Predictions by the H-P model are compared to MINERvA and MicrooBooNE (and T2K, see appendix A) data. As expected, the refined modeling in regions of higher W tend to improve the agreement with the data. Notably, in comparison to the TKI measurements, the overprediction in the pN Fermi-motion peak region, as previously seen from the ∆-P model, has now been significantly reduced, indicating an improvement at the pion production level — the remaining overestimate, however, could be ambiguously attributed to either the production cross section or FSIs [58], as events can migrate from the pNpeak to the tail by modifying the latter. The shapes of the δαT and pN distributions indicate some universality between pionless and pion production; however, capturing the shape of the pN Fermi-motion peak turns out to be a challenge for current nuclear models implemented in MC generators. This observation should motivate further development of the initial state modelling for pion production, beyond the local Fermi gas and the effective spectral function approaches. Modification of the pion production mechanism due to in-medium effects may also lead to important changes in the predicted cross sections [ 59 ]. Medium modification of the ∆ resonance [ 41 ] are not explicitly included here and can lead to a non-negligible reduction of the cross section for MINERvA kinematics [ 30 ]. This is the case of the modification of the ∆ decay width, which in the model of [ 41 ] leads to the opening of pionless ∆decay channels. This mechanism is partially included in the NuWro cascade as pion absorption; however, further and dedicated studies on this subject are needed. Apart from nucleonic effects, there are nuclear effects that should be addressed when comparing with neutrino-nucleus cross section data. In this paper, we have discussed FSI included in the intranuclear cascade and, briefly, the medium-modification of the resonance – 16 – JHEP12(2024)141 5.7 4.3 4.1 4.4 3.9 3.7 3.6 3.5 4.1 3.5 4.3 3.6 3.5 3.5 3.4 1.6 2.0 2.4 2.8 3.2 (GeV) min W 1.6 2.0 2.4 2.8 3.2 (GeV) max W 3.5 4 4.5 5 5.5 2 χ , ndf=9 T αδ, 0 π(a) LFG, MINERvA 138.5 125.2 110.9 112.3 103.1 100.3 100.0 103.5 91.7 93.4 99.0 94.8 91.7 90.7 92.0 1.6 2.0 2.4 2.8 3.2 (GeV) min W 1.6 2.0 2.4 2.8 3.2 (GeV) max W 95 100 105 110 115 120 125 130 135 2 χ , ndf=12 N p, 0 π(b) LFG, MINERvA Figure 11. The χ2 value of the H-P predictions on the MINERvA π0 measurement with varying transition window, (Wmin, Wmax). properties. However, a more comprehensive understanding should incorporate other quantum mechanical effects not addressed here, like the distortion of the hadrons (i.e., elastic FSI for the outgoing nucleon and pion), which cannot be addressed by classical or semiclassical approaches like intranuclear cascade or models based on a factorized cross section; alternatively, these effects could be incorporated into the model describing the elementary vertex [ 60 ]. Currently, work is in progress to improve the Hybrid model by unitarizing the amplitude and incorporating additional contributions to the SPP mechanism, such as ρ - and ω -meson exchanged and higher mass resonances. The present work paves the road so that these and other potential improvements could be easily incorporated in NuWro. In the data-model comparison, we also appreciate future measurements that expand the phase space to higher W which is crucial to differentiate pion production models in the transitional region. – 17 – JHEP12(2024)141 0 0.2 0.4 0.6 0.8 )c (GeV/ N p 0 1 2 3 4 39− 10× /nucleon)c/GeV/ 2 (cm N p/dσd 0 πMINERvA CC 2.0 GeV−: 1.6 W 2.4 GeV−: 2.0 W 2.8 GeV−: 2.4 W 3.2 GeV−: 2.8 W Figure 12. The H-P model prediction of pN distribution for MINERvA CC π0 with different transition regions. The 2.8–3.2 GeV window is the nominal choice of our model. Particle Momentum pAngle θ µ−0.25 GeV/cto 7 GeV/c<70◦ π+0.15 GeV/cto 1.2 GeV/c<70◦ p0.45 GeV/cto 1.2 GeV/c<70◦ Table 2. The phase space cuts definition for T2K CC-π+TKI measurement [56]. Acknowledgments We thank Clarence Wret for pointing out an update of the MINERvA data of ref. [ 15 ]. R.G.-J. is supported by Project No. PID2021-127098NA-I00 funded by MCIN/AEI/10.13039/ 501100011033/FEDER, UE. N.J. and K.N. are supported by the Fund for Scientific Research Flanders (FWO) and by Ghent University Special Research Fund. X.L. is supported by the STFC (UK) Grant No. ST/S003533/2. A.N. is supported by Fermi Research Alliance, LLC under Contract No. DE-AC02-07CH11359 with the U.S. Department of Energy, Office of Science, Office of High Energy Physics. J.T.S. is supported by Polish Ministery of Science grant UMOWA 2023/WK/04. Q.Y. and Y.Z. are supported by National Natural Science Foundation of China (NSFC) under contract 12221005. A Comparison with T2K π+TKI data With a neutrino beam with an energy of approximately 0.6 GeV scattering off a hydrocarbon target, T2K has measured TKI in the CC- π+ production [ 56 ]: νµ + A→µ− + p + π+ + X , where X is the hadronic system that can contain nucleons but no mesons. The phase space cuts are defined in table 2. As shown in figure 13, the W -distribution is dominated by the ∆ ++ resonance, as expected given the low energy of the neutrinos. Figure 14 depicts the ∆-P and H-P predictions, showing limited improvement by the latter due to the low- W dominance of the sample. As noted by ref. [ 56 ], the inherent challenge in describing the data appears to stem from the initial state of LFG. – 18 – JHEP12(2024)141 1 1.5 2 2.5 (GeV)W 0 0.5 1 1.5 39− 10× /GeV) 2 (cmW/d σ d ) 2 cm 40− /10σ, LFG ( + πTKI H-P (1.6) H-P: Hybrid (1.5) H-P: Hybrid: Hydrogen (0.4) -P (1.6)∆ (1.4)∆-P: ∆ : Hydrogen (0.4)∆-P: ∆ QE (0.01) Figure 13. Similar to figure 3but for the T2K π+TKI measurement [56]. 0 50 100 150 (degree) T α δ 0 0.5 1 1.5 2 42− 10× /degree/nucleon) 2 (cm T α δ /d σ d /ndf=3.6/3 2 χ-P LFG ∆(a) hydrogen + π1 πmulti- + πT2K 0 50 100 150 (degree) T α δ 0 0.5 1 1.5 2 42− 10× /degree/nucleon) 2 (cm T α δ /d σ d /ndf=3.5/3 2 χ(b) H-P LFG hydrogen + π1 πmulti- + πT2K 0 0.5 1 1.5 (GeV/c) N p 0 0.2 0.4 0.6 0.8 39− 10× /nucleon)c/GeV/ 2 (cm N p/d σ d 10× /ndf=25.6/4 2 χ-P LFG ∆(c) + π1 πmultihydrogen + πT2K 0 0.5 1 1.5 (GeV/c) N p 0 0.2 0.4 0.6 0.8 39− 10× /nucleon)c/GeV/ 2 (cm N p/d σ d 10× /ndf=22.4/4 2 χ(d) H-P LFG + π1 πmultihydrogen + πT2K Figure 14. Similar to figure 8lower panels but for the T2K π+TKI measurement [56]. B Comparison with MINERvA TKI data using ESF For completeness, the comparison with the MINERvA data using the initial state ESF is fully shown in figure 15. C Single pion production dσ dWwith 3 GeV neutrino To illustrate model behavior, a neutrino energy of 8 GeV is used in figure 2. For comparison, figure 16 presents the same analysis for a lower energy of 3 GeV . As can be seen, at 3 GeV , the contribution from the transition region is diminished. – 19 – JHEP12(2024)141 0 50 100 150 (degree) T αδ 0 5 10 15 20 42− 10× /degree/nucleon) 2 (cm T αδ/dσd /ndf=19/12 2 χ-P ESF ∆(a) QE (SF) 0 π1 + π1 2p2h (ESF) πMINERvA 0 0 50 100 150 (degree) T αδ 0 5 10 15 20 42− 10× /degree/nucleon) 2 (cm T αδ/dσd /ndf=19/12 2 χ(b) H-P ESF QE (SF) 0 π1 + π1 2p2h (ESF) πMINERvA 0 0 0.2 0.4 0.6 0.8 )c (GeV/ N p 0 2 4 6 39− 10× /nucleon)c/GeV/ 2 (cm N p/dσd /ndf=61/24 2 χ-P ESF ∆(c) QE (SF) 0 π1 + π1 2p2h (ESF) πMINERvA 0 0 0.2 0.4 0.6 0.8 )c (GeV/ N p 0 2 4 6 39− 10× /nucleon)c/GeV/ 2 (cm N p/dσd /ndf=56/24 2 χ(d) H-P ESF QE (SF) 0 π1 + π1 2p2h (ESF) πMINERvA 0 0 50 100 150 (degree) T αδ 0 5 10 42− 10× /degree/nucleon) 2 (cm T αδ/dσd /ndf=4.1/9 2 χ-P ESF ∆(e) QE (SF) 0 π1 + π1 πmulti- 0 πMINERvA 0 50 100 150 (degree) T αδ 0 5 10 42− 10× /degree/nucleon) 2 (cm T αδ/dσd /ndf=3.5/9 2 χ(f) H-P ESF QE (SF) 0 π1 + π1 πmulti- 0 πMINERvA 0 0.2 0.4 0.6 0.8 )c (GeV/ N p 0 2 4 39− 10× /nucleon)c/GeV/ 2 (cm N p/dσd /ndf=202/12 2 χ-P ESF ∆(g) QE (SF) 0 π1 + π1 πmulti- 0 πMINERvA 0 0.2 0.4 0.6 0.8 )c (GeV/ N p 0 2 4 39− 10× /nucleon)c/GeV/ 2 (cm N p/dσd /ndf=134/12 2 χ(h) H-P ESF QE (SF) 0 π1 + π1 πmulti- 0 πMINERvA Figure 15. Similar to figure 8but with the non-QE initial state modelled by ESF. – 20 – JHEP12(2024)141 1 2 3 4 (GeV)W 0 20 40 60 39− 10× /GeV) 2 (cmW/d σ d + π p − µ → p µ ν(a) 3 GeV, H-P Hybrid -P∆ ∆ PYTHIA 1 2 3 4 (GeV)W 0 5 10 15 39− 10× /GeV) 2 (cmW/d σ d 0 π p − µ → n µ ν(b) 3 GeV, H-P Hybrid -P∆ ∆ PYTHIA 1 2 3 4 (GeV)W 0 5 10 15 39− 10× /GeV) 2 (cmW/d σ d + π n − µ → n µ ν(c) 3 GeV, H-P Hybrid -P∆ ∆ PYTHIA Figure 16. ∆and Hybrid model predictions for SPP by an 3GeV neutrino off the nucleon: (a) νµp→µ−pπ+ , (b) νµn→µ−pπ0 , and (c) νµn→µ−nπ+ . The trend and difference between models are similar to the figure with 8GeV neutrino in figure 2. But the lower energy suppressed the contribution in the transition region. – 21 – JHEP12(2024)141 Data Availability Statement. This article has no associated data or the data will not be deposited. Code Availability Statement. This article has no associated code or the code will not be deposited. Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited. References [1] DUNE collaboration, Deep Underground Neutrino Experiment (DUNE), Far Detector Technical Design Report, Volume I Introduction to DUNE,2020 JINST 15 T08008 [arXiv:2002.02967] [INSPIRE]. [2] Hyper-Kamiokande Protocollaboration, Physics potential of a long-baseline neutrino oscillation experiment using a J-PARC neutrino beam and Hyper-Kamiokande,PTEP 2015 (2015) 053C02 [arXiv:1502.05199] [INSPIRE]. [3] NOvA collaboration, First measurement of muon-neutrino disappearance in NOvA,Phys. Rev. D93 (2016) 051104 [arXiv:1601.05037] [INSPIRE]. [4] T2K collaboration, The T2K Experiment,Nucl. Instrum. Meth. A 659 (2011) 106 [arXiv:1106.1238] [INSPIRE]. [5] T2K collaboration, Constraint on the matter-antimatter symmetry-violating phase in neutrino oscillations,Nature 580 (2020) 339 [Erratum ibid. 583 (2020) E16] [arXiv:1910.03887] [INSPIRE]. [6] Super-Kamiokande collaboration, Evidence for oscillation of atmospheric neutrinos,Phys. Rev. Lett. 81 (1998) 1562 [hep-ex/9807003] [INSPIRE]. [7] JUNO collaboration, Neutrino Physics with JUNO,J. Phys. G 43 (2016) 030401 [arXiv:1507.05613] [INSPIRE]. [8] K.M. Graczyk, D. Kiełczewska, P. Przewłocki and J.T. Sobczyk, CA 5axial form factor from bubble chamber experiments,Phys. Rev. D 80 (2009) 093001 [arXiv:0908.2175] [INSPIRE]. [9] T. Sjöstrand, S. Mrenna and P.Z. Skands, PYTHIA 6.4 Physics and Manual,JHEP 05 (2006) 026 [hep-ph/0603175] [INSPIRE]. [10] A. Bodek and U.K. Yang, Modeling neutrino and electron scattering inelastic crosssections in the few GeV region with effective LO PDFs TV Leading Order, in the proceedings of the 2nd International Workshop on Neutrino-Nucleus Interactions in the Few GeV Region, Irvine, U.S.A., December 12–15 (2002) [hep-ex/0308007] [INSPIRE]. [11] J.C. Collins, D.E. Soper and G.F. Sterman, Factorization of Hard Processes in QCD,Adv. Ser. Direct. High Energy Phys. 5(1989) 1 [hep-ph/0409313] [INSPIRE]. [12] K. Niewczas et al., Angular distributions in Monte Carlo event generation of weak single-pion production,Phys. Rev. D 103 (2021) 053003 [arXiv:2011.05269] [INSPIRE]. [13] R. González-Jiménez et al., Electroweak single-pion production off the nucleon: from threshold to high invariant masses,Phys. Rev. D 95 (2017) 113007 [arXiv:1612.05511] [INSPIRE]. – 22 – JHEP12(2024)141 [14] A. Nikolakopoulos, R. González-Jiménez, N. Jachowicz and J.M. Udías, Assessing the theory-data tension in neutrino-induced charged pion production: The effect of final-state nucleon distortion,Phys. Rev. D 107 (2023) 053007 [arXiv:2210.12144] [INSPIRE]. [15] MINERvA collaboration, Charged Pion Production in νµ Interactions on Hydrocarbon at ⟨Eν⟩ = 4.0 GeV,Phys. Rev. D 92 (2015) 092008 [arXiv:1406.6415] [INSPIRE]. [16] MicroBooNE collaboration, Measurement of the differential cross section for neutral pion production in charged-current muon neutrino interactions on argon with the MicroBooNE detector,Phys. Rev. D 110 (2024) 092014 [arXiv:2404.09949] [INSPIRE]. [17] MINERvA collaboration, Probing nuclear effects with neutrino-induced charged-current neutral pion production,Phys. Rev. D 102 (2020) 072007 [Erratum ibid. 110 (2024) 059903] [arXiv:2002.05812] [INSPIRE]. [18] E. Hernández, J. Nieves and M. Valverde, Weak Pion Production off the Nucleon,Phys. Rev. D 76 (2007) 033005 [hep-ph/0701149] [INSPIRE]. [19] D. Drechsel, S.S. Kamalov and L. Tiator, Unitary Isobar Model — MAID2007,Eur. Phys. J. A 34 (2007) 69 [arXiv:0710.0306] [INSPIRE]. [20] H. Kamano, S.X. Nakamura, T.S.H. Lee and T. Sato, Isospin decomposition of γN →N∗ transitions within a dynamical coupled-channels model,Phys. Rev. C 94 (2016) 015201 [arXiv:1605.00363] [INSPIRE]. [21] S.X. Nakamura, H. Kamano and T. Sato, Dynamical coupled-channels model for neutrino-induced meson productions in resonance region,Phys. Rev. D 92 (2015) 074024 [arXiv:1506.03403] [INSPIRE]. [22] J.E. Sobczyk et al., Angular distributions in electroweak pion production off nucleons: odd parity hadron terms, strong relative phases and model dependence,Phys. Rev. D 98 (2018) 073001 [arXiv:1807.11281] [INSPIRE]. [23] O. Lalakulich, E.A. Paschos and G. Piranishvili, Resonance production by neutrinos: The Second resonance region,Phys. Rev. D 74 (2006) 014009 [hep-ph/0602210] [INSPIRE]. [24] L. Alvarez-Ruso, E. Hernández, J. Nieves and M.J. Vicente Vacas, Watson’s theorem and the N∆(1232) axial transition,Phys. Rev. D 93 (2016) 014016 [arXiv:1510.06266] [INSPIRE]. [25] E. Hernández et al., nu induced threshold production of two pions and N*(1440) electroweak form factors,Phys. Rev. D 77 (2008) 053009 [arXiv:0710.3562] [INSPIRE]. [26] M. Guidal, J.M. Laget and M. Vanderhaeghen, Pion and kaon photoproduction at high-energies: Forward and intermediate angles,Nucl. Phys. A 627 (1997) 645 [INSPIRE]. [27] M.M. Kaskulov and U. Mosel, Deep exclusive charged πelectroproduction above the resonance region,Phys. Rev. C 81 (2010) 045202 [arXiv:1001.1952] [INSPIRE]. [28] T. Vrancx and J. Ryckebusch, Charged pion electroproduction above the resonance region,Phys. Rev. C 89 (2014) 025203 [arXiv:1310.7715] [INSPIRE]. [29] T. Vrancx, L. De Cruz, J. Ryckebusch and P. Vancraeyveld, Consistent interactions for high-spin fermion fields,Phys. Rev. C 84 (2011) 045201 [arXiv:1105.2688] [INSPIRE]. [30] R. González-Jiménez, K. Niewczas and N. Jachowicz, Pion production within the hybrid relativistic plane wave impulse approximation model at MiniBooNE and MINERvA kinematics, Phys. Rev. D 97 (2018) 013004 [arXiv:1710.08374] [INSPIRE]. [31] A. Nikolakopoulos et al., Modeling neutrino-induced charged pion production on water at T2K kinematics,Phys. Rev. D 97 (2018) 093008 [arXiv:1803.03163] [INSPIRE]. – 23 –