On the systematic uncertainties in DUNE and their role in New Physics studies

In the recent years experiments have established the existence of neutrino oscillations and most of the oscillation parameters have been measured with a good accuracy. The search for New Physics in neutrino oscillation will be an experimental concrete possibility in the next future. In this paper we investigate the ability of the DUNE facility to search for Non Standard Interaction (NSI) in neutrino propagation in matter, emphasizing the role of different assumptions on the shape and absolute normalization errors of both νe and νμ signals. We also study in detail the effects of NSI and systematics in the extraction of standard oscillation parameters.


Introduction
The standard neutrino mixing angles and mass differences have been determined with very good accuracy [1]- [3] apart from the yet unknown mass ordering, the sign of the atmospheric mass squared difference, the octant of the atmospheric mixing angle well as the CP-violating phase δ CP . Non-oscillation experiments have also provided important information on the neutrino cross sections on nuclei [4], which is triggering by itself an intense field of research. Excluding some hints for neutrino oscillation into additional sterile neutrino states [5]- [6], there are no experimental evidences that neutrinos possess non-standard properties beyond those described by the Standard Model (SM). Thus, the presence of new phenomena in the neutrino sector is still an open question which attracts a lot of attention in both theoretical and experimental communities, absorbed by the investigation of models able to produce small perturbations on top of the usual three-flavor oscillation and by the requirements of the experimental facilities needed to detect them. In such an intriguing panorama, a crucial role is played by the systematic uncertainties which affect the sensitivity of a given experiment to oscillation parameters and, even more seriously, to New Physics. For this reason, the relevant sources of uncertainty (like flux normalization, energy calibration, nuclear effects, detector efficiency and so on) must be contemplated in any simulation of future detectors aiming at testing the standard paradigm of neutrino oscillation [7]. It is well known that correlated systematics can be reduced by using the unoscillated event distribution at Near Detector to predict the distribution at the Far Detector; this strategy works very well for short-baseline experiments and has been used with great success in reactor experiments to measure θ 13 . Even though the Near to Far extrapolation at long baselines does not seem to provide the same advantages in reducing the systematic uncertainties (because, for example, the absolute cross section is known with an accuracy of 10-20%, the kinematics is more difficult to treat because of JHEP08(2018)028 several interaction mechanisms and in appearance measurements initial and final states are different), it is nevertheless of crucial importance to estimate what the most important sources of systematics are and, consequently, to understand the needs and requirements for the future generation of Near Detectors [8]. Since the final answer obviously depends on the chosen experimental facility and on what kind of new physics scenario we want to address, we start here this ambitious program checking which reasonable combinations of shape 1 and absolute normalization uncertainties of both ν e appearance (app) and ν µ disappearance (dis) signals can ensure the best performance in ameliorating the current bounds on the matter Non Standard Interaction (NSI) parameters (see [9] for a recent review) at the DUNE facility [10]. This will be the subject of section 2. At the same time, we update the existing bounds on NSI using the latest results from the global fits on standard oscillation parameters, section 3, as well as the impact of marginalizing over NSI parameters in the estimate of the octant, CP and mass hierarchy determination sensitivities, section 4. Notice that a similar interesting analysis was presented in [11] where, however, only the overall systematics for signal and background normalization were allowed to vary. In this work we do not attempt to take simultaneously into account possible NSI affecting production and/or detection [12] as we consider them as subleading effects compared to the matter effect modifications [13].
NSI affecting neutrino propagation are generated through four-fermion operators involving Dirac bilinears of neutrinos and of fermions in the Earth matter, which modify the neutrino evolution equation in the flavor basis according to [14]- [16]: where ∆ ij = ∆m 2 ij /2E, U P M N S is the neutrino mixing matrix, A ≡ 2 √ 2G F n e and αβ are effective parameters encoding the action of the four-fermion operators. Thus, beside the standard oscillation parameters, the parameter space is enriched by six more moduli | αβ | and three more phases φ αβ , on which no constraints have been derived so far. A reduction of the number of independent moduli can be achieved subtracting to the whole Hamiltonian the matrix τ τ × 1 and redefining ee − τ τ → ee and µµ − τ τ → µµ ; in this case τ τ is set to zero in our numerical simulations. The current 90% Confidence Level (CL) bounds on the real part of the NSI parameters used throughout this paper are reported in table 1 and extracted from [17]. No constraints on the new CP phases are available so far.
If θ 23 is different from π/4 [18] (but see also [19]- [20] for a detailed perturbative derivation of the transition probability formulae), the leading order dependence of P (ν µ → ν e ) and P (ν µ → ν µ ) on the NSI parameters are reported in table 2. Here we use the shorthand notation αβ = , with the meaning O(1) to identify a common order in the perturbative expansion of the probabilities. Thus, a change of systematics related to the ν µ → ν e transition most probably will affect the determination of eµ and eτ , with scarse effects on the other parameters µµ and µτ , while a variation in the shape and absolute JHEP08 (2018) Table 2. Leading order dependence of the oscillation probabilities P (ν µ → ν e ) and P (ν µ → ν µ ) on the NSI parameters.
normalizations of the ν µ signal will reflect on a different sensitivity to µµ and µτ . For ee the situation is a bit different; in fact, even though it is perturbatively suppressed at O( 3 ), nevertheless the current bounds allow it to assume O(1) values, thus making possible a significant dependence on the systematic choice.

Experimental setup and treatment of the systematics
In this paper we focus on the DUNE experiment as described in [10]; we consider a 40 Kton liquid argon detector at a baseline L = 1300 Km and on-axis with respect to the beam direction. The neutrino fluxes correspond to the Optimized Desing flux, from a proton beam energy of 80 GeV and a beam power of 1.07 MW. 3.5 years of data taking are assumed in both neutrino and antineutrino modes, for a total exposure of 300 kt·MW·year for both ν e appearance and ν µ disappearance. No near detector is assumed in our numerical simulations. At the far detector the ν µ disappearance sample is composed of ν µ CC interactions, with main backgrounds from neutral current (NC) interactions where charged pions are misidentified as muons, and ν τ +ν τ CC interactions in which the produced τ s decay to a muon and two neutrinos. The ν e appearance sample is composed of ν e CC interactions from ν µ → ν e oscillation and the relevant backgrounds originate from intrinsic ν e +ν e beam contamination, NC and ν µ +ν µ CC interactions in which a photon from electromagnetic neutral pion decays is misidentified as an electron, and from ν τ +ν τ interactions in which the resulting τ decays to an electron. It should be noticed that in both neutrino and antineutrino running modes we consider the sum ν µ → ν e ⊕ν µ →ν e as signal appearance events because the information coming from the wrong-sign events can be of some relevance. According to [22], the same approach is also used for the disappearance channels, for which the sum ν µ → ν µ ⊕ν µ →ν µ is considered.

JHEP08(2018)028
The Liquid Argon Time Projection Chamber (LArTPC) performance parameters used to get the sensitivity plots discussed in this paper have been generated using the DUNE Fast Monte Carlo simulation [23] and translated into GLoBES files [24,25] made publicly available as ancillary files in [22]. Beside cross section files describing neutrino charged and neutral current interactions in Argon (generated using GENIE 2.8.4), the DUNE Collaboration also provided true-to-reconstructed smearing matrices (not modified in this work) as well as selection efficiencies as a function of bin energies for the various signals and backgrounds.
Our implementation of the χ 2 and the method of determining the confidence regions is based on the pull method [7,26,27] and represents the standard implementation of systematic uncertainties in GLoBES. The χ 2 is obtained after the minimization over the nuisance parameters ξ [7]. For each transition channel c and energy bin i we use a Poissonian χ 2 of the form: is the predicted number of events in the i-th energy bin for a given channel c, for a set of oscillation parameters θ and nuisance parameters ξ. O c,i , instead, is the observed event rate, that is the one corresponding to assumed true values of the oscillation parameters. Both F c,i and O c,i receive contributions from different sources s, that tipically include signal and background rates specified by R c,s,i ( θ), so that for example: The auxiliary parameters a c,s have the form a c,s ≡ k w c,s,k ξ k , in which the coefficients w c,s,k can assume the values one or zero depending on whether a particular nuisance parameter ξ k affects the contribution from the source s to channel c or not, respectively. Thus, the total χ 2 is given by: where the last two contributions are the pull terms associated with the shape and the overall signal normalization, respectively. The ξ φ, i shape parameters are bin-to-bin uncorrelated whereas the normalization parameter ξ N is fully correlated between different energy bins i. In [22] only overall normalizations for signal and backgrounds have been taken into account; for ν e andν e signal modes they are fixed to 2% each, while a 5% is assumed for ν µ andν µ signals. For the backgrounds the normalization uncertainties range from 5% (for intrinsic ν e and misidentified ν µ ) to 20% (for misidentified τ ) with a 10% assumed for the neutral currents. While the background normalization uncertainties are not changed in our numerical simulations, we consider several options for the signal normalization uncertainties and, in addition, we also take into account various level of normalizations. In order to maintain the possible combinations of shape and normalization errors to a reasonable level, we decided to split the systematic uncertainties in two different classes: those in which only the ν e signal uncertainties are changed (and the ones related to µ's are fixed, case-I) and, JHEP08(2018)028 viceversa, those in which only the ν µ signal uncertainties are changed (case-II). Thus the following combinations for case-I are considered: ν e app = 2.5%, 5% shape normalization : ν µ dis = 5% ν e app = 2%, 7% , which give rise to four different scenarios: • A, very optimistic, where the pair (shape, absolute) = (2%, 2.5%); • B 1 , (shape, absolute) = (2%, 5%); • B 2 , (shape, absolute) = (7%, 2.5%); • C, very pessimistic, with (shape, absolute) = (7%, 5%).
It is worth mentioning that systematic uncertainties between 2% and 5% have also been assumed in [10] whereas the pessimistic 7% has been introduced by us (but also taken into account in [11]) to contemplate more conservative estimates on the physics reach of DUNE. For the sake of simplicity, while retaining the nominal binning for the spectrum, we adopt wider intervals for the systematic uncertainties: [0  Table 3. Central values and relative uncertainties of the standard mixing parameters extracted from [1]. As in [10], for non-Gaussian parameters the relative uncertainty is computed using 1/6 of the 3σ allowed range.

Numerical results
In this section we discuss in details how the bounds on the NSI parameters change if different assumptions on systematics are made according to the discussion of the previous section. In deriving the CL intervals for a given αβ , we marginalize over all other standard and non-standard parameters: for the NSI's, the phases are left free in the whole [0, 2π) interval while the moduli are marginalized taking into account the 90% CL ranges quoted in All numerical results of this section have been obtained using a modified version of GLoBES which includes non-standard interactions affecting the propagation processes for active neutrinos [29]. The event rates for signal and background in both neutrino and antineutrino modes for αβ = 0 are reported in

Results for case-I
According to the discussion after table 1, we start considering the bounds at the 90% CL (1 degree of freedom) that DUNE can set on the parameters | eµ | (left panel of figure 1) 2 Since the main goal of the paper is to study the effects of the systematics on the study of the NSI scenario more than deriving new bounds on the NSI parameters, we adopted here the simple approach to take external priors on the standard mixing parameters. 3 As shown in [28], bounds obtained with all the parameters free to vary are not importantly affected by   • it is clear that the behavior of the χ 2 function is quite similar for the pairs (A, B 1 ) and (B 2 , C), signalizing that the shape uncertainty (which is the same within the pairs) is by far more relevant than the overall signal normalization. In fact, even doubling it from 2.5% (cases A and B 2 ) to 5% (cases B 1 and C) we hardly see any significant difference. In fact, this has to be expected since the disappearance distribution fixes the normalisation of both muon neutrino and antineutrino event samples; JHEP08(2018)028  • in the most favorable cases (A, B 1 ), DUNE can set a bound on | eµ | at the 90% CL around b eµ = 0.130 while in the worst cases (B 2 , C) this is deteriorated to b eµ = 0.135, thus a ∼4% less stringent bound (not a significant difference given the simple treatment of the systematics adopted in this paper); • for | eτ | the bounds at the 90% CL are b eτ ∼ 0.32 and b eτ ∼ 0.33 for the cases (A, B 1 ) and (B 2 , C), respectively, which entails a worsening by ∼3%.
We have also checked that a further increase of the absolute normalization systematics to 7% does not produce relevant changes in the obtainable bounds on | eµ | and | eτ |. As already anticipated, for the other parameters µµ and | µτ | there is no important dependence on the assumption on the systematic errors, as we can see in the central and right panels of figure 2. However, some dependence on systematics is visible for ee (left panel of figure 2): considering the extreme cases A and C, the allowed 90% CL interval for case A is about 7% smaller than for case C.
For the sake of completeness, we summarize in table 5 the 90% CL bounds on the five real parts of the αβ parameters. 5 Comparing with the limits quoted in table 1, we observe that a strong improvement can be set on the bounds of | eµ |, which passes from b eµ < 0.36 to b eµ 0.13 for almost every choice of the systematics. For | eτ |, instead, in all cases a reduction of the upper bound by roughly ∼40% can be obtained. As expected, the JHEP08(2018)028   other off-diagonal parameters whose dependence in the oscillation probability P (ν µ → ν e ) is subleading, do not benefit of any improvements of the systematics, even in the more aggressive case A. A separate discussion should be devoted to ee ; in fact, while the current limit consists of two separate islands above and below ee = 0, DUNE alone will be able to exclude negative ee but not those positive values within the current bounds [17].

Results for case-II
For the case-II, the dependence of P (ν µ → ν µ ) on the NSI parameters suggests that a change in the systematics affects more µµ and µτ than the others; this is clearly visible in figure 3, where we reported the χ 2 function for µµ (left panel) and | µτ | (right panel), and in figure 4 where ee (left panel), | eµ | (central panel) and | eτ | (right panel) have been considered. The trend of having much better bounds in the cases (A µ , B 1µ ) than (B 1µ , C µ ) is confirmed by our numerical results although, if compared to the case-I, there is a more clear separation among the curves in the same pairs. This is due to the fact that the muon disappearance is now the channel under discussion, which provides a much better statistics than the appearance one. Going into details, we observe that the whole allowed range of µµ improves by ∼26% from case C µ to A µ while that of | µτ | by roughly 15%, see figure 3, thus making extremely relevant a reduction of the shape uncertainty for the ν µ signal. No significant changes (below 3%) are observed for the other moduli, see figure 4, not even for ee for which the 3 -dependence in P (ν µ → ν µ ) cannot be compensated and made comparable with the -dependence of µµ and µτ by O(1) moduli.

JHEP08(2018)028
In table 6 we report the obtained 90% CL bounds also for the case-II. The most relevant feature is that, like for case-I, more stringent bounds by a factor of ∼3 can be set on | eµ | and, to a lesser extent, on | eτ |, while it is clear that DUNE cannot improve the current limits for all other parameters, even in the most aggressive scenarios for the systematic uncertainties.

Effects of NSI parameters on the extraction of standard parameters
It is interesting to estimate how the presence of matter NSI can worsen the sensitivity to the octant of θ 23 , to mass ordering and to CP violation at DUNE. To make things simpler, we limit ourselves to two extreme scenarios, built from the two cases analyzed above: the (A, A µ ) scenario (OPTimistic), which means that the shape normalization error is 2% and the absolute normalization error is 2.5% for both ν e and ν µ signals, and the (C, C µ ) scenario (PESsimistic) in which the shape normalization error is 7% and the absolute normalization error is 5%, again for ν e and ν µ signals. In this way we are confident that our results will include a large class of intermediate assumptions on systematic errors. In both situations, the errors and correlations among the several sources of backgrounds are the same as above. Unless stated otherwise, in all plots were NSI parameters are taken into account, they have been marginalized over assuming the bounds reported in table 1.
Determination of the octant of θ 23 . The uncertainty in the determination of whether θ 23 is larger or smaller than maximal mixing stems from the fact that its measurements are mainly due to the ν µ disappearance channel which depends on sin 2 2θ 23 . In this respect, the JHEP08(2018)028 fact the DUNE can have simultaneous access to both ν µ disappearance and ν e appearance channels provides an useful synergy to probe the octant hypothesis. To estimate the sensitivity, we adopt the following metric [31]: without imposing any priors on the atmospheric angle. Our results (for NO only) are reported in figure 5. In the right panel we show the situation when all NSI parameters are fixed and set to zero, thus this corresponds to pure SM (all NSI vanishing) results. The dashed horizontal lines indicate the 3σ and 5σ CL; the area in the light green region represents the range in sensitivity due to potential variations in the true value of δ CP in the case of the PES scenario while the region in blue is obtained for the OPT case. We clearly see that DUNE can determine the octant at high confidence level as soon as θ 23 is away from maximal mixing by ∼ 3 • for the OPT case and by ∼ 4 • for the PES case; to be more precise, in the most favorable cases, a 5σ discovery of the octant would be guaranteed outside the intervals θ 23 ∈ [42.0 • , 49.1 • ] for the OPT assumption and θ 23 ∈ [40.7 • , 50.2 • ] for PES. This also means that for the θ 23 central value reported in table 1 the significance with which DUNE can resolve the octant is scarse [32]. Notice also that the difference between the PES and OPT cases is more evident for smaller values of sin 2 2θ 23 , where ∆χ 2 can be as large as four units. The situation gets sensibly worse when marginalization over the NSI parameters is taken into account, as we can see in the right panel of figure 5. In JHEP08(2018)028 particular, the discovery of the octant drops below 2σ in the PES case even for values of θ 23 at the extremes of the range analyzed here while a 2σ sensitivity can only be reached in the OPT case and for mixing angles outside the interval θ 23 ∈ [38 • , 52 • ], as one would have expected [33], and for favorable values of the standard CP phase.
CP violation sensitivity. In the standard framework of three neutrino oscillation, a signal indicating violation of the CP symmetry in the lepton sector will be observable if the true value of δ CP is sufficiently different from the CP conserving values δ CP = 0 and δ CP = π. In order to estimate the capability of DUNE to determine leptonic CP violation, we make use of the following indicator: . Given the fact that the smallest the value of sin 2 θ 23 the largest the sensitivity to CP violation is, we generally expect a smaller significance with respect to the results quoted in [10] even in the case of the SM only, as our sin 2 θ 23 ∼ 0.59 while in [10] they used sin 2 θ 23 = 0.45 (and similar values for the other mixing angles). In the fit procedure we marginalized over all undisplayed parameters; the true value of θ 23 is set in the second octant, according to the central values quoted in table 3.
Our results are reported in figure 6 where we displayed the CP discovery potential as a function of the true value of the leptonic CP phase, for both NO (left panel) and IO (right panel). The bands represents the range in sensitivity obtained under the two different assumptions for the systematics, with the implicit meaning that the upper curves in each bands correspond to the OPT case while the lower ones to the PES case. The results for the SM are in blue/solid line, those considering the marginalization over the NSI parameters in red/dashed.
In the SM case, the DUNE setup adopted in this paper is enough to guarantee a 5σ discovery potential for the leptonic phase only for the inverted ordering of the neutrino mass eigenstates, in ∼31% of the cases in the OPT scenario; if PES systematics are considered, only a 3σ discovery potential can be reached in ∼50% of the true δ CP . For NO and for the OPT case, an encouraging 4σ can be reached in the ∼39% of the δ CP values, roughly the same percentage obtainable in the PES case at 3σ. The marginalization over the NSI parameters strongly reduces the CP discovery potential [34] such that not even 3σ discovery can be claimed.
In this context the dependence on the choice of systematics is quite relevant for both SM and NSI-marginalized cases: the sensitivities for the OPT and PES assumptions differ by ∼20% and ∼15% for |δ CP | ∼ π/2, respectively, almost independently on the considered hierarchy.
Sensitivity to mass hierarchy. The long-baseline and high neutrino energies of the DUNE flux are particularly suited to explore the sensitivity to mass hierarchy, that is to determine the sign of the atmospheric mass difference ∆m 2 31 . As it is well known [35]- [36], the marginalization over all NSI parameters produces a loss in the sensitivity to mass hierarchy (both for NO and IO); in our numerical simulations, this barely reaches values JHEP08(2018)028 larger than ∆χ 2 ∼ 1.4 when δ CP ∼ π/2. This conclusion can be traced back to the intrinsic degeneracy involving ee and δ CP . Since the primary goal of this paper is to illustrate the effects of systematics, we decided to marginalize ee in the positive interval reported in table 1, with a true value ee = 0.7 and to keep fixed all other NSI parameters but eτ , whose correlation with ee has been shown to be important in decreasing the sensitivity of long-baseline experiments to several physical quantities [37]. We quantify the experimental sensitivity to the mass hierarchy using: 6 true normal ordering , Our results are reported in figure 7 where we show ∆χ 2 M H as a function of the true value of the CP phase δ CP , for the two cases where the true ordering is the NO (left panel) or the IO one (right panel). The plots in the boxes depict the same variable ∆χ 2 M H computed in the SM. The bands represent the range in sensitivity due to the different assumptions on the systematics, with the meaning that the largest ∆χ 2 M H corresponds to the OPT case. Within our setups, the mass ordering in the SM case can be determined at ∆χ 2 M H = 5 for every value of the true δ CP for both NO and IO while the same is not true when the NSI marginalization is taken into account: in this case, the mass hierarchy discovery potential 6 As shown in refs. [38,39], the ∆χ 2 metric does not follow the χ 2 function; it is used here as a representative of the mean or the most likely value of the true ∆χ 2 that would have been obtained in an ensemble of experiments. is larger than ∆χ 2 M H = 3 for every δ CP and occasionally also larger than 4. While in the standard case (no NSI) the systematics change the sensitivity by 20 to 28% (for both assumptions on the true neutrino mass ordering), when NSI are taken into account the hypothesis on the systematic effects change the results by a ∼12% on average.

Conclusions
In this paper we tried to characterize the role of different assumptions on shape and normalization uncertainties for the ν µ and ν e signals at the DUNE facilities in the determination of the bounds on the NSI parameters. The main results of our study can be summarized in the following points: • shape uncertainty is much more relevant than the signal normalization uncertainty in setting the bounds on NSI parameters; • for | eµ |, reasonable assumptions on the systematics (for both appearance and disappearance channels) are enough to lower by a factor of ∼3 the existing upper limits; • on the other hand, for | eτ | only a 30-40% reduction of the current upper limits can be reached in DUNE, depending on whether the assumptions of case-I or II are taken into account, respectively; • µµ and | µτ | strongly benefit of smaller systematic uncertainties in the ν µ sector but even in the most favorable cases, A and A µ , the current limits cannot be ameliorated; • for ee , the configuration of DUNE used in this paper helps in disfavoring negative values, independently on the assumptions on systematics.

JHEP08(2018)028
In addition to the previous considerations, we have studied the impact of the marginalization on the NSI parameters on the sensitivity to the octant of θ 23 , to mass ordering and to CP violation at DUNE. In the case of the octant of θ 23 , a 2σ sensitivity in the NSI-marginalized case can be reached for very optimistic assumptions on the systematics and for mixing angles outside the interval θ 23 ∈ [38 • , 52 • ] and for favorable values of the standard CP phase. For the sensitivity to CP violation, the choice of systematics can alter the results at |δ CP | ∼ π/2 by 15-20%, where the largest impact is seen in the SM case; in this respect not much difference has been found for the different hierarchies. Finally, we have observed that the choice of the systematics is not crucial to get a sensitivity ∆χ 2 M H above 5σ in the SM case and well above 3σ in the NSI-marginalized case for every value of δ CP , even though the effects of changing from the PES to the OPT assumptions are well above 20% and 10%, respectively.