Long-baseline neutrino oscillation physics potential of the DUNE experiment

The sensitivity of the Deep Underground Neutrino Experiment (DUNE) to neutrino oscillation is determined, based on a full simulation, reconstruction, and event selection of the far detector and a full simulation and parameterized analysis of the near detector. Detailed uncertainties due to the flux prediction, neutrino interaction model, and detector effects are included. DUNE will resolve the neutrino mass ordering to a precision of 5$\sigma$, for all $\delta_{\mathrm{CP}}$ values, after 2 years of running with the nominal detector design and beam configuration. It has the potential to observe charge-parity violation in the neutrino sector to a precision of 3$\sigma$ (5$\sigma$) after an exposure of 5 (10) years, for 50\% of all $\delta_{\mathrm{CP}}$ values. It will also make precise measurements of other parameters governing long-baseline neutrino oscillation, and after an exposure of 15 years will achieve a similar sensitivity to $\sin^{2} 2\theta_{13}$ to current reactor experiments.

Abstract The sensitivity of the Deep Underground Neutrino Experiment (DUNE) to neutrino oscillation is determined, based on a full simulation, reconstruction, and event selection of the far detector and a full simulation and parameterized analysis of the near detector. Detailed uncertainties due to the flux prediction, neutrino interaction model, and detector effects are included. DUNE will resolve the neutrino mass ordering to a precision of 5σ, for all δ CP values, after 2 years of running with the nominal detector design and beam configuration. It has the potential to observe chargeparity violation in the neutrino sector to a precision of 3σ (5σ) after an exposure of 5 (10) years, for 50% of all δ CP values. It will also make precise measurements of other parameters governing long-baseline neutrino oscillation, and after an exposure of 15 years will achieve a similar sensitivity to sin 2 2θ 13 to current reactor experiments.

Introduction
The Deep Underground Neutrino Experiment (DUNE) is a next-generation, long-baseline neutrino oscillation experiment which will carry out a detailed study of neutrino mixing utilizing high-intensity ν µ andν µ beams measured over a long baseline. DUNE is designed to make significant contributions to the completion of the standard three-flavor picture by measuring all the parameters governing ν 1 -ν 3 and ν 2 -ν 3 mixing in a single experiment. Its main scientific goals are the definitive determination of the neutrino mass ordering, the definitive observation of charge-parity symmetry violation (CPV) for more than 50% of possible true values of the charge-parity violating phase, δ CP , and precise measurement of oscillation parameters, particularly δ CP , sin 2 2θ 13 , and the octant of θ 23 . These measurements will help guide theory in understanding if there are new symmetries in the neutrino sector and whether there is a relationship between the generational structure of quarks and leptons [1]. Observation of CPV in neutrinos would be an important step in understanding the origin of the baryon asymmetry of the universe [2,3].
The DUNE experiment will observe neutrinos from a high-power neutrino beam peaked at ∼2.5 GeV but with a broad range of neutrino energies, a near detector (ND) located at Fermi National Accelerator Laboratory, in Batavia, Illinois, USA, and a large liquid argon time-projection chamber (LArTPC) far detector (FD) located at the 4850 ft level of Sanford Underground Research Facility (SURF), in Lead, South Dakota, USA, 1285 km from the neutrino production point. The neutrino beam provided by Long-Baseline Neutrino Facility (LBNF) [4] is produced using protons from Fermilab's Main Injector, which are guided onto a graphite target, and a traditional horn-focusing system to select and focus particles produced in the target [5]. The polarity of the focusing magnets can be reversed to produce a beam dominated by either muon neutrinos or muon antineutrinos. A highly capable ND will constrain many systematic uncertainties for the oscillation analysis. The 40-kt (fiducial) FD is composed of four 10 kt (fiducial) LArTPC modules [6][7][8]. The deep underground location of the FD reduces cosmogenic and atmospheric sources of background, which also provides sensitivity to nucleon decay and low-energy neutrino detection, for example, the possible observation of neutrinos from a core-collapse supernova [5].
The entire complement of neutrino oscillation experiments to date has measured five of the neutrino mixing parameters [9][10][11]: the three mixing angles θ 12 , θ 23 , and θ 13 , and the two squared-mass differences ∆m 2 21 and |∆m 2 31 |, where ∆m 2 ij = m 2 i − m 2 j is the difference between the squares of the neutrino mass states in eV 2 . The neutrino mass ordering (i.e., the sign of ∆m 2 31 ) is unknown, though recent results show a weak preference for the normal ordering [12][13][14]. The value of δ CP is not well known, though neutrino oscillation data are beginning to provide some information on its value [12,15].
The oscillation probability of ν µ → ν e through matter in the standard three-flavor model and a constant density approximation is, to first order [16]: where a = ± G F N e √ 2 ≈ ± 1 3500 km ρ 3.0 g/cm 3 , G F is the Fermi constant, N e is the number density of electrons in the Earth's crust, ∆ ij = 1.267∆m 2 ij L/E ν , L is the baseline in km, and E ν is the neutrino energy in GeV. Both δ CP and a terms are positive for ν µ → ν e and negative forν µ →ν e oscillations; i.e., a neutrino-antineutrino asymmetry is introduced both by CPV (δ CP ) and the matter effect (a). The origin of the matter effect asymmetry is simply the presence of electrons and absence of positrons in the Earth [17,18].
The (anti-)electron neutrino appearance probability is shown in Figure 1 at the DUNE baseline of 1285 km as a function of neutrino energy for several values of δ CP . Fig. 1 The appearance probability at a baseline of 1285 km, as a function of neutrino energy, for δ CP = −π/2 (blue), 0 (red), and π/2 (green), for neutrinos (top) and antineutrinos (bottom), for normal ordering.
DUNE has a number of features that give it unique physics reach, complementary to other existing and planned experiments [19][20][21]. Its broad-band beam makes it sensitive to the shape of the oscillation spectrum for a range of neutrino energies. DUNE's relatively high energy neutrino beam enhances the size of the matter effect and will allow DUNE to measure δ CP and the mass ordering simultaneously. The unique LArTPC detector technology will enhance the resolution on DUNE's measurement of the value of δ CP , and along with the increased neutrino energy, gives DUNE a different set of systematic uncertainties to other experiments, making DUNE complementary with them.
This paper describes studies that quantify DUNE's expected sensitivity to long-baseline neutrino oscillation, using the accelerator neutrino beam. Note that atmospheric neutrino samples would provide additional sensitivity to some of the same physics, but are not included in this work. The flux simulation and associated uncertainties are described in Section 2. Section 3 describes the neutrino interaction model and systematic variations. The near and far detector simulation, reconstruction, and event selections are described in Sections 4 and 5, respectively, with a nominal set of event rate predictions given in Section 6. Detector uncertainties are described in Section 7. The methods used to extract oscillation sensitivities are described in Section 8. The primary sensitivity results are presented in Section 9. We present our conclusions in Section 10.

Neutrino Beam Flux and Uncertainties
The expected neutrino flux is generated using G4LBNF [5,22], a Geant4-based [23] simulation of the LBNF neutrino beam. The simulation uses a detailed description of the LBNF optimized beam design [5], which includes a target and horns designed to maximize sensitivity to CPV given the physical constraints on the beamline design.
Neutrino fluxes for neutrino-enhanced, forward horn current (FHC), and antineutrino-enhanced, reverse horn current (RHC), configurations of LBNF are shown in Figure 2. Uncertainties on the neutrino fluxes arise primarily from uncertainties in hadrons produced off the target and uncertainties in the design parameters of the beamline, such as horn currents and horn and target positioning (commonly called "focusing uncertainties") [5]. Given current measurements of hadron production and LBNF estimates of alignment tolerances, flux uncertainties are approximately 8% at the first oscillation maximum and 12% at the second. These uncertainties are highly correlated across energy bins and neutrino flavors. The unoscillated fluxes at the ND and FD are similar, but not identical. The relationship is well understood, and flux uncertainties mostly cancel for the ratio of fluxes between the two detectors. Uncertainties on the ratio are dominated by focusing uncertainties and are ∼1% or smaller except at the falling edge of the focusing peak (∼4 GeV), where they rise to 2%. The rise is due to the presence of many particles which are not strongly focused by the horns in this energy region, which are particularly sensitive to focusing and alignment uncertainties. The near-to-far flux ratio and uncertainties on this ratio are shown in Fig. 3.
Beam-focusing and hadron-production uncertainties on the flux prediction are evaluated by reproducing the full beamline simulation many times with variations of the input model according to those uncertainties. The resultant uncertainty on the neutrino flux prediction is described through a covariance matrix, where each bin corresponds to an energy range of a particular beam mode and neutrino species, separated by flux at the ND and FD. The output covariance matrix has 208 × 208 bins, despite having only ∼30 input uncertainties. To reduce the number of parameters used in the fit, the covariance matrix is diagonalized, and each principal component is treated as an uncorrelated nuisance parameter. The 208 principal components are ordered by the magnitude of their corresponding eigenvalues, which is the variance along the principal component (eigenvector) direction, and only the first ∼30 are large enough that they need to be included. This was validated by including more flux parameters and checking that there was no significant change to the sensitivity for a small number of test cases. By the 10th principal component, the eigenvalue is 1% of the largest eigenvalue. As may be expected, the largest uncertainties correspond to the largest principal components as shown in Figure 4. The largest principal component (component 0) matches the hadron production uncertainty on nucleonnucleus interactions in a phase space region not covered by data. Components 3 and 7 correspond to the dataconstrained uncertainty on proton interactions in the target producing pions and kaons, respectively. Components 5 and 11 correspond to two of the largest focusing uncertainties, the density of the target and the horn current, respectively. Other components not shown either do not fit a single uncertain parameter or may represent two or more degenerate systematics or ones that produce anti-correlations in neighboring energy bins.  Fig. 4 Select flux principal components are compared to specific underlying uncertainties from the hadron production and beam focusing models. Note that while these are shown as positive shifts, the absolute sign is arbitrary.
Future hadron production measurements are expected to improve the quality of, and the resulting constraints on, these flux uncertainty estimates. Approximately 40% of the interactions that produce neutrinos in the LBNF beam simulation have no direct data constraints. Large uncertainties are assumed for these interactions. The largest unconstrained sources of uncertainty are proton quasielastic interactions and pion and kaon rescattering in beamline materials. The proposed EMPHATIC experiment [24] at Fermilab will be able to constrain quasielastic and low-energy interactions that dominate the lowest neutrino energy bins. The NA61 experiment at CERN has taken data that will constrain many higher energy interactions, and also plans to measure hadrons produced on a replica LBNF target, which would provide tight constraints on all interactions occurring in the target. A similar program at NA61 has reduced flux uncertainties for the T2K experiment from ∼10% to ∼5% [25]. Another proposed experiment, the LBNF spectrometer [26], would measure hadrons after both production and focusing in the horns to further constrain the hadron production uncertainties, and could also be used to experimentally assess the impact of shifted alignment parameters on the focused hadrons (rather than relying solely on simulation).

Neutrino interaction model and uncertainties
A framework for considering the impact of neutrino interaction model uncertainties on the oscillation analysis has been developed. The default interaction model is implemented in v2.12.10 of the GENIE generator [27,28]. Variations in the cross sections are implemented in various ways: using GENIE reweighting parameters (sometimes referred to as "GENIE knobs"); with ad hoc weights of events that are designed to parameterize uncertainties or cross-section corrections currently not implemented within GENIE; or through discrete alternative model comparisons. The latter are achieved through alternative generators, alternative GENIE configurations, or custom weightings, which made extensive use of the NUISANCE package [29] in their development.
The interaction model components and uncertainties can be divided into seven groups: (1) initial state, (2) hard scattering and nuclear modifications to the quasielastic, or one-particle one-hole (1p1h) process, (3) multinucleon, or two-particle two-hole (2p2h), hard scattering processes, (4) hard scattering in pion production processes, (5) higher invariant mass (W ) and neutral current (NC) processes, (6) final-state interactions (FSI), (7) neutrino flavor dependent differences. Uncertainties are intended to reflect current theoretical freedom, deficiencies in implementation, and/or current experimental knowledge. The default nuclear model in GENIE describing the initial state of nucleons in the nucleus is the Bodek-Ritchie global Fermi gas model [30]. There are significant deficiencies that are known in global Fermi gas models: these include a lack of consistent incorporation of the high-momentum tails in the nucleon momentum distribution that result from correlations among nucleons; the lack of correlation between location within the nucleus and momentum of the nucleon; and an incorrect relationship between momentum and energy of the off-shell, bound nucleon within the nucleus. They have also been shown to agree poorly with neutrinonucleus scattering data [31]. GENIE modifies the nucleon momentum distribution empirically to account for short-range correlation effects, which populates the high-momentum tail above the Fermi cutoff, but the other deficiencies persist. Alternative initial state models, such as spectral functions [32,33], the mean field model of GiBUU [34], or continuum random phase approximation (CRPA) calculations [35] may provide better descriptions of the nuclear initial state [36], but are not considered further here.
The primary uncertainties considered in 1p1h interactions (ν l + n → l − + p,ν l + p → l + + n) are the axial form factor of the nucleon and the nuclear screeningfrom the so-called random phase approximation (RPA) calculations-of low momentum transfer reactions. The Valencia group's [37,38] description of RPA comes from summation of W ± self-energy terms. In practice, this modifies the 1p1h (quasielastic) cross section in a nontrivial way, with associated uncertainties presented in Ref. [39], which were evaluated as a function of Q 2 . Here we use T2K's 2017/8 parameterization of the Valencia RPA effect [12]. The shape of the correction and error is parameterized with a third-order Bernstein polynomial up to Q 2 = 1.2 GeV 2 where the form transitions to a decaying exponential. The BeRPA (Bernstein RPA) function has three parameters controlling the behavior at increasing Q 2 (A, B and D), a fourth parameter (E) that controls the high-Q 2 tail, and a fifth (U), which changes the position at which the behaviour changes from polynomial to exponential. The BeRPA parameterization modifies the central value of the model prediction, as decribed in Table 3. BeRPA parameters E and U are not varied in the analysis described here, the parameters A and B have a prefit uncertainty of 20%, and D has a prefit unertainty of 15%. The axial form factor parameterization we use, a dipole, is known to be inadequate [40]. However, the convolution of BeRPA uncertainties with the limited axial form factor uncertainties do provide more freedom as a function of Q 2 , and the two effects combined likely provide adequate freedom for the Q 2 shape in quasielastic events. BBBA05 vector form factors are used [41].
The 2p2h contribution to the cross section comes from the Valencia model [37,38], the implementation in GENIE is described in Ref. [42]. However, MIN-ERvA [43] and NOvA [44] have shown that this model underpredicts observed event rates on carbon. The extra strength from the "MINERvA tune" to 2p2h is applied as a two-dimensional Gaussian in (q 0 , q 3 ) space, where q 0 is the energy transfer from the leptonic system, and q 3 is the magnitude of the three momentum transfer) to fit reconstructed MINERvA CC-inclusive data [43]. Reasonable predictions of MINERvA's data are found by attributing the missing strength to any of 2p2h from np initial state pairs, 2p2h from nn initial state pairs, or 1p1h (quasielastic) processes. The default tune uses an enhancement of the np and nn initial strengths in the ratio predicted by the Valencia model, and alternative systematic variation tunes ("MnvTune" 1-3) attribute the missing strength to the individual interaction processes above. We add uncertainties for the energy dependence of this missing strength based on the MINERvA results [43], and assume a generic form for the energy dependence of the cross section using the "A" and "B" terms taken from Ref. [45]. These uncertainties are labeled E 2p2h and are separated for neutrinos and antineutrinos. We add uncertainties on scaling the 2p2h prediction from carbon to argon on electron-scattering measurements of short-range correlated (SRC) pairs taken on multiple targets [46], separately for neutrinos (ArC2p2h ν) and antineutrinos (ArC2p2hν).
GENIE uses a modified version of the Rein-Sehgal (R-S) model for pion production [47], including only the 16 resonances recommended by the Particle Data group [48], and excluding interferences between resonances. The cross section is cut off at invariant masses, W ≥ 1.7 GeV (2 GeV in the original R-S model). No inmedium modifications to the resonances are included, and by default they decay isotropically in their rest frame, although there is a parameter denoted here as "θ π from ∆-decay", for changing the angular distribution of pions produced through ∆ resonance decays to match the experimentally observed distributions used in the original R-S paper [47]. Resonance decays to η and γ (plus a nucleon) are included from Ref. [48]. We use a tuning of the GENIE model to reanalyzed neutrino-deuterium bubble chamber data [49,50] as our base model, as noted in Table 3. We note that an improved Rein-Sehgal-like resonance model has been developed [51], and has been implemented in Monte Carlo generators, although is not used as the default model in the present work.
The deep inelastic scattering (DIS) model implemented in GENIE uses the Bodek-Yang parametrization [52], using GRV98 parton distribution functions [53]. Hadronization is described by the AKGY model [54], which uses the KNO scaling model [55] for invariant masses W ≤ 2.3 GeV and PYTHIA6 [56] for invariant masses W ≥ 3 GeV, with a smooth transition between the two for intermediate invariant masses. A number of variable parameters affecting DIS processes are included in GENIE, as listed in Table 3, and described in Ref. [52]. In GENIE, the DIS model is extrapolated to all values of invariant mass, and replaces the non-resonant background to pion production in the R-S model.
The NOvA experiment [57] developed uncertainties beyond those provided by GENIE to describe their single pion to DIS transition region data. We follow their findings, and implement separate, uncorrelated uncertainties for all perturbations of 1, 2, and ≥ 3 pion final states, CC/NC, neutrinos/antineutrinos, and interactions on protons/neutrons, with the exception of CC neutrino 1-pion production, where interactions on protons and neutrons are merged, following [50], which modifies the central value of the model prediction, as listed in Table 3. This leads to 23 distinct uncertainty channels with a label to denote the process it affects: NR [ν,ν] [CC,NC] [n,p] [1π,2π,3π]. Each channel has an uncertainty of 50% for W ≤ 3 GeV, and an uncertainty which drops linearly above W = 3 GeV until it reaches a flat value of 5% at W = 5 GeV, where external measurements better constrain this process.
GENIE includes a large number of final state uncertainties on its final state cascade model [58][59][60], which are summarized in Table 2. A recent comparison of the underlying interaction probabilities used by GENIE is compared with other available simulation packages in Ref. [61].
The cross sections include terms proportional to the lepton mass, which are significant contributors at low energies where quasielastic processes dominate. Some of the form factors in these terms have significant uncertainties in the nuclear environment. Ref. [62] ascribes the largest possible effect to the presence of poorly constrained second-class current vector form factors in the nuclear environment, and proposes a variation in the cross section ratio of σ µ /σ e of ±0.01/Max(0.2 GeV, E ν ) for neutrinos and ∓0.018/Max(0.2 GeV, E ν ) for antineutrinos. Note the anticorrelation of the effect in neutrinos and antineutrinos. This parameter is labeled ν e /ν e norm.
An additional normalization uncertainty (NC norm.) of 20% is applied to all NC events at the ND in this analysis to investigate whether the small contamination of NC events which passed the simple selection cuts had an effect on the analysis. Although a similar systematic could have been included (uncorrelated) at the FD, it was not in this analysis.
Finally, some electron-neutrino interactions occur at four-momentum transfers where a corresponding muonneutrino interaction is kinematically forbidden, therefore the nuclear response has not been constrained by muon-neutrino cross-section measurements. This region at lower neutrino energies has a significant overlap with the Bodek-Ritchie tail of the nucleon momentum distribution in the Fermi gas model [30]. There are signif- icant uncertainties in this region, both from the form of the tail itself and from the lack of knowledge about the effect of RPA and 2p2h in this region. Here, a 100% uncertainty is applied in the phase space present for ν e but absent for ν µ (labeled ν e phase space (PS)).
The complete set of interaction model uncertainties includes GENIE implemented uncertainties (Tables 1  and 2), and new uncertainties developed for this effort (Table 4) which represent uncertainties beyond those implemented in the GENIE generator.
Tunes which are applied to the default model, using the dials described, which represent known deficiencies in GENIE's description of neutrino data, are listed in Table 3.
The way model parameters are treated in the analysis is described by three categories: -Category 1: expected to be constrained with on-axis data; uncertainties are implemented in the same way for ND and FD. N. PROD, Nucleon π-production probability ±20% π CEX, π charge exchange probability ±50% π EL, π elastic reaction probability ±10% π INEL, π inelastic reaction probability ±40% π ABS, π absorption probability ±20% π PROD, π π-production probability ±20% Table 2 The intra-nuclear hadron transport systematic parameters implemented in GENIE with associated uncertainties considered in this work. Note that the 'mean free path' parameters are omitted for both N-N and π-N interactions as they produced unphysical variations in observable analysis variables.  Table 3 Neutrino interaction cross-section systematic parameters that receive a central-value tune and modify the nominal event rate predictions.
-Category 2: implemented in the same way for ND and FD, but on-axis ND data alone is not sufficient to constrain these parameters. They may be constrained by additional ND samples in future analyses, such as off-axis measurements. -Category 3: implemented only in the FD. Examples are parameters which only affect ν e and ν e rates which are small and difficult to precisely isolate from background at the ND.
All GENIE uncertainties (original or modified), given in Tables 1 and 2, are all treated as Category 1. Table 4, which describes the uncertainties beyond those available within GENIE, includes a column identifying which of these categories describes the treatment of each additional uncertainty.  Table 4 List of extra interaction model uncertainties in addition to those provided by GENIE, and the category to which they belong in the analysis. Note that in this analysis, the NC norm. systematic is not applied at the FD, as described in the text.

The Near Detector Simulation and Reconstruction
The ND hall will be located at Fermi National Accelerator Laboratory (Fermilab), 574 m from where the protons hit the beam target, and 60 m underground. The baseline design for the DUNE ND system consists of a LArTPC with a downstream magnetized multi-purpose detector (MPD), and an on-axis beam monitor. Additionally, it is planned for the LArTPC and MPD to be movable perpendicular to the beam axis, to take measurements at a number of off-axis angles. The use of off-axis angles is complementary to the on-axis analysis described in this work through the DUNE-PRISM concept, originally developed in the context of the J-PARC neutrino beamline in Ref.
[63]. We note that there are many possible ND samples which are not included in the current analysis, but which may either help improve the sensitivity in future, or will help control uncertainties to the level assumed here. These include: neutrinoelectron scattering studies, which can independently constrain the flux normalization to ∼2% [64]; additional flux constraints from the low-ν method, which exploits the fact that the low energy transfer (lowν) cross section is roughly constant with neutrino energy [65][66][67][68][69][70]; and using interactions on the gaseous argon (GAr) in the MPD. There is also the potential to include events where the muon does not pass through the MPD, e.g. using multiple Coulomb scattering to estimate the muon momentum [71]. The LArTPC is a modular detector based on the ArgonCube design [72], with fully-3D pixelated readout [73] and optical segmentation [74]. These features greatly reduce reconstruction ambiguities that hamper monolithic, projective-readout time projection chambers (TPCs), and enable the ND to function in the high-intensity environment of the DUNE ND site. Each module is itself a LArTPC with two anode planes and a shared central cathode. The active dimensions of each module are 1×3×1 m (x×y×z), where the z direction is along the neutrino beam axis, and the y direction points upward. Charge drifts in the ±x direction, with a maximum drift distance of 50 cm for ionization electrons. The full liquid argon (LAr) detector consists of an array of modules in a single cryostat. The minimum active size required for full containment of hadronic showers in the LBNF beam is 3 × 4 × 5 m. High-angle muons can also be contained by extending the width to 7 m. For this analysis, 35 modules are arranged in an array 5 modules deep in the z direction and 7 modules across in x so that the total active dimensions are 7 × 3 × 5 m. The total active LAr volume is 105 m 3 , corresponding to a mass of 147 tons.
The MPD used in the analysis consists of a high-pressure gaseous argon time-projection chamber (GArTPC) in a cylindrical pressure vessel at 10 bar, surrounded by a granular, high-performance electromagnetic calorimeter, which sits immediately downstream of the LAr cryostat. The pressure vessel is 5 m in diameter and 5 m long. The TPC is divided into two drift regions by a central cathode, and filled with a 90%:10% Ar:CH 4 gas mixture, such that 97% of neutrino interactions will occur on the Ar target. The GArTPC is described in detail in Ref.
[5]. The electromagnetic calorimeter (ECAL) is composed of a series of absorber layers followed by arrays of scintillator and is described in Ref. [75]. The entire MPD sits inside a magnetic field, which allows the MPD to precisely measure the momentum and discriminate the sign of particles passing through it.
Neutrino interactions are simulated in the active volume of the LArTPC. The propagation of neutrino interaction products through the LArTPC and MPD detector volumes is simulated using a Geant4-based model [23]. Pattern recognition and reconstruction software has not yet been developed for the ND. Instead, we Bottom: Acceptance as a function of hadronic energy; the black line is for the full Fiducial Volume (FV) while the red line is for a 1 × 1 × 1 m 3 volume in the center, where the acceptance is higher due to the better hadron containment. The blue curve shows the expected distribution of true hadronic energy in the DUNE ND flux normalized to unity; 56% of events have hadronic energy below 1 GeV where the acceptance is high.
perform a parameterized reconstruction based on true energy deposits in active detector volumes as simulated by Geant4.
The analysis described here uses events originating in the LAr component, within a fiducial volume (FV) that excludes 50 cm from the sides and upstream edge, and 150 cm from the downstream edge of the active region, for a total of 6 × 2 × 3 m 2 . Muons with kinetic energy greater than ∼1 GeV typically exit the LAr. An energetic forward-going muon will pass through the ECAL and into the gaseous TPC, where its momentum and charge are reconstructed by curvature. For these events, it is possible to differentiate between µ + and µ − event by event. Muons that stop in the LAr or ECAL are reconstructed by range. Events with wideangle muons that exit the LAr and do not match to the GArTPC are rejected, as the muon momentum is not reconstructed. The asymmetric transverse dimensions of the LAr volume make it possible to reconstruct wide-angle muons with some efficiency.
The charge of muons stopping in the LAr volume cannot be determined event by event. However, the wrong-sign flux is predominantly concentrated in the high-energy tail, where leptons are more likely to be forward and energetic. In FHC beam running, the wrongsign background in the focusing peak is negligibly small, and µ − is assumed for all stopping muon tracks. In RHC beam running, the wrong-sign background is larger in the peak region. Furthermore, high-angle leptons are generally at higher inelasticity, which enhances the wrong-sign contamination in the contained muon subsample. To mitigate this, a Michel electron is required at the end of the track. The wrong-sign µ − captures on Ar with 75% probability, effectively suppressing the relative µ − component by a factor of four.
True muons and charged pions are evaluated as potential muon candidates. The track length is determined by following the true particle trajectory until it undergoes a hard scatter or ranges out. The particle is classified as a muon if its track length is at least 1 m, and the mean energy deposit per centimeter of track length is less than 3 MeV. The mean energy cut rejects tracks with detectable hadronic interactions. The minimum length requirement imposes an effective threshold on true muons of about 200 MeV kinetic energy, but greatly suppresses potential NC backgrounds with lowenergy, non-interacting charged pions. Charged-current events are required to have exactly one muon, and if the charge is reconstructed, it must be of the appropriate charge.
As in the FD reconstruction described in Section 5, hadronic energy in the ND is reconstructed by summing all charge deposits in the LAr active volume that are not associated with the muon. To reject events where the hadronic energy is poorly reconstructed due to particles exiting the detector, a veto region is defined as the outer 30 cm of the active volume on all sides. Events with more than 30 MeV total energy deposit in the veto region are excluded from the analysis. This leads to an acceptance that decreases as a function of hadronic energy, as shown in the bottom panel of Figure 5. Neutron energy is typically not observed, resulting in poor reconstruction of events with energetic neutrons, as well as in events where neutrons are produced in secondary interactions inside the detector. The reconstructed neu-trino energy is the sum of the reconstructed hadronic energy and the reconstructed muon energy.
The oscillation analysis presented here includes samples of ν µ andν µ charged-current interactions originating in the LAr portion of the ND, as shown in Figure 6. These samples are binned in two dimensions as a function of reconstructed neutrino energy and inelasticity, and E rec ν are the reconstructed muon and neutrino energies, respectively. Backgrounds to ν ( -) µ CC arise from NC π ± production where the pion leaves a long track and does not shower. Muons below about 400 MeV kinetic energy have a significant background from charged pions, so these CC events are excluded from the selected sample. Wrong-sign contamination in the beam is an additional background, particularly at low reconstructed neutrino energies in RHC.

The Far Detector Simulation and Reconstruction
The 40-kt DUNE FD consists of four separate LArTPC detector modules, each with a FV of at least 10 kt, installed ∼1.5 km underground at the Sanford Underground Research Facility (SURF) [76]. DUNE is committed to deploying both single-phase [77] and dual-phase [78] LArTPC technologies, and is investigating advanced detector designs for the fourth detector module. As such, the exact order of construction and number of modules of each design is unknown. In this work, the FD reconstruction performance is assessed assuming a single-phase design for all four modules, which does not fully exploit the benefits of different technologies with independent systematics in the sensitivity studies. A full simulation chain has been developed, from the generation of neutrino events in a Geant4 model of the FD geometry, to efficiencies and reconstructed neutrino energy estimators of all samples used in the sensitivity analysis.
The total active LAr volume of each single-phase DUNE FD detector module is 13.9 m long, 12.0 m high and 13.3 m wide, with the 13.3 m width in the drift direction subdivided into four independent drift regions, with two shared cathodes. Full details of the single-phase detector module design can be found in Ref. [79]. The total active volume of each module is ∼13 kt, the FV of at least 10 kt is defined by studies of neutrino energy resolution, using the neutrino energy estimators described below. At the anode, there are two wrapped-wire readout induction planes, which are offset by ±35.7 • to the vertical, and a vertical collection plane. Neutrino interactions of all flavors are simulated such that weights can be applied to produce samples for any set of oscillation parameters. The interaction model described in Section 3 was used to model the neutrino-argon interactions in the volume of the cryostat, and the final-state particles are propagated in the detector through Geant4. The electronics response to the ionization electrons and scintillation light is simulated to produce digitized signals in the wire planes and photon detectors (PDs) respectively.
Raw detector signals are processed using algorithms to remove the impact of the LArTPC electric field and electronics response from the measured signal, to identify hits, and to cluster hits that may be grouped together due to proximity in time and space. Clusters from different wire planes are matched to form highlevel objects such as tracks and showers. These high level objects are used as inputs to the neutrino energy reconstruction algorithm.
The energy of the incoming neutrino in CC events is estimated by adding the lepton and hadronic energies reconstructed using the Pandora toolkit [80,81]. If the event is selected as ν µ CC, the neutrino energy is estimated as the sum of the energy of the longest reconstructed track and the hadronic energy. The energy of the longest reconstructed track is estimated from its range if the track is contained in the detector. If the longest track exits the detector, its energy is estimated from multiple Coulomb scattering. The hadronic energy is estimated from the charge of reconstructed hits that are not in the longest track, and corrections are applied to each hit charge for recombination and the electron lifetime. An additional correction is made to the hadronic energy to account for missing energy due to neutral particles and final-state interactions.
If the event is selected as ν e CC, the energy of the neutrino is estimated as the sum of the energy of the reconstructed electromagnetic (EM) shower with the highest energy and the hadronic energy. The former is estimated from the charges of the reconstructed hits in the shower, and the latter from the hits not in the shower; the recombination and electron lifetime corrections are applied to the charge of each hit. The same hadronic shower energy calibration is used for both ν andν based on a sample of ν andν events.
In the energy range of 0.5-4 GeV that is relevant for oscillation measurements, the observed neutrino energy resolution is ∼15-20%, depending on lepton flavor and reconstruction method. The muon energy resolution is 4% for contained tracks and 18% for exiting tracks. The electron energy resolution is approximately 4% ⊕ 9%/ √ E, with some shower leakage that gives rise to a non-Gaussian tail that is anticorrelated with the hadronic energy measurement. The hadronic en-ergy resolution is 34%, which could be further improved by identifying individual hadrons, adding masses of charged pions, and applying particle-specific recombination corrections. It may also be possible to identify final-state neutrons by looking for neutron-nucleus scatters, and use event kinematics to further inform the energy estimate. These improvements are under investigation and are not included in this analysis.
Event classification is carried out through image recognition techniques using a convolutional neural network, named convolutional visual network (CVN). Detailed descriptions of the CVN architecture can be found in Ref. [82]. The primary goal of the CVN is to efficiently and accurately produce event selections of the following interactions: ν µ CC and ν e CC in FHC, andν µ CC andν e CC in RHC. In order to build the training input to the DUNE CVN three images of the neutrino interactions are produced, one for each of the three readout views of the LArTPC, using the reconstructed hits on individual wire planes. Each pixel contains information about the integrated charge in that region. An example of a simulated 2.2 GeV ν e CC interaction is shown in a single view in Figure 7 demonstrating the fine-grained detail available from the LArTPC technology.
The CVN is trained using approximately three million simulated neutrino interactions. A statistically independent sample is used to generate the physics measurement sensitivities. The training sample is chosen to ensure similar numbers of training examples from the different neutrino flavors. Validation is performed to ensure that similar classification performance is obtained for the training and test samples to ensure that the CVN is not overtrained. For the analysis presented here, we use the CVN score for each interaction to belong to one of the following classes: ν µ CC, ν e CC, ν τ CC and NC. The ν e CC score distribution, P (ν e CC), and the ν µ CC score distribution, P (ν µ CC), are shown in Figure 8. Excellent separation between the signal and background interactions is seen in both cases. The event selection require-ment for an interaction to be included in the ν e CC (ν µ CC) is P (ν e CC) > 0.85 (P (ν µ CC) > 0.5), optimized to produce the best sensitivity to charge parity (CP) violation. Since all of the flavor classification scores must sum to unity, the interactions selected in the two event selections are completely independent. The same selection criteria are used for both FHC and RHC beam running. Figure 9 shows the efficiency as a function of reconstructed energy (under the electron neutrino hypothesis) for the ν e event selection, and the corresponding selection efficiency for the ν µ event selection. The ν e and ν µ efficiencies in both FHC and RHC beam modes all exceed 90% in the neutrino flux peak.
The ability of the CVN to identify neutrino flavor is dependent on its ability to resolve and identify the charged lepton. Backgrounds originate from the misidentification of charged pions for ν µ disappearance, and photons for ν e appearance. The probability for these backgrounds to be introduced varies with the momentum and isolation of the energy depositions from the pions and photons. The efficiency was also observed to drop as a function of track/shower angle (with respect to the incoming neutrino beam direction) when energy depositions aligned with wire planes. The shapes of the efficiency functions in lepton momentum, lepton angle, and hadronic energy fraction (inelasticity) are all observed to be consistent with results from previous studies, including hand scans of LArTPC simulations. The CVN is susceptible to bias if there are features in the data that are not present in the simulation, so before its use on data, it will be important to comprehensively demonstrate that the selection is not sensitive to the choice of reference models. A discussion of the bias studies performed so far, and those planned in future, can be found in Ref. [82].

Expected Far Detector Event Rate and Oscillation Parameters
In this work, FD event rates are calculated assuming the following nominal deployment plan, which is based on a technically limited schedule: -Start of beam run: two FD module volumes for total fiducial mass of 20 kt, 1.2 MW beam -After one year: add one FD module volume for total fiducial mass of 30 kt -After three years: add one FD module volume for total fiducial mass of 40 kt -After six years: upgrade to 2.4 MW beam Table 5 shows the conversion between number of years under the nominal staging plan, and kt-MW-years, which are used to indicate the exposure in this analysis. For all studies shown in this work, a 50%/50% ratio of FHC to RHC data was assumed, based on studies which showed a roughly equal mix of running produced a nearly optimal δ CP and mass ordering sensitivity. The exact details of the run plan are not included in the staging plan.  Table 5 Conversion between number of years in the nominal staging plan, and kt-MW-years, the two quantities used to indicate exposure in this analysis.
Event rates are calculated with the assumption of 1.1 ×10 21 protons on target (POT) per year, which assumes a combined uptime and efficiency of the Fermilab accelerator complex and the LBNF beamline of 57% [5]. Figures 10 and 11 show the expected rate of selected events for ν e appearance and ν µ disappearance, respectively, including expected flux, cross section, and oscillation probabilities, as a function of reconstructed neutrino energy at a baseline of 1285 km. The spectra are shown for a 3.5 year (staged) exposure each for FHC and RHC beam modes, for a total run time of seven years. The rates shown are scaled to obtain different exposures. Tables 6 and 7 give the integrated rate for the ν e appearance and ν µ disappearance spectra, respectively. Note that the total rates are integrated over the range of reconstructed neutrino energies used in the analysis, 0.5-10 GeV. The nominal neutrino oscillation parameters used in Figures 10 and 11 and the uncertainty on those parameters (used later in the analysis) are taken from the NuFIT 4.0 [9,83] global fit to neutrino data, and their values are given in Table 8. See also Refs. [10] and [11] for other recent global fits.
As can be seen in Figure 10, the background to ν e appearance is composed of: (1) CC interactions of ν e andν e intrinsic to the beam; (2) misidentified NC interactions; (3) misidentified ν µ andν µ CC interactions; and (4) ν τ andν τ CC interactions in which the τ 's decay leptonically into electrons/positrons. NC and ν τ backgrounds emanate from interactions of higherenergy neutrinos that feed down to lower reconstructed neutrino energies due to missing energy in unreconstructed final-state neutrinos. The selected NC and CC ν µ generally include an asymmetric decay of a relatively high energy π 0 coupled with a prompt photon conversion. As can be seen in Figure 11, the backgrounds to the ν µ disappearance are due to wrong-sign ν µ interactions, which cannot easily be distinguished in the unmagnetized DUNE FD, and NC interactions, where a pion has been misidentified as the primary muon. As expected, the ν µ background in RHC is much larger than theν µ background in FHC.

Detector Uncertainties
Detector effects impact the event selection efficiency as well as the reconstruction of quantities used in the oscillation fit, such as neutrino energy. The main sources of detector systematic uncertainties are limitations of the expected calibration and modeling of particles in the detector.
The ND LArTPC uses similar technology to the FD, but important differences lead to uncertainties that do not fully correlate between the two detectors. First, the readout technology is different, as the ND LArTPC uses pixels as well as a different, modular photon detector. Therefore, the charge response will be different between near and far detectors due to differences in electronics readout, noise, and local effects like alignment. Second, the high-intensity environment of the ND complicates associating detached energy deposits to events, a problem which is not present in the FD. Third, the calibration strategies will be different. For example, the ND has a high-statistics calibration sam-  Fig. 10 νe andνe appearance spectra: reconstructed energy distribution of selected νe CC-like events assuming 3.5 years (staged) running in the neutrino-beam mode (top) and antineutrino-beam mode (bottom), for a total of seven years (staged) exposure. Statistical uncertainties are shown on the datapoints. The plots assume normal mass ordering and include curves for δ CP = −π/2, 0, and π/2. ple of through-going, momentum-analyzed muons from neutrino interactions in the upstream rock, which is not available with high statistics for the FD. Finally, the reconstruction efficiency will be inherently different due to the relatively small size of the ND. Containment of charged hadrons will be significantly worse at the ND, especially for events with energetic hadronic showers or with vertices near the edges of the FV. An uncertainty on the overall energy scale is included in the analysis presented here, as well as particle response uncertainties that are separate and uncorrelated between four species: muons, charged hadrons, neutrons, and electromagnetic showers. In the ND, muons reconstructed by range in LAr and by curvature in the MPD are treated separately. The energy scale and particle response uncertainties are allowed to vary with energy; each term is described by three free parameters:

Sample
Expected Events  Table 6 νe andνe appearance rates: integrated rate of selected νe CC-like events between 0.5 and 10.0 GeV assuming a 3.5-year (staged) exposure in the neutrino-beam mode and antineutrino-beam mode. The rates are shown for both NO and IO, and signal events are shown for both δ CP = 0 and δ CP = −π/2.
where E rec is the nominal reconstructed energy, E rec is the shifted energy, and p 0 , p 1 , and p 2 are free fit parameters that are allowed to vary within a priori constraints. Note that the parameters produce a shift to the kinematic variables in an event, as opposed to simply assigning a weight to each simulated event. The energy scale and resolution parameters are conservatively treated as uncorrelated between the ND and FD. With a better understanding of the relationship between ND and FD calibration and reconstruction techniques, it   Table 8 Central value and relative uncertainty of neutrino oscillation parameters from a global fit [9,83] to neutrino oscillation data. The matter density is taken from Ref. [84]. Because the probability distributions are somewhat non-Gaussian (particularly for θ 23 ), the relative uncertainty is computed using 1/6 of the 3σ allowed range from the fit, rather than 1/2 of the 1σ range. For θ 23 , θ 13 , and ∆m 2 31 , the best-fit values and uncertainties depend on whether normal mass ordering (NO) or inverted mass ordering (IO) is assumed. may be possible to correlate some portion of the energy response. The full list of assumed energy scale uncertainties is given as Table 9. In addition to the uncertainties on the energy scale, uncertainties on energy resolutions are also included. These are treated as fully uncorrelated between the near and far detectors and are taken to be 2% for muons, charged hadrons, and EM showers and 40% for neutrons.

Particle type
Allowed variation p, π ± 5% 5% 5% e, γ, π 0 2.5% 2.5% 2.5% n 20% 30% 30% Table 9 Uncertainties applied to the energy response of various particles. p 0 , p 1 , and p 2 correspond to the constant, square root, and inverse square root terms in the energy response parameterization given in Equation 2. All are treated as uncorrelated between the ND and FD.
The scale of these assumed uncertainties is motivated by what has been achieved in recent experiments, including calorimetric based approaches (NOvA, MINERvA) and LArTPCs (LArIAT, MicroBooNE, Ar-goNeuT). The DUNE performance is expected to significantly exceed the performance of these current surfacebased experiments. NOvA [44] has achieved < 1% (5%) uncertainties on the energy scale of muons (protons). Uncertainties associated to the pion and proton reinteractions in the detector medium are expected to be controlled from ProtoDUNE and LArIAT data, as well as the combined analysis of low density (gaseous) and high density (LAr) NDs. Uncertainties in the E field also contribute to the energy scale uncertainty [85], and calibration is needed (with cosmics at ND, laser system at FD) to constrain the overall energy scale. The recombination model will continue to be validated by the suite of LAr experiments and is not expected to be an issue for nominal field provided minimal E field distortions. Uncertainties in the electronics response are controlled with a dedicated charge injection system and validated with intrinsic sources, Michel electrons and 39 Ar.
The response of the detector to neutrons is a source of active study and will couple strongly to detector technology. The validation of neutron interactions in LAr will continue to be characterized by dedicated measurements (e.g., CAPTAIN [86,87]) and the LAr program (e.g., ArgoNeuT [88]). However, the association of the identification of a neutron scatter or capture to the neutron's true energy has not been demonstrated, and significant reconstruction issues exist, so a large uncertainty (20%) is assigned comparable to the observations made by MINERvA [89] assuming they are attributed entirely to the detector model. Selection of photon candidates from π 0 is also a significant reconstruction challenge, but a recent measurement from MicroBooNE indicates this is possible and the reconstructed π 0 invariant mass has an uncertainty of 5%, although with some bias [90].
The p 1 and p 2 terms in Equation 2 allow the energy response to vary as a function of energy. The energy dependence is conservatively assumed to be of the same order as the absolute scale uncertainties given by the p 0 terms.
In addition to impacting energy reconstruction, the E field model also affects the definition of the FD fiducial volume, which is sensitive to electron drift. An additional 1% uncertainty is assumed on the total fiducial mass, which is conservatively treated as uncorrelated between the ν µ and ν e samples due to the potential distortion caused by large electromagnetic showers in the electron sample. These uncertainties affect only the overall normalization, and are called FV numu FD and FV nue FD in Figure 12.
The ND and FD have different acceptance to CC events due to the very different detector sizes. The FD is sufficiently large that acceptance is not expected to vary significantly as a function of event kinematics. However, the ND selection requires that hadronic showers be well contained in LAr to ensure a good energy resolution, resulting in a loss of acceptance for events with energetic hadronic showers. The ND also has regions of muon phase space with lower acceptance due to tracks exiting the side of the TPC but failing to match to the MPD, which are currently not used in the sensitivity analysis.
Uncertainties are evaluated on the muon and hadron acceptance of the ND. The detector acceptance for muons and hadrons is shown in Figure 5. Inefficiency at very low lepton energy is due to events being misreconstructed as neutral current. For high energy, forward muons, the inefficiency is only due to events near the edge of the FV where the muon happens to miss the MPD. At high transverse momentum, muons begin to exit the side of the LAr active volume, except when they happen to go along the 7 m axis. The acceptance is sensitive to the modeling of muons in the detector.
An uncertainty is estimated based on the change in the acceptance as a function of muon kinematics.
Inefficiency at high hadronic energy is due to the veto on more than 30 MeV deposited in the outer 30 cm of the LAr active volume. Rejected events are typically poorly reconstructed due to low containment, and the acceptance is expected to decrease at high hadronic energy. Similar to the muon reconstruction, this acceptance is sensitive to detector modeling, and an uncertainty is evaluated based on the change in the acceptance as a function of true hadronic energy.

Sensitivity Methods
Previous DUNE sensitivity predictions have used the GLoBES framework [7, 91,92]. In this work, fits are performed using the CAFAna [93] analysis framework, developed originally for the NOvA experiment. Systematics are implemented using one-dimensional response functions for each analysis bin, and oscillation weights are calculated exactly, in fine (50 MeV) bins of true neutrino energy. For a given set of inputs, flux, oscillation parameters, cross sections, detector energy response matrices, and detector efficiency, an expected event rate can be produced. Minimization is performed using the minuit [94] package.
Oscillation sensitivities are obtained by simultaneously fitting the ν µ → ν µ ,ν µ →ν µ (Figure 11), ν µ → ν e , andν µ →ν e (Figure 10) FD spectra along with the ν µ FHC andν µ RHC samples from the ND (Figure 6). In the studies, all oscillation parameters shown in Table 8 are allowed to vary. Gaussian penalty terms (taken from Table 8) are applied to θ 12 and ∆m 2 12 and the matter density, ρ, of the Earth along the DUNE baseline [84]. Unless otherwise stated, studies presented here include a Gaussian penalty term on θ 13 (also taken from Table 8), which is precisely measured by experiments sensitive to reactor antineutrino disappearance [95][96][97]. The remaining parameters, sin 2 θ 23 , ∆m 2 32 , and δ CP are allowed to vary freely, with no penalty term. Note that the penalty terms are treated as uncorrelated with each other, or other parameters, which is a simplification. In particular, the reactor experiments that drive the constraint on θ 13 in the NuFIT 4.0 analysis are also sensitive to ∆m 2 32 , so the constraint on θ 13 should be correlated with ∆m 2 32 . We do not expect this to have a significant impact on the fits, and this effect only matters for those results with the θ 13 Gaussian penalty term included.
Flux, cross section, and FD detector parameters are allowed to vary in the fit, but constrained by a penalty term proportional to the pre-fit uncertainty. ND detector parameters are not allowed to vary in the fit, but their effect is included via a covariance matrix based on the shape difference between ND prediction and the "data" (which comes from the simulation in this sensitivity study). The covariance matrix is constructed with a throwing technique. For each "throw", all ND energy scale, resolution, and acceptance parameters are simultaneously thrown according to their respective uncertainties, and the modified prediction is produced by varying the relevant quantities away from the nominal prediction according to the thrown parameter values. The bin-to-bin covariance is determined by comparing the resulting spectra with the nominal prediction, in the same binning as is used in the oscillation sensitivity analysis. This choice protects against overconstraining that could occur given the limitations of the parameterized ND reconstruction described in Section 4 taken together with the high statistical power at the ND, but is also a simplification.
The compatibility of a particular oscillation hypothesis with both ND and FD data is evaluated using a negative log-likelihood ratio, which converges to a χ 2 at high-statistics [48]: where ϑ and x are the vector of oscillation parameter and nuisance parameter values respectively; M i ( ϑ, x) and D i are the Monte Carlo (MC) expectation and fake data in the ith reconstructed bin (summed over all selected samples), with the oscillation parameters neglected for the ND; ∆x j and σ j are the difference between the nominal and current value, and the prior uncertainty on the jth nuisance parameter with uncertainties evaluated and described in Sections 2, 3 and 7; and V kl is the covariance matrix between ND bins described previously. In order to avoid falling into a false minimum, all fits are repeated for four different δ CP values (-π, -π/2, 0, π/2), both mass orderings, and in both octants, and the lowest χ 2 value is taken as the minimum. Two approaches are used for the sensitivity studies presented in this work. First, Asimov studies [98] are carried out in which the fake (Asimov) dataset is the same as the nominal MC. In these, the true value of Treatment of the oscillation parameters for the simulated data set studies. Note that for some studies θ 13 has a Gaussian penalty term applied based on the NuFIT 4.0 value, and for others it is thrown uniformly within a range determined from the NuFIT 4.0 3σ allowed range.
all systematic uncertainties and oscillation parameters except those of interest (which are fixed at a test point) remain unchanged, and can vary in the fit, but are constrained by their pre-fit uncertainty. Second, studies are performed where many statistical and systematic throws are made according to their pre-fit Gaussian uncertainties, and fits of all parameters are carried out for each throw. A distribution of post-fit values is built up for the parameter of interest. In these, the expected resolution for oscillation parameters is determined from the spread in best-fit values obtained from an ensemble of throws that vary according to both the statistical and systematic uncertainties. For each throw, the true value of each nuisance parameter is chosen randomly from a distribution determined by the a priori uncertainty on the parameter. For some studies, oscillation parameters are also randomly chosen as described in Table 10. Poisson fluctuations are then applied to all analysis bins, based on the mean event count for each bin after the systematic adjustments have been applied. For each throw in the ensemble, the test statistic is minimized, and the best-fit value of all parameters is determined. The median throw and central 68% of throws derived from these ensembles are shown. Sensitivity calculations for CPV, neutrino mass ordering, and octant are performed, in addition to studies of oscillation parameter resolution in one and two dimensions. In these cases, the experimental sensitivity is quantified using a likelihood ratio as the test statistic: where χ 2 B and χ 2 A are both obtained from Equation 3, using a coherent systematic and statistical throw. The size of ∆χ 2 is a measure of how well the data can exclude model B in favor of model A, given the uncertainty in the model. For example, the sensitivity for excluding the IO in favor of the NO would be given as  Fig. 12 The ratio of post-fit to pre-fit uncertainties for various systematic parameters for a 15-year staged exposure. The red band shows the constraint from the FD only in 15 years, while the green shows the ND+FD constraints. Flux parameters are named "Flux #i" representing the ith principal flux component, cross-section parameter names are given in Section 3, and detector systematics are described in Section 7, where the p 0 , p 1 and p 2 parameters are described in Table 9.
Note that the ∆χ 2 for the mass ordering may be negative, depending on how the test is set up. The sensitivity for discovering CPV is the preference for the CP violating hypothesis over the CP conserving hypothesis, χ 2 0,π − χ 2 CPV . Post-Fit uncertainties on systematic parameters are shown for Asimov fits at the NuFIT 4.0 best-fit point to both the ND+FD samples, and the FD-only samples in Figure 12, as a fraction of the pre-fit systematic uncertainties described in Sections 2, 3, and 7. The FD alone can only weakly constrain the flux and cross-section parameters, which are much more strongly constrained when the ND is included. The ND is, however, unable to strongly constrain the FD detector systematics as they are treated as uncorrelated, and due to the treatment of ND detector systematics in a covariance matrix in Equation 3. Adding the ND does slightly increase the constraint on detector parameters as it breaks degeneracies with other parameters. Several important crosssection uncertainties are also not constrained by the ND. In particular, an uncertainty on the ratio of ν µ to ν e cross sections is totally unconstrained, which is not surprising given the lack of ND ν e samples in the current analysis. The most significant flux terms are constrained at the level of 20% of their a priori values. Less significant principal components have little impact on the observed distributions at either detector, and receive weaker constraints.  Figure 13 shows the pre-and post-fit systematic uncertainties on the FD FHC samples for Asimov fits at the NuFIT 4.0 best-fit point including both ND and FD samples with a 15 year exposure. It shows how the parameter constraints seen in Figure 12 translate to a constraint on the event rate. Similar results are seen for the RHC samples. The large reduction in the systematic uncertainties is largely due to the ND constraint on the systematic uncertainties apparent from Figure 12.

Sensitivities
In this section, various sensitivity results are presented. For the sake of simplicity, unless otherwise stated, only true normal ordering is shown. Possible variations of sensitivity are presented in two ways. Results produced using Asimovs are shown as lines, and differences between two Asimov scenarios are shown with a colored band. Note that the band in the Asimov case is purely to guide the eye, and does not denote a confidence interval. For results produced using many throws of oscillation parameters, systematic and statistical uncertainties, ∼300,000 throws were used to calculate the sensitivity for each scenario. The median sensitivity is shown with a solid line, and a transparent filled area indicates the region containing the central 68% of throws, which can be interpreted as the 1σ uncertainty on the sensitivity. Figure 14 shows the significance with which CPV (δ CP = [0, ±π]) can be observed in both NO and IO as a function of the true value of δ CP for exposures corresponding to seven and ten years of data, using the staging scenario described in Section 6, and using the toy throwing method described in Section 8 to investigate their effect on the sensitivity. This sensitivity has a characteristic double peak structure because the significance of a CPV measurement necessarily decreases around CP-conserving values of δ CP . The median CPV sensitivity reaches 5σ for a small range of values after an exposure of seven years in NO, and a broad range of values after a ten year exposure. In IO, DUNE has slightly stronger sensitivity to CPV, and reaches 5σ for a broad range of values after a seven year exposure. Note that with statistical and systematic throws, the median sensitivity never reaches exactly zero. Figure 15 shows the DUNE Asimov sensitivity to CPV in NO when the true values of θ 23 , θ 13 , and ∆m 2 32 vary within the 3σ range allowed by NuFIT 4.0. The largest effect is the variation in sensitivity with the true value of θ 23 , where degeneracy with δ CP and matter effects are significant. Values of θ 23 in the lower octant lead to the best sensitivity to CPV. The true values of θ 13 and ∆m 2 32 are highly constrained by global data and, within these constraints, do not have a dramatic impact on the sensitivity. Note that in the Asimov cases shown in Figure 15, the median sensitivity reaches 0 at CP-conserving values of δ CP (unlike the case with the throws as in Figure 14), but in regions far from CP-conserving values, the Asimov sensitivity and the median sensitivity from the throws agree well. Figure 16 shows the result of Asimov studies investigating the significance with which CPV can be determined in NO for 75% and 50% of δ CP values, and when  δ CP = −π/2, as a function of exposure in kt-MW-years, which can be converted to years using the staging scenario described in Section 6. The width of the bands show the impact of applying an external constraint on θ 13 . CP violation can be observed with 5σ significance after about seven years (336 kt-MW-years) if δ CP = −π/2 and after about ten years (624 kt-MW-years) for 50% of δ CP values. CP violation can be observed with 3σ significance for 75% of δ CP values after about 13 years of running. In the bottom plot of Figure 16, the width of the bands shows the impact of applying an external constraint on θ 13 , while in the bottom plot, the width of the bands is the result of varying the true value of sin 2 θ 23 within the NuFIT 4.0 90% C.L. allowed region. Fig. 17 Significance of the DUNE determination of the neutrino mass ordering, as a function of the true value of δ CP , for seven (blue) and ten (orange) years of exposure. The width of the transparent bands cover 68% of fits in which random throws are used to simulate statistical variations and select true values of the oscillation and systematic uncertainty parameters, constrained by pre-fit uncertainties. The solid lines show the median sensitivity. Figure 17 shows the significance with which the neutrino mass ordering can be determined in both NO and IO as a function of the true value of δ CP , for both seven and ten year exposures, including the effect of all other oscillation and systematic parameters using the toy throwing method described in Section 8. The characteristic shape results from near degeneracy between matter and CPV effects that occurs near δ CP = π/2 (−δ CP = π/2) for true normal (inverted) ordering. Studies have indicated that special attention must be paid to the statistical interpretation of neutrino mass ordering sensitivities [99][100][101] because the ∆χ 2 metric does not follow the expected chi-square function for one degree of freedom, so the interpretation of the ∆χ 2 as the sensitivity is complicated. However, it is clear from Figure 17 that DUNE is able to distinguish the mass ordering for both true NO and IO, and using the corrections from, for example, Ref. [99], DUNE would still achieve 5σ significance for the central 68% of all throws shown in Figure 17. We note that for both seven and ten years (it was not checked for lower exposures), there were no parameter throws used in generating the plots (∼300,000 each) for which the incorrect mass ordering was preferred. Figure 18 shows the DUNE Asimov sensitivity to the neutrino mass ordering when the true values of θ 23 , θ 13 , and ∆m 2 32 vary within the 3σ range allowed by NuFIT 4.0. As for CPV (in Figure 15), the largest variation in sensitivity is with the true value of θ 23 , but in this case, the upper octant leads to the best sensitivity. Again, the true values of θ 13 and ∆m 2 32 do not have a dramatic impact on the sensitivity. The median Asimov sensitivity tracks the median throws shown in Figure 17 well for the reasonably high exposures tested -this was not checked for exposures below seven years (336 kt-MW-years). Figure 19 shows the result of Asimov studies assessing the significance with which the neutrino mass ordering can be determined for 100% of δ CP values, and when δ CP = −π/2, as a function of exposure in kt-MWyears, for true NO. The width of the bands show the impact of applying an external constraint on θ 13 . The bottom plot shows the impact of varying the true value of sin 2 θ 23 within the NuFIT 4.0 90% C.L. region. As DUNE will be able to establish the neutrino mass ordering at the 5σ level for 100% of δ CP values after a relatively short period, these plots only extend to 500 kt-MW-years.
The measurement of ν µ → ν µ oscillations depends on sin 2 2θ 23 , whereas the measurement of ν µ → ν e oscillations depends on sin 2 θ 23 . A combination of both ν e appearance and ν µ disappearance measurements can probe both maximal mixing and the θ 23 octant. Figure 20 shows the sensitivity to determining the octant as a function of the true value of sin 2 θ 23 , in both NO and IO. We note that the octant sensitivity strongly depends on the use of the external θ 13 constraint.
In addition to the discovery potential for neutrino the mass ordering and CPV, and sensitivity to the θ 23 octant, DUNE will improve the precision on key param-  eters that govern neutrino oscillations, including δ CP , sin 2 2θ 13 , ∆m 2 31 , and sin 2 θ 23 . Figure 21 shows the resolution, in degrees, of DUNE's measurement of δ CP , as a function of the true value of δ CP , for true NO. The resolution on a parameter is produced from the central 68% of post-fit parameter values using many throws of the systematic and remaining oscillation parameters, and statistical throws. The resolution of this measurement is significantly better near CP-conserving values of δ CP , compared to maximally CP-violating values. For fifteen years of exposure, resolutions between 5 • -15 • are possible, depend- ing on the true value of δ CP . A smoothing algorithm has been applied to interpolate between values of δ CP at which the full analysis has been performed. Figures 22 and 23 show the resolution of DUNE's measurements of δ CP and sin 2 2θ 13 and of sin 2 2θ 23 and ∆m 2 32 , respectively, as a function of exposure in kt-MW-years. The resolution on a parameter is produced from the central 68% of post-fit parameter values using Fig. 21 Resolution in degrees for the DUNE measurement of δ CP , as a function of the true value of δ CP , for seven (blue), ten (orange), and fifteen (green) years of exposure. The width of the band shows the impact of applying an external constraint on θ 13 . many throws of the systematic other oscillation parameters, and statistical throws. As seen in Figure 21, the δ CP resolution varies significantly with the true value of δ CP , but for favorable values, resolutions near five degrees are possible for large exposure. The DUNE measurement of sin 2 2θ 13 approaches the precision of reactor experiments for high exposure, allowing a comparison between the two results, which is of interest as a test of the unitarity of the PMNS matrix.
One of the primary physics goals for DUNE is the simultaneous measurement of all oscillation parameters governing long-baseline neutrino oscillation, without a need for external constraints. Figure 24 shows the 90% constant ∆χ 2 allowed regions in the sin 2 2θ 13 -δ CP and sin 2 θ 23 -∆m 2 32 planes for seven, ten, and fifteen years of running, when no external constraints are applied, compared to the current measurements from world data. An additional degenerate lobe visible at higher values of sin 2 2θ 13 and in the wrong sin 2 θ 23 octant is present in the seven and ten year exposures, but is resolved after long exposures. The time to resolve the degeneracy with DUNE data alone depends on the true oscillation parameter values. For shorter exposures, the degeneracy observed in Figure 24 can be resolved by introducing an external constraint on the value of θ 13 . Figure 25 shows two-dimensional 90% constant ∆χ 2 allowed regions in the sin 2 θ 23 -δ CP plane with an external constraint on θ 13 applied. In this case, the degenerate octant solution has disappeared for all exposures shown. Resolution of DUNE measurements of δ CP (top) and sin 2 2θ 13 (bottom), as a function of exposure in kt-MW-years. As seen in Figure 21, the δ CP resolution has a significant dependence on the true value of δ CP , so curves for δ CP = −π/2 (red) and δ CP = 0 (green) are shown. For δ CP , the width of the band shows the impact of applying an external constraint on θ 13 . No constraint is applied when calculating the sin 2 2θ 13 resolution.  δ CP = -π/2, 0, π/2 were used in the sin 2 θ 23 -δ CP plane. It can be observed that the resolution in the value of sin 2 θ 23 is worse at sin 2 θ 23 = 0.5, and improves for values away from maximal in either octant. As was seen in Figure 21, the resolution of δ CP is smaller near the CP-conserving value of δ CP = 0, and increases towards the maximally CP-violating values δ CP = ±π/2.
The exposures required to reach selected sensitivity milestones for the nominal analysis are summarized in Table 11. Note that the sensitivity to CPV and for determining the neutrino mass ordering was shown to be dependent on the value of θ 23 in Figures 15 and 18, so these milestones should be treated as approximate. δ CP = −π/2 is taken as a reference value of maximal CPV close to the current global best fit. Similarly, a resolution of 0.004 on sin 2 2θ 13 is used as a reference as the current resolution obtained by reactor experiments.

Conclusion
The analyses presented here are based on full, endto-end simulation, reconstruction, and event selection of FD Monte Carlo and parameterized analysis of ND Monte Carlo of the DUNE experiment. Detailed uncertainties from flux, the neutrino interaction model, and detector effects have been included in the analysis. Sensitivity results are obtained using a sophisticated, custom fitting framework. These studies demonstrate that DUNE will be able to measure δ CP to high precision, unequivocally determine the neutrino mass ordering, and make precise measurements of the parameters governing long-baseline neutrino oscillation.
We note that further improvements are expected once the full potential of the DUNE ND is included in the analysis. In addition to the samples used explicitly in this analysis, the LArTPC is expected to measure numerous exclusive final-state CC channels, as well as ν e and NC events. Additionally, neutrino-electron elastic scattering [64] and the low-ν technique [65][66][67][68][69][70] may be used to constrain the flux. Additional samples of events from other detectors in the DUNE ND complex are not explicitly included here, but there is an assumption that we will be able to control the uncertainties to the level used in the analysis, and it should be understood that that implicitly relies on having a highly capable ND.   11 Exposure in years, assuming true normal ordering and equal running in neutrino and antineutrino mode, required to reach selected physics milestones in the nominal analysis, using the NuFIT 4.0 best-fit values for the oscillation parameters. The staging scenario described in Section 6 is assumed. Exposures are rounded to the nearest year.
DUNE will be able to establish the neutrino mass ordering at the 5σ level for 100% of δ CP values between two and three years. CP violation can be observed with 5σ significance after ∼7 years if δ CP = −π/2 and after ∼10 years for 50% of δ CP values. CP violation can be observed with 3σ significance for 75% of δ CP values after ∼13 years of running. For 15 years of exposure, δ CP resolution between five and fifteen degrees are possible, depending on the true value of δ CP . The DUNE measurement of sin 2 2θ 13 approaches the precision of reactor experiments for high exposure, allowing measurements that do not rely on an external sin 2 2θ 13 constraint and facilitating a comparison between the Fig. 26 Two-dimensional 90% constant ∆χ 2 confidence regions in the sin 2 θ 23 -δ CP (left) and sin 2 θ 23 -∆m 2 32 (right) planes for different oscillation parameter values and seven, ten, and fifteen years of exposure, with equal running in neutrino and antineutrino mode. The 90% C.L. region for the NuFIT 4.0 global fit is included in yellow for comparison. In all cases, an external constraint on the value of θ 13 is applied. The true oscillation parameter values used are denoted by stars, and the NuFIT 4.0 best fit values are used as the true value of all those not explicitly shown. Test values of sin 2 θ 23 = 0.42, 0.5, 0.58 were used for both top and bottom plots. In the top plot, test values of δ CP = -π/2, 0, π/2 were used. DUNE and reactor sin 2 2θ 13 results, which is of interest as a potential signature for physics beyond the standard model. DUNE will have significant sensitivity to the θ 23 octant for values of sin 2 θ 23 less than about 0.47 and greater than about 0.55. We note that the results found are broadly consistent with those found in Ref.
The measurements made by DUNE will make significant contributions to completion of the standard threeflavor mixing picture, and provide invaluable inputs to theory work understanding whether there are new symmetries in the neutrino sector and the relationship between the generational structure of quarks and leptons. The observation of CPV in neutrinos would be an important step in understanding the origin of the baryon asymmetry of the universe. The precise measurements of the three-flavor mixing parameters that DUNE will provide may also yield inconsistencies that point us to physics beyond the standard three-flavor model.