Revisiting the T2K data using different models for the neutrino-nucleus cross sections

We present a three-flavour fit to the recent \nu\mu -->\nu e and \nu\mu -->\nu\mu T2K oscillation data with different models for the neutrino-nucleus cross section. We show that, even for a limited statistics, the allowed regions and best fit points in the (\theta_{13},\delta_{CP}) and (\theta_{23},\Delta m^2_{atm}) planes are affected if, instead of using the Fermi Gas model to describe the quasielastic cross section, we employ a model including the multinucleon emission channel.


Introduction
Recently the T2K collaboration has released data in both ν µ → ν e appearance [1] and ν µ → ν µ disappearance [2] modes; in the first case, six events passed all the selection criteria, implying (under the assumption of a normal ordering of the neutrino mass eigenstates): with the CP phase δ CP undetermined. In the disappearance channel, the 31 events collected by T2K are fitted with: (sin 2 2θ 23 ) T 2K = 0.98 |∆m 2 atm | T 2K = 2.65 · 10 −3 eV 2 .
The aim of this work is to reanalyse the T2K data to assess the impact of different models for the ν-nucleus cross sections on the determination of oscillation parameters. This work can be considered as a generalization of Ref. [3], where the impact of different modelizations of quasielatic cross sections in the low-gamma beta-beam regime was analyzed. In the present case we consider two different models involving not only quasielastic but also pion production and inclusive cross sections. On one hand, we choose a model as similar as possible to the one used by the T2K collaboration. They simulate the neutrino-nucleus interaction using the NEUT Monte Carlo Generator [4]. Even if we do not know the details of the last tunings performed by the collaboration to take into account for the recent measurements of K2K [5,6], MiniBooNE [7,8] and SciBooNE [9,10], we treat the several exclusive channels using the same models implemented in NEUT. As a consequence, we consider the Fermi Gas [11] for the quasielastic channel and the Rein and Sehgal model [12] for pion production. The second model considered in our analysis is the one of Martini, Ericson, Chanfray and Marteau [13], in the following called "MECM model". It is based on the nuclear response functions calculated in random phase approximation and allows an unified treatment of the quasielastic, the multinucleon emission channel and the coherent and incoherent pion production. The agreement with the experimental data in the pion production channels [6,7,9] has been proved. Nevertheless the main feature of this MECM model is the treatment of the multinucleon emission channel in connection with the quasielastic. In fact, as suggested in [13,14], the inclusion of this channel in the quasielastic cross section is a possible explanation of the MiniBooNE quasielastic total cross section [8], apparently too large with respect to many theoretical predictions [15] employing the standard value of the axial mass. Since the MiniBooNE experiment, as well as many others involving Cherenkov detectors, defines a "quasielastic" event as the one in which only a final charged lepton is detected, the ejection of a single nucleon (a genuine quasielastic event) is only one possibility, and one must in addition consider events involving a correlated nucleon pair from which the partner nucleon is also ejected. This leads to the excitation of 2 particle-2 hole (2p-2h) states; 3p-3h excitations are also possible. Nowadays other models [16,17,18,19] have included the multinucleon contribution in the computation of the cross sections relevant for the MiniBooNE quasielastic kinematics, improving the agreement with the experimental data. For a brief review see for example [20]. Recently, it has been shown [21] that the MECM model can also reproduce the MiniBooNE flux averaged double differential cross section [8] which is a directly measured quantity and hence free from the model-dependent uncertainties in the neutrino energy reconstruction, and the total inclusive cross section [20] (also employed by T2K as described below) measured by SciBooNE [10]. In the following we will use the cross sections obtained in the two different approaches described above in several exclusive channels (quasielastic and pion production), as well as in the inclusive one, for both charged current (CC) and neutral current (NC) interactions on carbon and oxygen (the targets used in near and far T2K detectors, respectively) and for two neutrino flavours ν µ and ν e . Although all exclusive channels are involved in the analysis, we will refer to the first model as "the Fermi Gas model" and to the second approach as "the MECM model".
In order to perform our comparison among the above-mentioned models, we first need to correctly normalize the Fermi gas to the T2K event rates, at both near (ND) and far (FD) detectors; we use the following algorithm: 1-normalization of the cross section with the ν µ inclusive CC at the ND; according to [1], we have to reproduce 1529 ν µ inclusive events, collected using 2.9 × 10 19 POT, in the energy range [0 − 5] GeV, with an active detector mass of 1529 Kg 1 at a distance of 280 m from the ν source and half a year of data taking (Run 1). Notice that only the muon neutrino cross sections can be correctly normalized; we assume that the same normalization also applies for the ν e cross section, although they could differ at the µ production threshold (in any case away from the peak of the neutrino flux); 2-computation of the expected events (and energy distributions) at the far detector in the appropriate two-parameter plane ((sin 2 2θ 13 , δ CP ) for appearance and (θ 23 , ∆m 2 atm ) for disappearance); 3-normalization to the T2K spectral distributions.
Step #3 is needed to get rid of the experimental efficiencies applied by the T2K collaboration to the signal and background events. This means that the bin contents of our simulated distributions (obtained at point #2) are corrected by coefficients, generally of O(1) that we consider as a detector property, and then not further modified. For a different model, we repeat step #1 and then go to step #2, using the same normalization coefficients extracted in step #3 with the Fermi gas. We make use of the GloBES [22] and MonteCUBES [23] softwares for the computation of event rates (and related χ 2 functions) expected at the T2K ND and FD detectors. The fluxes of ν µ , ν e and their CP-conjugate counterparts predicted at the FD in absence of oscillations have been extracted directly from Fig.1 of [1], whereas the ν µ flux at the ND has been obtained from [2]. Such fluxes (the relevant ones summarized in Fig.1) are given for 10 21 POT. As already stressed, for the relevant cross sections we assumed  Figure 1: Fluxes at near (ν µ only) and far (ν µ and ν e ) detectors. 1 We thank Scott Oser for providing such a number to us. that the T2K collaboration uses some "sophisticated" version of the Fermi gas model [11]. In Fig.2 we show the inclusive and QE cross sections in the FG model (dashed lines) and in the MECM model (solid line) used in our simulation, after having correctly normalized the inclusive cross sections to the event rate at the ND. Especially for the MECM model, this procedure involves a degree of extrapolation of the inclusive cross sections towards neutrino energies beyond the validity of model itself. However, neutrino fluxes above O(1) GeV drop very fast and we checked that different kind of extrapolations do not alter our conclusions.  The important feature here is that, even after the normalization procedure, the MECM CCQE cross section is still larger than the FG predictions, in the energy range relevant for appearance studies. This is due to the inclusion of the multinucleon component and will be the main reason of the differences between the results obtained in the two models. Note on the contrary that the inclusive cross sections are not really different.
We clearly see that CP violating effects are encoded in the interference term proportional to the product of the solar mass splitting and the baseline, implying a scarce dependence of this facility on δ CP when only the ν µ → ν e channel (and the current luminosity) is considered.

Extracting the T2K data
Events in the far detector (obtained with 2.9 × 10 20 POT) are ν e CCQE from ν µ → ν e oscillation, with main backgrounds given by ν e contamination in the beam and neutral current events with a misidentified π 0 . The experimental data have been grouped in 5 reconstructed-energy bins, from 0 to 1.25 GeV and they are summarized in Tab.1. The expectations for signal and backgrounds have been computed by the T2K collaboration from Monte Carlo simulations, for fixed value of the oscillation parameters, namely sin 2 2θ 12 = 0.8794, sin 2 2θ 13 = 0.1, sin 2 2θ 23 = 1 and ∆m 2 sol = 7.5 × 10 −5 eV 2 , ∆m 2 atm = +2.4 × 10 −3 eV 2 . In order to normalize our event rates to the T2K Monte Carlo expectations, we extracted these numbers from  Notice that we used the central bin energy as a reference value for the neutrino energy in a given bin; this could be different from the reconstructed neutrino energies used by the T2K collaboration. To mimic possible uncertainties associated to the neutrino energy reconstruction, we apply an energy smearing function to distribute the rates in the various energy bins. Other choiches, more related to microscopical calculations [25,26,27] are also possible. In particular, an analysis of the validity of the approximation contained in the identification of the reconstructed neutrino energy via a two-body kinematics with the real neutrino energy has been done in [27], where the MECM model has been employed. The role of several nuclear effects such as Pauli blocking, Fermi motion, RPA correlations and multinucleon component has been studied in details. This analysis was performed, among others, considering the T2K conditions at near and far detectors, paying a particular attention to the ν µ → ν e appearance mode. The ratios among our computation and the T2K data (energy dependent efficiencies) are summarized in Tab.2. This procedure (corresponding to step #3 of the previous paragraph) allows us to take into account all the detection efficiencies to different neutrino flavours in the Super Kamiokande detector. Once computed, these corrective factors are used in the simulations done with a different cross section, since we assume here that they are features of the detector and not of the neutrino interactions. We see that for ν e,µ → ν e transitions these numbers are just O(1) coefficients, which makes channel bin 1 bin 2 bin 3 bin 4 bin 5 ν µ → ν e 1.  us confident that the normalization procedure correctly accounts for the main experimental features. The same is not true for the NC events which, however, have not been normalized to the ND as for the CC interactions. As a check, we also computed the expected events for sin 2 2θ 13 =0, obtaining 0.1 ν µ → ν e events and 0.72 ν e → ν e events (and the same neutral current rate), in good agreement with the T2K expectations [1].

Fit to the data
Equipped with these results, we performed a χ 2 analysis to reproduce the confidence level regions in the (sin 2 2θ 13 , δ CP )-plane shown in Fig.6 of [1]. Contrary to what has been done in the official T2K paper, we make a complete three-neutrino analysis of the experimental data, marginalizing over all parameters not shown in the confidence regions. As external input errors, we used 3% on θ 12 and ∆m 2 sol , 8% on θ 23 and 6% on ∆m 2 atm . We use a constant energy resolution function σ(E ν ) = 0.085 and, for simplicity, we adopt a 7% normalization error for the signal and 30% for the backgrounds. We also used energy calibration errors fixed to 10 −4 for the signal and 5 · 10 −2 for the backgrounds; normalization and energy calibration errors take into account the impact of systematic errors in the χ 2 computation. Assuming a normal hierarchy spectrum, the best-fit point from the fit procedure is (obviously): with χ 2 min = 1.69; the related contour plot is shown in Fig.3. Compared to the official release, the plot is in quite good agreement, although the allowed values of θ 13 around maximal CP violation δ CP = π/2 are a bit larger (this is the effect of including the errors of the atmospheric parameters into the fit procedure). We now apply the same procedure to determine θ 13 using the MECM cross sections described in [13]. In doing that, we normalize the cross sections to the ND events and then compute the number of oscillated events (and related backgrounds), to be compared with the experimental T2K data. We assume that the efficiencies reported in Tab.2 are exactly the same, since they are a property of the SK detectors and then independent on the cross section model. With these assumptions, we get the following number of expected rates for sin 2 2θ 13 =0.1: channel bin 1 bin2 bin3 bin4 bin5 total estimates ν µ → ν e 0.  It is clear that larger rates need smaller θ 13 to reproduce the data (the effect of the CP phase δ is negligible with such a statistics). The best fit point is: with χ 2 min = 1.53, and the contour plot is shown in Fig.4. We can appreciate a substantial improvement in the determination of the reactor angle, whose largest value is 0.24, to be compared with 0.35 obtained with the Fermi Gas. To make a more direct comparison on θ 13 between the FG and MECM results, in Fig.5 we show the χ 2 − χ 2 min function, computed marginalizing over all other oscillation parameters (including δ CP ). At 1σ, we get: sin 2 2θ F G 13 = 0.12 ( +0.08 −0.09 ) .
They are clearly compatible although, as expected, θ M ECM 13 < θ F G 13 .

The disappearance channel
We extend the previous analysis to include the first disappearance ν µ → ν µ data [2]. In the two-flavour limit, (the one where both θ 13 and ∆m 2 sol are vanishing) the ν µ → ν µ probability reads [28]: Effects related to θ 13 are clearly sub-dominant, so that this channel is particularly useful to extract information on the atmospheric parameters. The T2K collaboration collected 31 data events, grouped in 13 energy bins, as one can see from Fig.3 of [2]. The sample extend up to 6 GeV and it is mainly given by ν µ CCQE, ν µ CC non-QE, ν e CC and NC. Differently from the appearance channel, we cannot normalize their energy distribution to the channelby-channel T2K Monte Carlo expectation since, as far as we know, such information has not been released. We can only normalize our FG cross section to the total rates shown in Tab.I of [2], which amounts to 17.3, 9.2, 1.8 and <0.1 events for ν µ CCQE, ν µ CC non-QE, NC and ν e CC, respectively. Such numbers refer to sin 2 (2θ 23 ) = 1.0 and |∆m 2 atm | = 2.4 × 10 −3 eV 2 , with all other neutrino mixing parameters vanishing. For the sake of completeness, we summarize in Tab.4 the T2K data as well as the energy distributions of the ν µ CCQE and ν µ CC non-QE as obtained using the MECM cross sections. In evaluating such numbers, we assume a variable bin size, centered in the neutrino energy corresponding to the T2K data. In our fit procedure we have assumed a total normalization for the NC as given in [2], but with appropriate energy distributions according to the FG and MECM cross sections. We have also adopted a conservative 15% normalization error and energy calibration error at the level of 10 −3 for both signal and background. The results of our fit procedure can  Table 4: T2K events and bin distributions for the ν µ CCQE and ν µ CC non-QE rates in the MECM model. be appreciated in Fig.6, where we show the 90% CL for the Fermi Gas (dashed line) and the MECM model (solid line), in the case of normal hierarchy. We plot the 2 degrees of freedom (dof) confidence levels in the (θ 23 , ∆m 2 atm ) (left panel) and (sin 2 2θ 23 , ∆m 2 atm ) (right panel, to facilitate the comparison with the official T2K results) planes. Again, the plots have been obtained marginalizing over the not shown parameters (a full three-flavour analysis); we considered a 50% error on sin 2 2θ 13 (with best fit at sin 2 2θ 13 = 0.0059) and δ CP undetermined. We obtained:  Some comments are in order; first of all, we observe that, for both models, the best fit point is different from the T2K ones, which is (sin 2 2θ 23 ) T 2K = 0.98 |∆m 2 atm | T 2K = 2.65 · 10 −3 eV 2 ; this is somehow obvious since we normalized our events to the MC predictions obtained for a different set of atmospheric parameters. The MECM cross section gives a better determination of both θ 23 and ∆m 2 atm , mainly due to the larger statistics than the FG; at the same time, the disappearance probability in Eq.(7), for negligible solar mass difference and reactor angle, is smaller if the atmospheric mass difference is larger, for fixed sin 2 2θ 23 . This is what happens here, where a smaller P (ν µ → ν µ ) (and then a larger ∆m 2 atm ) is needed in the MECM model to partially compensate for the larger cross section.

Future perspectives
The statistics used by the T2K collaboration to make the disappearance study is only a 2% of the rates expected at the end of the experiment. It makes sense to ask how the previous results would modify if the accumulated statistics would be larger than the current one. We limit ourselves to consider a number of events with the same energy distribution as the experimental ones but bin contents larger by factor of 10, in both appearance and disappearance channels. In the analysis of the appearance channel, the (weak) information on θ 13 contained in the disappearance sample should not be neglected (as we did previously); at the same time, the dependence on the atmospheric parameters from the appearance sample cannot in principle be neglected when studying the disappearance data. For this reason, we prefer to combine both ν µ → ν e and ν µ → ν µ oscillation data, and study the sensitivity to the reactor and atmospheric parameters as we did in the previous sections, marginalizing over the parameters not expressly shown. Notice that such an approach would not give any additional information on the mixing parameters if adopted with the current T2K statistics: in fact, we see from Fig.6 that the uncertainties on θ 23 and ∆m 2 atm obtained from the T2K data are larger than the adopted external errors on these parameters in the appearance channel, so that adding the ν µ → ν µ data will not improve the sensitivity to θ 13 ; similarly, the dependence on the reactor angle in P (ν µ → ν µ ) is sub-leading and the impact of the disappearance channel in the appearance measurement is also negligible. We stress that extracting information on the mixing parameters combining appearance and disappearance channels is also mandatory to avoid some inconsistencies emerged in the official T2K fits, where |∆m 2 atm | is fixed to 2.4 × 10 −3 eV 2 in the appearance analysis whereas the best fit point obtained from the disappearance data is 2.6 × 10 −3 eV 2 . The results of our procedure are shown in Fig.7, where we display the 90% CL in the (sin 2 2θ 13 , δ CP )-plane (left panel) and (θ 23 , ∆m 2 atm )-plane (right panel) for the MECM (solid line) and FG (dashed one) models in the case the current T2K statistics is increased by a factor of 10. The minimum of the χ 2 still gets reasonable values: we obtain χ 2 min ∼ 20 in the appearance analysis and χ 2 min ∼ 30 in disappearance. In both panels we can appreciate the effects of the increased statistics, as  expected: for a given model of cross section, the allowed regions are strongly restricted with respect to the current situation. The best fit values for δ CP are somehow different in the two models (δ CP ∼ 0 and δ CP ∼ 144 • for the FG and MECM models, respectively), although statistically non very significant. Such intervals for sin 2 2θ 13 , θ 23 and ∆m 2 atm are summarized in Tab.5 (δ CP is obviously still unconstrained). We have checked that, if we only use the appearance channel to extract θ 13 , the results are sligthly different: although the best fit value is practically indistinguishable from the one quoted in Tab.5, the confidence regions are a bit larger, with significant overlap with the above mentioned analysis. To see stronger effects due to the θ 13 dependence in the ν µ → ν µ transition, we need a more accurate spectral information [29]. Similar conclusions can also be drawn for the disappearance channel: with only a factor of 10 more statistics and no appearance contribution, the best fit for the atmospheric parameters remain almost the same whereas the 90% CL region for θ 23 shows a smaller lower limit (from 40.1 • to 39.8 • ) in the FG model.
Finally, we observe that such an increased statistics is necessary to make marginally  Table 5: 90% intervals for sin 2 2θ 13 , θ 23 and ∆m 2 atm , for the MECM and FG models in the case the current T2K statistics is increased by a factor of 10. In parenthesis, the best fit points.
incompatible the FG and MECM sin 2 2θ 13 results, see Fig.8, obtained marginalizing over δ CP also. In fact, at 1σ we get:

The inverted hierarchy case
For the sake of completeness, we have repeated the same computations as above under the hypothesis that the neutrino mass spectrum is of inverted type (IH). With the current T2K statistics, we cannot appreciate huge differences in the results obtained using the two different models for the cross section. Then, we limit ourselves here to the case where the stastistics is larger by a factor of 10. Our results are summarized in Fig.9 and Tab.6. Comparing the left panel of Fig.9 with the corresponding one in Fig.7, we clearly see that an inverted spectrum prefers larger values for θ 13 , in both models. The best fit of the CP phases is different among the two mass orderings but not really significant. In the atmospheric plane, right panel of Fig.9, we observe the same pattern as in the normal hierarchy case, that is the MECM tends to give a better resolution for both ∆m 2 atm and θ 23 than the FG model.    Table 6: 90% intervals for sin 2 2θ 13 , θ 23 and ∆m 2 atm , for the MECM and FG models in the case the current T2K statistics is increased by a factor of 10 and the "data" are fitted with the inverted hierarchy. In parenthesis, the best fit points.

Conclusions
In this paper we have studied the impact of using different models for the neutrino-nucleus cross section in the determination of the θ 13,23 mixing angles and the atmospheric mass difference ∆m 2 atm using the recent T2K data, for both appearance and disappearance channels. Although the statistics is not large enough to draw definite conclusions, we have seen that a more refined treatments of nuclear effects in neutrino interactions can have some impact in the achievable precision on the mixing parameters. In particular, the MECM model predicts a large CCQE cross section, compared to the FG model, which results in a small θ 13 needed to fit the data in the ν µ → ν e channel. At the same time, a larger ∆m 2 atm is required to fit the data in the ν µ disappearance channel, since a smaller disappearance probability is needed to compensate for the larger cross sections. Interestingly enough, with 10 times more statistics the two models tend to give substantial different results in terms of best fit points and parameter uncertainties (of course, better than before) but their predictions are still compatible to each other.