Measurement of atmospheric neutrino oscillations with very large volume neutrino telescopes

Neutrino oscillations have been probed during the last few decades using multiple neutrino sources and experimental set-ups. In the recent years, very large volume neutrino telescopes have started contributing to the field. First ANTARES and then IceCube have relied on large and sparsely instrumented volumes to observe atmospheric neutrinos for combinations of baselines and energies inaccessible to other experiments. Using this advantage, the latest result from IceCube starts approaching the precision of other established technologies, and is paving the way for future detectors, such as ORCA and PINGU. These new projects seek to provide better measurements of neutrino oscillation parameters, and eventually determine the neutrino mass ordering. The results from running experiments and the potential from proposed projects are discussed in this review, emphasizing the experimental challenges involved in the measurements.


I. INTRODUCTION
Massive, mixed neutrinos, inferred from the phenomenon of oscillations, remain until this day the only physics found beyond the original formulation of the Standard Model. While the Standard Model can be extended to account for these experimental facts, precise measurements of the parameters involved in the phenomenon are necessary to constrain the different theories that attempt to explain it.
The current knowledge favors the existence of three active neutrinos (flavor eigenstates ν e , ν µ , ν τ ) whose mixing can be fully determined by the PMNS 1 matrix U . The matrix is often parameterized as the product of three rotation matrices, related to the mixing angles θ 12 , θ 13 and θ 23 , and a complex CP phase δ: where c ij ≡ cos θ ij and s ij ≡ sin θ ij . The last matrix in the multiplication does not affect neutrino oscillations, and only exists if neutrinos are Majorana particles [1]. The three angles that determine the mixing matrix are known to a precision of 10% or better [2][3][4]. The absolute differences of the square of the masses (mass splittings ∆m 2 ij = m 2 i − m 2 j with i, j = 1, 2, 3), which play a role * juan.pablo.yanez@desy.de † kouchner@apc.univ-paris7.fr 1 Pontecorvo-Maki-Nakagawa-Sakata in oscillations, are known to a precision better than 5%. While the sign of the mass splitting between the states 1-2 (m 2 > m 1 ) is known from matter effects in solar neutrino oscillations, the relative difference between m 3 and m 1 remains unknown. The relative values of the neutrino masses are commonly referred to as the neutrino mass ordering (NMO) 2 , which has two possible options: the normal ordering (NO), with m 1 < m 2 < m 3 , and the inverted ordering (IO), with m 3 < m 1 < m 2 . The determination of the NMO is important as the parameter can discriminate between flavor symmetry models [6]. Also, the sensitivity of experiments attempting to determine the neutrino nature depend on the NMO [7,8]. Finally, knowing the NMO would help to measure the value of the δ phase, which in turn would be an important step forward towards solving the fundamental question of the prevalence of matter over anti-matter in the Universe. Better measurements of all the parameters involved in neutrino oscillations are therefore necessary to understand if the current model is correct and how to incorporate it to the Standard Model.
In this view, atmospheric neutrinos remain a promising tool for studying oscillations: they cover a wide energy range, from MeV to TeV, and can reach a detector after traveling distances from a few to about 12700 km when they cross the Earth. No man-made beam covers a similar parameter space. However, the flux strongly decreases with energy. Detecting it implies building large detectors, such as very large volume neutrino telescopes (VLVNTs) to study atmospheric neutrino oscillations.
In the recent years first ANTARES and then Ice- 2 The determination of the mass ordering is sometimes confused with the determination of the neutrino "mass hierarchy", which requires additional information on the absolute scale of the neutrino masses [5].
Typeset by REVT E X arXiv:1509.08404v2 [hep-ex] 3 Nov 2015 Cube/DeepCore have proven that these studies are feasible by analyzing interactions of neutrinos with energy as low as 15 GeV. The current result from IceCube on sin 2 θ 23 and ∆m 2 32 reaches a precision which is only a factor of three to four times less stringent than global fits which combine all available data [2,3,9]. Building upon the success of these studies, proposals of new, more densely instrumented telescopes, as extensions or part of new projects, have appeared. PINGU and ORCA aim to improve the precision of these measurements and reduce the energy threshold to a few GeV, where matter effects are strong, and use these effects to measure the NMO.
This review begins by covering, in Section II, the current knowledge on atmospheric neutrinos, and how oscillations affect their flux. Section III describes the design and operation of VLVNTs. Special attention is paid to relevant sources of uncertainty. The neutrino oscillation results produced by these experiments until this date are covered in Section IV. Section V discusses the possible studies with future detectors, and a short summary is given in Section VI.

II. ATMOSPHERIC NEUTRINO OSCILLATIONS
The flux of atmospheric neutrinos in the energies relevant for a VLVNT, together with how the flux is modified by neutrino oscillations for those neutrinos that cross the Earth, are the topics covered in this Section.

A. A neutrino beam from cosmic rays
Cosmic rays (CR) continuously arrive at the Earth from all directions and interact with nuclei in the atmosphere at altitudes of about 25 km above sea level and initiate showers of particles. During the shower development, charged mesons are produced that eventually decay in comparable numbers of muons and neutrinos: µ − → e − +ν e + ν µ , The atmospheric muons produced in air showers can travel long distances before they decay. They are able to penetrate deep into the Earth, depending on their energy and the material that they are crossing [1], and constitute the dominant background in the measurement of atmospheric neutrinos.
which changes by up to 20% both for electron and muon (anti)neutrinos. Apart from that, the models agree that atmospheric neutrinos follow a power-law energy spectrum with a spectral index close to 3 in the energy range E ν = [3 − 100] GeV. Measurements from [13] estimate an uncertainty of ±0.04 on the spectral index of atmospheric neutrinos. A similar uncertainty has been also derived from varying the underlying cosmic ray model in neutrino flux calculations [23]. Muon neutrinos dominate the flux, and also have the hardest spectral index (see Fig. 1). Since electron (anti)neutrinos mainly come from muons that lose energy before decaying (see Eq. 2) their spectral index is softer, and their relative contribution to the total neutrino flux depends on energy and direction. The direction-averaged flux ofν µ is between 1.1 and 1.3 times smaller than that of ν µ , depending on the energy. The electron (anti)neutrino direction-averaged flux at a few GeV is about 2.5 times smaller than its muon (anti)neutrino counterpart. Figure 2 shows iso-contours of the neutrino flux flavor ratio as a function of energy and zenith angle for neutrinos that cross the Earth. The flux difference between ν µ and ν e grows both with energy and | cos θ z |. Already at 40 GeV the direction-averaged ratio is close to four.
Above 10 GeV the neutrino flux as a function of zenith angle is almost symmetric around cos θ z = 0. The angular dependence of the flux at the detection site is influenced by hadronization processes, the local atmospheric density and geomagnetic effects at the interaction point. The uncertainties associated to these processes result in energy-dependent modeling errors of the arrival zenith angle by up to 20 % on the ratio ν/ν for muon neutrinos and 8 % for electron neutrinos [23].

B. Neutrino oscillations at Eν ≥ 5 GeV
The flux of atmospheric neutrinos at a detection site is modified by oscillations. The oscillations that act over the L/E parameter space accessible with atmospheric neutrinos that cross the Earth (L/E ∼ 10 1 − 10 3 km/GeV) are mainly driven by the large mass splitting, ∆m 2 32 ∆m 2 31 and the mixing angles θ 13 , θ 23 . These are therefore parameters that VLVNTs are sensitive to.
Neutrinos propagating in matter are subject to a potential due to coherent forward scattering with the particles in the medium [24]. For explanatory purposes, we consider the case of neutrinos traveling through matter with constant electron density that results in a potential A = ±2 √ 2 G F n e (x) E ν , where G F is the Fermi constant, and the plus (minus) sign corresponds to neutrinos (antineutrinos). Computation of neutrino oscillation probabilities for the relevant energies has been done in [25], from where we take the approximations for the ν µ to ν e transition, given by while the survival probability of ν µ is a somewhat more complicated expression, and the transitions to ν τ are simply In these expressions ∆ ≡ ∆m 2 31 and ∆ M is the effective mass splitting in matter, given by The superscript M also accompanies θ 13 , whose effective value in matter is The mixing angle θ 23 is known to be close to maximal (∼ π/4), and |∆m 2 31 | is of the order of 10 −3 eV 2 [1]. The angle θ 13 has been recently measured and found to be small but non-zero [26][27][28]. It is then the case that θ M 13 can acquire any value, depending on the neutrino energy and the electron density of the material being crossed, as shown in Eq. 7. For a low electron density or neutrino energy, the parameters (and equations) in vacuum are recovered. A particularly interesting case appears when A = ∆m 2 31 cos 2θ 13 , which gives θ M 13 = π/4, maximizing the mixing between states 1-3, i.e. a resonance appears [29]. The effective mass splitting acquires its minimum value under this condition, and is reduced by a factor sin 2θ 13 .
The resonance that leads to maximal 1-3 mixing can only happen if the potential A and the mass difference ∆m 2 31 have the same sign, so for neutrinos in the case of NO and anti-neutrinos in the case of IO. Identifying whether the resonance takes place in neutrinos or antineutrinos is a way to identify the NMO.
For A ∆m 2 31 cos 2θ 13 a saturation effect occurs, where the effective angle in matter goes to π/2 and the effective mass splitting is then well approximated by A. In the saturated regime transitions of the type ν e → ν µ , given in Eq. 3, are suppressed by the factor sin 2 2θ M 13 . The effective matter parameters also modify P µµ by making the last two terms in Eq. 4 go to zero, resulting in the simpler expression with all the oscillated ν µ turning into ν τ . The CP-violating phase δ is not present in the approximate formulas shown. The reason is that the parameter δ always appears in oscillation probabilities accompanied by a factor ∆m 2 21 /∆m 2 31 , which suppresses its contribution [30]. Note, however, that the approximations presented here serve the purpose of explaining the main features of neutrino oscillations in matter. Figures contained in this review, as well as the latest data analyses discussed, use numerical calculations of oscillation probabilities that do not rely on simplified analytical expressions.
C. An oscillating atmospheric neutrino flux The atmospheric neutrinos under consideration, of a few GeV, are mostly ν µ +ν µ produced around a height of 25 km in the atmosphere, where the matter density is low enough to be approximated as vacuum. For most production angles the neutrinos proceed to cross the Earth, which has a non-negligible matter density. Earth's matter profile can be well explained as concentric shells, each one with a constant density [31]. To study the transitions that take place, consider the oscillation parameters from [3] and the electron number density of the mantle, n e = 2.5 cm −3 N A , where N A is Avogadro's number. Neutrinos crossing the mantle experience the resonance around E ν 6 GeV (see Eq. 7), while the saturation condition A ∆m 2 31 cos 2θ 13 is fulfilled already at E ν ∼ 12 GeV. Neutrinos measured by VLVNTs then experience oscillations either in the resonant or saturated regime, depending on the energy threshold of the detector.
Another interesting effect takes place on neutrinos that cross the Earth's core. These neutrinos experience a symmetric electron density profile that changes abruptly. For the right combination of neutrino energy and electron densities, a so-called parametric resonance can appear [32][33][34][35]. The effect, however, is not the dominant one at the energies at which future projects (Section V) will be sensitive to.
In the saturated regime atmospheric neutrino oscillations are independent of the mass ordering, dominated by ν µ → ν τ transitions and well described by Eq. 8. Near the resonance condition transitions involving electron (anti)neutrinos also play a role, and patterns become complex. Figure 3 shows the survival probabilities of ν e and ν µ for neutrinos and normal mass ordering. The original electron neutrino flux is expected to fully disappear due to matter effects over E ν = [5,8] GeV and cos θ z = [−0.9, −0.5]. The suppression of these oscillations due to saturation can be observed at about 10 GeV. The survival probability of ν µ shows abrupt changes that are due to the effects of matter. Muon neutrinos oscillate even if the resonance conditions are not fulfilled, which makes the effects of the resonance less obvious than for electron neutrinos. Resonant matter effects appear in the ν µ survival probability as modifications on the otherwise smooth and periodic disappearance pattern, as shown in Fig. 3. Saturation is reached above 15 GeV and the survival probability becomes smooth. Figure 4 shows the transition probabilities of ν e and ν µ into different flavors for the arrival direction cos θ z = −0.7 assuming a normal mass ordering. They correspond to a one-dimensional projection of Fig. 3 along the dashed line. The bands demonstrate how the uncertainties on the oscillation parameters impact the expected probabilities. For ν e it is easy to observe the same disappearance as in Fig. 3, with neutrinos oscillating equally into ν µ and ν τ . Transitions of ν µ to other flavors are complicated by matter effects, which open the ν µ ↔ ν e channel and thus modify the survival probability of ν µ .
Measurements of neutrino fluxes above the saturation energy of about 15 GeV are largely independent of θ 13 , the neutrino/antineutrino admixture of the sample and the ordering of neutrino masses. They provide excellent data for determining sin 2 θ 23 as well as |∆m 2 31 |. The NMO can only be accessed with neutrinos below 15 GeV, where matter induced resonances occur either for neutrinos or antineutrinos. The survival probability of muon (anti)neutrinos, the main component of atmospheric neutrinos, is modified by matter effects by about 20%. As will be discussed in Section III, VLVNT cannot separate neutrinos from antineutrinos event-wise, and instead rely on the ν/ν flux ratio and the difference in cross sections to identify whether oscillation probabilities of neutrinos or antineutrinos are modified by matter effects.
An interesting feature introduced by matter effects is that, instead of oscillating fully into ν τ , muon neutrinos also change into ν e . Transitions of these type are almost symmetric between the two flavors (see Fig. 4), but since the flux of ν µ is several times that of ν e at the energy and zenith angle of interest (see Fig. 2), the net effect is a significant excess of electron neutrinos with respect to the original ν e flux. In the NO, the ν e flux is enhanced, while for an IO the enhancement is realized forν e . Because of the initial ν µ /ν µ flux ratio and the differences in the ν e /ν e cross sections, different orderings result in a different number of detected events. Figure 5 shows the ratio between expected interaction rates of ν e +ν e for normal and inverted orderings, including all of the oscillation channels. A factor of 2.1 is applied to neutrinos to account for the difference in cross sections. The normal mass ordering predicts up to 30% more events in the region E ν = [5,8]  Expected interaction rate of electron neutrinos and antineutrinos predicted by a NO over the rate predicted assuming an IO. Using the oscillation parameters in [3]. Because of the flux ratio νµ/νµ and the cross section difference, estimated to be 2.1 times larger for neutrinos than antineutrinos, more electron neutrino interactions are expected for a NO.
surements of the flux of atmospheric electron neutrinos thus provide suitable data for determining the NMO. The VLVNTs currently in operation are presented in details in the next section. With an energy threshold of about 15 GeV they operate in the saturated regime. They can measure muon neutrino disappearance as well as tau neutrino appearance, and thus θ 23 and |∆m 2 31 |. Measuring the sign of ∆m 2 31 , on the other hand, requires measuring differences in oscillation probabilities below this threshold (see Figs. 3 and 4). This is the main goal of the next-generation detectors discussed in Section V.

III. VERY LARGE VOLUME NEUTRINO TELESCOPES
A generic VLVNT is a three-dimensional array of photo-sensors detecting the Cherenkov light of charged particles produced after a neutrino interaction. The secondaries of neutrino interactions above a few GeV produce enough light so that they can be observed by sensors several meters apart. The spacing between the optical sensors define the energy threshold of VLVNTs, which is approximately 15 GeV in currently operating detectors.

A. VLVNTs in operation
The optical sensors of VLVNTs are deployed at depths of 1 km or more, in an optically transparent, naturally occurring medium. Sensors are laid out in lines or strings that are operationally independent. The spacing between sensors is uneven, being considerably larger in the x − y plane (in between lines/strings) than in the z plane. The sensors also have a preferred acceptance for light coming from below, although this might change for future detectors.
The neutrino telescopes currently in operation are Ice-Cube in Antarctica [36], ANTARES in the Mediterranean Sea [37], and the prototype of the Gigaton Volume Detector in Lake Baikal [38]. Both ANTARES and Ice-Cube have published studies of neutrino oscillations, and are therefore the only ones discussed in this review.

Detector design and layout
ANTARES is located between depths of 2025-2475 m, 20 km away from Toulon (French Riviera), in the Mediterranean Sea [37]. It comprises 885 optical modules (OMs) [39], distributed along 12 flexible lines. OMs are grouped in triplets, with 25 triplets per line. The distance between triplets is 14.5 m, and the separation between lines ranges from 60 to 70 m, as sketched in Fig. 6. Acoustic devices, tiltmeters and compasses are used to monitor the shape of the detector, which is influenced by sea currents.
IceCube is located at depths between 1450-2450 m at the geographic South Pole [36]. The in-ice part of Ice-Cube consists of 5160 downward-facing digital optical modules (DOMs) [40]. The detector has 86 strings, each holding 60 DOMs. Of these, 78 strings are arranged in a hexagonal grid with a typical distance of 125 m (horizontal spacing) and 17 m (vertical spacing) between DOMs. A sketch of the detector layout is shown in Fig. 7.
The lower center region of IceCube, from 1760 m down to 2450 m, houses DeepCore [41], a region of denser instrumentation (7 m DOM vertical spacing), where eight strings are separated by 40 − 70 m. Some 50 % of the PMTs in this region have 35 % higher quantum efficiency than the standard IceCube PMTs. The DeepCore fiducial volume used for data analysis is defined by a cylinder with a height of 350 m and a radius of approximately 150 m that starts below a dust layer, where the light transparency is reduced, as shown in Fig. 7. This volume, which corresponds to roughly 2.5 times that of ANTARES, encloses about 550 DOMs with reduced spacing and results in a threshold for detection and reconstruction of neutrinos of about 15 GeV.
The optical modules of both IceCube and ANTARES are glass spheres enclosing a ten-inch PMT, optical coupling gel, and a µ-metal cage for magnetic shielding. The IceCube OM digitizes the waveforms detected by the PMT inside the module before transmission [40], while the ANTARES OM keeps the readout to a minimum and only transmits the time and amplitude of a signal above threshold [42]. ANTARES optical modules have a baseline noise rate of 70 kHz at single photon level [43], while for IceCube (DeepCore) OMs the noise is 0.45 kHz (0.65 kHz) [44].

Optical medium and calibration
The optical properties of the medium affect the time of arrival and the number of detected Cherenkov photons. At the ANTARES site (salt water) the absorption length, which is 60 m for blue light (λ 470 nm) and 26 m for UV light (λ 375 nm), reduces the number of photons observed. The effective scattering length, which is 256 m for blue light and 122 m for UV light, is considerably larger than the spacing between sensors [45]. In the clear ice in which DeepCore is located the absorption length of UV light (λ 400 nm) is of the order of 200 m, which is larger than the spacing between sensors. The effective scattering length in the deep Antarctic ice is approximately 50 m, comparable to the string distance of DeepCore, thus significantly modifying the expected time of arrival of photons [46,47].
Water offers the advantage of being a homogeneous medium. Nonetheless, sea currents can deviate the detector lines so the position of the lines needs to be monitored constantly. This is achieved by combining acoustic triangulations with tilt and compass measurements yielding a precision better than 10 cm, which does not affect the angular resolution [48]. High sea currents can also trigger bioluminescence bursts that must be accounted for in the optical background simulation, in addition to the stable optical noise arising for 40 K decays. The latter can be used for determining the absolute detection efficiency of the optical modules.
In ice, the positions of the optical modules are fixed and known to within a few cm. Noise levels are constant and a hundred times lower than in salt water after the detector has stabilized. A disadvantage of using ice is that the medium is not homogeneous and its structure has to be modeled. This is particularly challenging in the immediate surroundings of the optical modules. Columns of the original glacier are melted to deploy the instrumentation. The refreezing process leaves behind clear ice near to the boundaries of the hole, and a cylinder of ice of about 10 cm in diameter with a high concentration of bubbles towards the center of the column. These changes in ice properties modify the DOM angular acceptance measured in the laboratory. Future detectors in ice will consider the possibility of de-gassing the water to avoid trapping air bubbles inside the hole ice, and with that reduce the impact of the medium.
The absolute optical efficiency of the optical modules as well as their angular acceptance must be determined in situ after deployment. ANTARES and IceCube use both controlled light sources and minimum ionizing muons to calibrate the efficiency and timing accuracy of their optical modules [49][50][51]. Relative arrival times are known with a precision better than 3 ns and 1.5 ns for IceCube [40] and ANTARES, respectively.

B. Neutrino interactions
The dominant neutrino interaction for most of the energy range that VLVNTs can access is neutrino-nucleon deep inelastic scattering (DIS), with other processes being only a sub-dominant contribution. Nonetheless, below 15 GeV, the region of interest to search for matter effects in neutrino oscillations and the NMO, quasi-elastic scattering and production of resonances compete with DIS processes. Figure 8 shows a calculation of the competing νN cross sections around the GeV region, together with the data available.
Most of the knowledge of neutrino-nucleon cross sections between 1−15 GeV comes from bubble chambers or spark chamber detectors which collected comparatively small data samples. Thus, the constraints on the models that describe them are rather weak [52]. The uncertainty with the largest impact on the neutrino cross sections for quasi-elastic and resonant interactions, which changes them by up to 40%, is the value of the axial mass that effectively describes the nucleon form factor and has an estimated error of 15%-25% [52,53]. DIS interactions in the crossover region have a small momentum transfer. Non-perturbative QCD calculations are required [54], and the estimated errors are as well of the order of 20% [55].
Deep inelastic scattering accounts for 90% or more of the total cross section of neutrinos and antineutrinos above an energy of roughly 12 GeV, as shown in Fig. 8. DIS in the perturbative regime is comparatively better understood than the processes discussed so far, with uncertainties coming mainly from the determination of the parton distribution functions (PDFs) of the nucleons. The uncertainties on the PDFs change the total cross section by 5% or less [55].
At these energies the neutrino-nucleon DIS charged current (CC) cross section is quasi-independent of the inelasticity y (y = 1 − E lepton /E ν ) of the interaction, while for antineutrinos the cross section is accompanied by a factor (1−y 2 ), which suppresses kinematic configurations where the hadronic part of the interaction takes most of the energy. The inelasticity dependence makes the total νN cross section about one half of that of νN . PDF uncertainties impact the inelasticity dependent cross section by only a few percent [56].
While the neutrino-nucleon DIS CC cross sections for ν e and ν µ are equal, the ν τ N one is suppressed due to the mass of the tau lepton. It is only at E ν ∼ 40 GeV that the cross section reaches half of the value of the other neutrino flavors [57].
In neutral current interactions (NC) one or several  [52] for details on the experiments contributing to the data points and [58] for a description of the model used). The contributing processes in this energy region include quasi-elastic (QE) scattering, resonance production (RES), and deep inelastic scattering (DIS). Taken from [52].
hadrons are produced, initiating a hadronic shower. In charged current (CC) interactions a hadronic shower is also present, but now the neutrino transforms into a charged lepton. Electrons and taus also initiate a shower of particles after they are produced 3 . Muons, on the other hand, travel practically undisturbed and lose energy at a quasi-constant rate. For muons passing through water dE/dx ∼ 0.25 GeV/m up to a few hundred GeV [1]. Once they travel distances comparable to the detector spacing they can be identified, and charged current ν µ interactions can be tagged.

C. Event reconstruction
Neutrino interactions are reconstructed using the number of photons recorded by the optical module (or time over threshold), as well as the time at which they arrive. For the energy range under consideration, the most general hypothesis is an interaction which produces a hadronic shower (all interaction types), an electromagnetic shower (ν e CC) or a long range muon (ν µ CC). The direction in which these particles are produced is reconstructed from the arrival times of the emitted photons. The Cherenkov light of muons is produced almost perfectly in a cone. The light coming from the cascade is also beamed in the Cherenkov angle, but the smearing due to multiple particle contributions to it is larger, which degrades the achievable precision of directional reconstructions. This smearing effect is stronger for hadronic showers.
The energy reconstruction of showers is primarily given by the number of photons detected from a given interaction, and its accuracy depends mainly on the reconstructed position of the interaction vertex. To estimate the energy an assumption has to be made on whether the shower is hadronic or electromagnetic. The energy of muons can be estimated by the observed range in the detector.
In principle it is possible to fit the directions of both cascade and track components in an interaction. The sparse instrumentations of the detectors, however, make it challenging. In the simplest approach, tracks and cascades are assumed to be collinear.

D. Simulation tools
The measurement of diffuse fluxes in VLVNTs, such as the one required to determine oscillation parameters, rely fully on the correct modeling of the experimental set-up. Atmospheric muons, the leading source of background, are simulated in IceCube using full showers and parameterizations obtained from CORSIKA [63]. ANTARES uses the MUPAGE program, which produces muons based on a parameterization tuned to MACRO data [64].
Neutrino interactions in IceCube are simulated using the GENIE package [65] (E ν ≤ 200 GeV) and NuGen/ANIS [66] (E ν ≥ 50 GeV). Besides GENIE, ANTARES uses an in-house neutrino generator based on LEPTO [67] for the full energy range, with the PYTHIA package [68] handling the hadronization processes. The neutrinos produced are then weighted to match the flux predictions of the Honda and/or Bartol groups [21,69].
The propagation of short-ranged particles produced in the interaction is done, both in IceCube and ANTARES, using the Geant software [70] as basis. Parameterizations of the light yield of these particles are produced by both experiments, and used to obtain the detector response to high energy hadrons, electrons and photons [71], while   [4,9,37,[59][60][61][62]. Expected neutrino events quoted from published results of νµ disappearance at analysis level (note that for VLVNTs this number can vary significantly depending on the studied range in energy, zenith angle and topology). COH refers to coherent pion production. For details on the other interaction channels and energy ranges see Fig. 8.
low energy hadrons (E ≤ 30 GeV) are propagated individually. Muons are propagated using code optimized for simulation of long ranged leptons, namely MUSIC [72] and MUM [73] in ANTARES, and MMC [74] in IceCube. The Cherenkov photons produced during the propagation of charged particles are individually traced through the ice in IceCube/DeepCore, while ANTARES uses expectation from tables. Low energy future projects (Section V) foresee using only individual tracing to assure that the optical properties of the medium are included in detail. After the photons are propagated, the response of the optical module is recreated, and events from simulation and experimental data are treated equally.

E. Large statistics vs. precise reconstruction
The current knowledge on the oscillation parameters in the atmospheric sector comes from experiments which differ from VLVNTs substantially: they are Super-Kamiokande [75], T2K [76], MINOS [62], which is no longer in operation, and the recently commissioned NOvA [77]. Table I contains a qualitative comparison of the detectors and neutrino sources used by these experiments and VLVNTs.
Super-Kamiokande, which also measures atmospheric neutrinos, has about twenty (ten) times the number of optical sensors as DeepCore (ANTARES), separated by a few cm, placed on a cylindrical tank with a diameter similar to the inter-string distance in ANTARES/DeepCore. Neutrinos are detected using the rings produced after the Cherenkov light of the charged products of the interaction hits the walls of the detector. Muons, elec-trons and pions can be identified by the differences in the ring pattern they produce. Because of its considerable smaller size and the steepness of the spectrum of atmospheric neutrinos, its operating energy is lower than that of VLVNTs.
Long baseline experiments, such as T2K, MINOS and NOvA, use neutrinos from particle accelerators and have near and far detectors. While T2K uses Super-Kamiokande as a far detector, MINOS and NOvA follow an experimental set-up where the far detector is smaller than Super-Kamiokande, but is more densely instrumented, can be magnetized, and observes the path of individual particles coming from a neutrino interaction. These set-ups benefit from their controlled neutrino source and detailed event reconstruction. Unlike the case of atmospheric neutrino experiments, long baseline experiments have a unique baseline and cover a narrow energy range, allowing for better precision but also limiting the L/E region that they can access. It should also be noted that, as stated in Section III B, the poor knowledge of neutrino interactions at energies of a few GeV introduces significant uncertainties in the data analysis of long baseline oscillation experiments.
VLVNTs have become competitive with accelerator based experiments thanks to the possibility of observing multiple combinations of baseline and energy (L/E), and with Super-Kamiokande because VLVNTs can collect large event samples and in an energy range where most events are DIS which can be modeled with high accuracy. The sparse instrumentation does not permit observation of small details of the interaction, but in the same way reduces the impact from uncertainties in the hadronization processes, one of the leading systematic uncertainties for MINOS [78] and T2K [4]. Reconstruction accuracy and proper handling of systematic uncertainties are the most important points to consider for precision measurements with VLVNT.

IV. NEUTRINO OSCILLATION MEASUREMENTS FROM RUNNING VLVNTS
The ANTARES and IceCube collaborations have published measurements of oscillations studying the muon neutrino disappearance channel. Above 15 GeV, where these detectors operate, muon neutrinos oscillate into tau neutrinos, following Eq. 8. Signal neutrinos, i.e. ν µ interacting via CC with E ν ∼ 25 GeV, are typically recorded by a handful of optical modules both for ANTARES and IceCube's DeepCore. The events develop over a distance of order of 100 m, and thus can be fully contained in both detectors.
The measurement of neutrino oscillations in VLVNTs follows a general strategy which begins with the reduction of the dominant sources of background, i.e. atmospheric muons and pure noise. Straight cuts are applied on variables of which the distribution for neutrinos differs from that of background sources. They generally aim for a neutrino purity higher than 95%.
For the currently published results of both experiments, the presence of a muon in a neutrino interaction is required for an event to be selected for analysis. The analyses are done by comparing the histograms of data and simulation as a function of the reconstructed variable(s) used. The simulation is modified by the physics parameters of interest, θ 23 and ∆m 2 32 , and by nuisance parameters which absorb the systematic uncertainties involved in the measurement. Errors are derived from a scan of the likelihood landscape, and/or directly using a χ 2 approximation.
The results of ANTARES and IceCube that have been made public until now use only events coming below the horizon. ANTARES removes the down-going region because it is dominated by atmospheric muons. IceCube uses the instrumentation outside DeepCore to veto atmospheric muons, nevertheless the contribution of these muons in the down-going region is still significant, so the region is also removed from analysis. This situation is different for Super-Kamiokande, where events from the entire zenith range are used in oscillation studies and top-down ratios are used to reduce uncertainties. Ongoing studies within IceCube are exploring the possibility of using neutrinos coming from above the horizon in future results [79].

A. First measurements of oscillations from ANTARES
The ANTARES collaboration presented the first results on the study of neutrino oscillations from a VLVNTs [59]. The analysis relied on the muon track reconstruction described in [80], which fits the depth at which the Cherenkov cone of light arrives at the OMs as a function of time. This corresponds to a hyperbola of which the orientation of the asymptotes depends on the zenith angle. An algorithm that searches for these patterns, without assuming any knowledge on the arrival angle of the emitter, was implemented. The algorithm is capable of rejecting noise hits and keeping events down to energies of 20 GeV (R µ =100 m) with photons in a single line, and 50 GeV (R µ = 250 m) in multiple lines. Misreconstructed muons that appear up-going are removed by selecting only events which have a good fit quality. This cut also effectively reduces the contribution of NC interactions from all flavors, and ν e CC interactions. The median zenith angle resolution with respect to the neutrino direction of single-line events is 3.0 • , and it reduces to 0.8 • for multi-line events. The energy of the neutrino is estimated solely by the muon range, resulting in a lower limit to the neutrino energy, where The analysis is done by comparing data and simulation as a function of E reco / cos θ reco by means of a χ 2 , combining single-and multi-line selections. Only events below the horizon (cos θ reco < −0.15) are considered. Systematic uncertainties are implemented using two normalization coefficients, for single-and multi-line events, as pull factors in the χ 2 following the method presented in [81]. These factors absorb the effects of changes in the average quantum efficiency (±10%), optical properties of sea water (±10%), the spectral index of atmospheric neutrinos (±0.03), and disagreements between data and simulation during the selection (varying cut values). The overall normalization of the ν µ flux and detector efficiency are left unconstrained.
The data analyzed were taken between March 2007 and December 2010, corresponding to a detector live time of 863 days. A total of 2126 neutrino candidates were selected. The measured oscillation parameters, which were found to be compatible with the world's average, are indicated in Fig 13. Data and simulation were in good agreement, as it can be seen in Fig. 9, which results in a χ 2 /NDF = 17.1/21. The case of no oscillations could be rejected at the 3σ confidence level. The ANTARES collaboration will proceed to an updated analysis of this kind with the full data sample collected until the end of the data taking, circa 2017.

B. First measurements from IceCube DeepCore
To this date, IceCube has reported results of four neutrino oscillation analyses of the low-energy DeepCore data. The selection, reconstruction and analysis methods have been refined in each step. The low-energy data for all studies comes from the DeepCore filter and trigger [41]. The main source of background at this stage are triggers due to sensor self-noise and atmospheric muons. The instrumentation outside the fiducial volume of Deep-Core (see Sec. II-A and Fig. 7) is used to tag atmospheric muons. Low energy neutrino interactions are required to start within the DeepCore fiducial volume, while no requirement is imposed for full containment.
Systematic uncertainties are accounted for using additional parameters which modify the expected number of events. An energy dependent term (E −γ , γ ± 0.05) and a free overall normalization absorb total cross section uncertainties, and the uncertainties on the spectral index of the neutrino flux. The electron neutrino flux is varied by ±20% around the predicted value. The cosmic ray models which predict the cosmic muon contamination are varied to obtain a robust estimate. The effects of changing the optical description of the pristine ice, as well as the re-frozen ice around the DOMs, are studied by producing multiple simulation sets.
The initial three oscillation studies from DeepCore, presented first herein, were restricted to a single year of detector live time. Two used a partial configuration (IC79, two DeepCore strings missing) and one used the full detector (IC86). The first analysis [82], from here on IC79-A, used a DeepCore low-energy sample where the effect of oscillations is expected (E ν < 100 GeV, 719 events) and an IceCube high-energy sample, where oscillations play no role, to constrain flux and detection uncertainties (E ν 1 TeV, 39638 events). The measurement was done by analyzing the distribution of events as a function of zenith angle in the low-energy sample (see Fig 10). The zenith angle of both samples was estimated using the muon track reconstruction described in [83]. Atmospheric muons were mainly removed by reconstructing all events as up-going, and making cuts on parameters related to the quality of the reconstruction (without muon tagging).
The data were analyzed using a χ 2 optimization with pulls, also following the method in [81]. The results ob-tained for the atmospheric oscillation parameters were compatible with contemporary global fits [84], although the errors were a factor 4 to 9 larger (see Fig. 13).
FIG. 10. Data and simulation expectation at world average oscillation parameters (in black) and the case of no oscillations (in red) for the low energy sample of IceCube's IC79-A analysis. Systematic uncertainties are split into a fully correlated part (hatched bands) and uncorrelated part (shaded bands). From [82].
Two subsequent analyses of the data, from here on IC79-B and IC86-A, created new event selections based on the rejection of atmospheric muons by using the veto, separating the background rejection from the reconstruction of events [85,86]. In both cases only the low-energy DeepCore data were analyzed.
The data used for IC79-B were acquired during the same period of time as for IC79-A, however, due to the change in the selection of events the final sample studied was a factor 10 larger. The zenith angle of events was reconstructed with a similar method as in IC79-A [85]. A second observable, the reconstructed muon range L muon [87] was used as an energy proxy and the data were analyzed as a function of both observables. The ratio of events with respect to the no oscillation scenario, together with the best fit, are shown as a function of reconstructed L osc /L reco in Fig. 11, where L osc is the distance the neutrino traveled and L reco is the reconstructed length of the muon produced in the interaction. The best fit and estimated errors of this method were similar to those of IC79-A.
The first analysis of data from the full detector configuration [86], IC86-A, was performed using a selection of photons and event reconstruction based on the method published by ANTARES [80]. The selection of photons was modified to remove multiply scattered photons instead of noise. Unscattered, or direct, photons were identified by restricting their possible arrival times to those given by the hyperbolic pattern that Cherenkov light produces as a function of time as it crosses a string. About 70% of the neutrino interactions which trigger the detector do not have a clear core of direct photons, and thus are removed.  11. Ratio of the distribution of oscillation length over reconstructed track length to the no-oscillation hypothesis from simulation in the IC79-B analysis. The best fit is also shown. From [85].
The direct photons found are used to fit track and cascade hypotheses. The zenith angle from the track fit was used as an observable, and the ratio of the χ 2 of the track and cascade fits was used to separate track-like from cascade-like events. An estimator of the total energy of the neutrino was also implemented, which takes the muon range estimator from IC79-A and also fits a hadronic cascade at the vertex.
In IC79-B and IC86-A the data were analyzed using a likelihood optimization with nuisance parameters to account for systematic uncertainties. For IC86-A, uncertainties related to the detector were also included as nuisance parameters. Simulation sets with varied detector settings were produced and interpolated at the final level of the analysis, allowing the fitter to make arbitrary modifications to them.
In similar livetime as IC79-A and IC79-B, IC86-A selected 1487 neutrino events for analysis. While the best fit obtained was in agreement with the other results, the error in ∆m 2 32 was reduced by about 20% with respect to IC79-A, while maintaining a similar precision on sin 2 2θ 23 . Figure 12 shows a comparison of data and best-fit simulation in projections in energy of the twodimensional histogram used in the analysis. A comparison of the confidence regions in sin 2 θ 23 and ∆m 2 32 of the single year analyses of IceCube DeepCore, together with the result from ANTARES, is shown in Fig. 13.

C. Precision measurements with IceCube DeepCore
The latest result from IceCube DeepCore [9] is an update to the IC86-A analysis introduced before, now with almost a thousand days of detector live time. The measurement demonstrates the potential for VLVNTs to become relevant experiments in the field of neutrino oscillations.
While the analysis strategy is still to focus on the selection on clear tracks, for which a core of direct photons can be identified, three large improvements are introduced, namely: • An optimization of the event selection, which results in 40% more events.
• The cosmic muon background is derived from data (tagged muons), avoiding the need of computationally expensive model-dependent simulation.
• An improved estimator of the energy deposited at the interaction point, which reduces the error on the total neutrino energy by more than 30% at 20 GeV.
A demonstration of how the data-derived background is used can be seen in Fig. 14, where the distribution of events as a function of reconstructed zenith angle at the final level and two earlier stages of the event selection is shown. At each step the cosmic muon background is more strongly suppressed. The contribution of atmospheric muons in the down-going region can be seen at all steps, including the final sample to be analyzed.
For their IC86-B result, the IceCube collaboration has expanded the list of possible sources of uncertainties considered. Non-DIS events are a non-negligible fraction of the sample at E reco ≤ 20 GeV, and additional cross sections uncertainties on these interactions (about 20%) were also included. A possible shift of 5% in the energy scale of hadronic showers was also taken into account.
In 950 days of live time, a total of 5174 events were observed, while 6830 were expected without oscillations. Note that the energy range of the search was reduced in comparison with IC86-A to E reco = [7,56] GeV. The data were analyzed in a full three-neutrino oscillation formalism, including the effects induced by matter as neutrinos cross the Earth. The parameters that best describe the data, assuming a normal mass ordering, are sin 2 θ 23 = 0.53 +0.09 −0.12 and ∆m 2 32 = 2.72 +0.19 −0.20 × 10 −3 eV 2 . No significant preference was found for either the normal or inverted mass orderings. Purely statistical uncertainties are +0.06 −0.08 for sin 2 θ 23 , and +0.14 −0.15 × 10 −3 eV 2 for ∆m 2 32 , from which it is deduced that statistical and systematic uncertainties have an almost equal impact on the result.
Data and simulation are in good agreement, with a χ 2 /NDF= 54.9/56 for the energy-zenith angle histogram used in the fit. Figure 15 compares the L reco /E reco distributions of data and best fit simulation, where the agreement can be observed 4 . The 90 % confidence contours on the atmospheric oscillation parameters obtained are shown in Fig. 16, together with the results from the other experiments leading the field.  [59] and Ice-Cube's single year measurements [82,85,86], compared to earlier measurements by MINOS [61], T2K [89] and Super-Kamiokande [90]. Taken from [91].
The results from VLVNTs will be further improved by adding statistics to the analyzed data sample and refining the reconstruction methods. However, the most decisive improvements will come with the construction of the next-generation VLVNTs presented in the next section.

V. NEUTRINO OSCILLATIONS WITH THE NEXT GENERATION OF VLVNTS
After the measurements from ANTARES and Ice-Cube/DeepCore in the atmospheric sector, the next goal of VLVNTs is to further decrease the energy threshold below the 15 GeV domain in order to improve the sensitivity to the PMNS matrix elements and determine the NMO. Measuring the neutrino mass ordering is the main objective of the forthcoming ORCA (Oscillation Research with Cosmics in the Abyss) [94] and PINGU (Precision IceCube Next Generation Upgrade) [95,96] detectors as part of the KM3NeT [97] and IceCube Gen-2 [98] infrastructures, respectively.

A. Design of future detectors
Both ORCA and PINGU will be more densely equipped than the currently operating detectors and should reach several megatons in instrumented volume. Their concepts are similar in many ways, with the most significant differences coming from the detection medium, the proposed detector layout and the (default) optical module design.

Hardware and detector geometry
The PINGU optical module will most likely be a simplified and modernized version of that of IceCube, which has demonstrated its stability and reliability over almost ten years of operation. The PINGU DOM design removes components that are no longer required, such as the local coincidence logic and the multiple amplification modes, while providing a larger dynamic range than the original IceCube DOM and improved time resolution of 2 ns [96]. A schematic view of the IceCube and PINGU (Gen2) DOMs is shown in Fig. 17. By maintaining the basic Ice-Cube design, the PINGU DOM minimizes risk and cost. The ORCA optical module will follow the KM3NeT design [97] with each DOM housing 31 small (3") PMTs arranged in a 17" glass sphere together with the associated electronics, as can be see from Fig. 18. This design offers the possibility of creating coincidences within the OM to suppress the large 40 K decay background as well as the thermal noise of the PMTs. The orientation of the PMTs within the OM is also used in the reconstruction of events, although not yet at its full potential. A single sphere houses three to four times the photo cathode Ratio to no osc. Data are compared to the best fit and expectation with no oscillations (top) and the ratio of data and best fit to the expectation without oscillations is also shown (bottom). Bands indicate estimated systematic uncertainties. From [9].
area of an ANTARES OM with an almost uniform angular coverage, improving the cost effectiveness by a factor four. Several prototypes of such a multi-PMT OM have been successfully tested in-situ [99]. The final layouts of ORCA and PINGU are still under optimization 5 . The current benchmark geometries used for establishing the detector performances consist  [61,92,93]. The log-likelihood profiles for individual oscillation parameters are also shown (right and top). A normal mass ordering is assumed. Updated from [9]. of 40 (115) strings with a horizontal spacing of ∼20 m for PINGU (ORCA). The vertical spacing is set to 6 m for ORCA and 3 m for PINGU. While a PINGU string will hold up to 96 DOMs, there are 18 DOMs in a default ORCA string. The maximum number of DOMs that a PINGU string can hold is given by the mechanical constraints of the down-hole cable and the appearance of shadowing effects, while for ORCA the constraint comes from the launcher vehicle (a large spherical frame in which the DOMs slot into dedicated cavities) used for string deployments. The separation between the sensors of both detectors is smaller than the absorption and scattering lengths of their respective media, making the optical properties of ice and salt water less relevant than for ANTARES and IceCube/DeepCore. The footprints of the ORCA and PINGU detectors are shown in Figure 19. The instrumented mass of both detectors is of order 3.5 to 4 Mt, and their effective masses reach the same value for neutrinos of energy above 10 GeV. While the PINGU extension is foreseen to be embedded inside the current IceCube/DeepCore detector (which will be used for background vetoing), the ORCA detector will be located around 10 km west from the ANTARES site, at a depth of 2475m.

Costs and timescale
PINGU estimates a cost of 48M$ for hardware and 23M$ for logistics [100]. The estimated cost of ORCA is 40Me. Funding request processes are currently driving the possible time line of the projects.
PINGU will be built as part of the IceCube Gen-2 project. From a technical point of view the installation of the detector at South Pole could start by the end of 2020 [101]. Based on the experience gained with the IceCube, the deployment is expected to take only three years. The first construction phase of ORCA, a demonstrator array of 6-7 strings (already funded), started in late 2014 with the deployment of the main electro-optical cable, followed by the deployment of a junction box in April 2015. The demonstrator is expected to be deployed by the end of 2016, and will be used to carry out studies of detectorrelated systematic effects and event reconstructions. In an optimistic case, the deployment of the full detector case could happen by 2020. Both PINGU and ORCA plan to take data during their construction phase.

B. Projected performance
The determination of the NMO, the main physics goal of these projects, relies on a detailed analysis of deviations of the order of ∼10% and ∼ 30% in the rates of detected atmospheric muon and electron neutrinos (see Figs. 3,4,5) as a function of energy and arrival zenith angle. Therefore, the key parameters that characterize the potential of a detector are its effective mass, the energy and zenith angle resolutions achievable, and its particle (mis-)identification capabilities. In the following discussion the latest, preliminary, studies from ORCA [94,102] and PINGU [96,103] are presented.
These studies are based on full Monte Carlo simulations adapted from IceCube and ANTARES. All ORCA results account for an optical background induced by 40 K decays of 5-10 kHz per PMT and a time-correlated hit rate of 500 Hz per OM (two coincident hits in different PMTs inside the same OM). Since PINGU DOMs will follow closely the design used for IceCube, the typical in situ behaviour of the IceCube/DeepCore DOMs, with a noise rate of 650 Hz, is used in the simulations.
The published results of ANTARES and IceCube have so far focused on ν µ disappearance, and therefore only selected events where a muon was observed. The sensitivity to the NMO, on the other hand, also comes from oscillations that involve ν e . It is therefore useful to detect all neutrino flavors, placing them in two categories depending on their topology: tracks and cascades (see Sec. V B 2).

Reconstruction of tracks and cascades
Track-like events are those where a muon is observed coming out of the interaction vertex. Track-like topologies are CC ν µ interactions as well as the ν τ CC interactions when the decay of the tau lepton produces a muon. The cascade-like topologies are CC ν e interactions, CC ν τ interactions without a muon in the final state, and NC interactions from all flavors. Independent studies indicate that after accounting for reasonable detector resolution effects, the cascade channel provides more sensitivity to the effects of the NMO. Note, however, that the two channels are complementary as track-like events can provide better precision in sin 2 θ 23 . It is consequently important to be able to distinguish the two topologies with high efficiency and purity.
The event reconstruction in PINGU is a simultaneous global likelihood fit of the interaction vertex position and time, the zenithal and azimuthal angles, the energy of the cascade at the vertex, and the length of the daughter muon track. The event hypothesis assumes that tracks and cascades are collinear. The likelihood is calculated using the time of arrival of single photons and the expected noise in the time windows analyzed. The expectations for minimum ionizing muon tracks and electromagnetic cascades needed for the likelihood are stored in tables, obtained from direct simulation of particle and photon propagation, as it is already done for IceCube [51]. An event is reconstructed by comparing photon expectation for a given event hypothesis to the photons observed. All the DOMs in PINGU, as well as those in IceCube/DeepCore, are used in the reconstruction [96].
Fitting eight parameters at once while simultaneously looking up expectations from tables makes the reconstruction CPU intensive, but in return it provides robust results and similar resolutions for track-like and cascadelike topologies. While it would be possible to use the information provided by this reconstruction to obtain an estimate of the inelasticity of the event, this has not been explored so far. Energy and zenith angle resolutions for different interactions are shown in Figs. 20 and 21, together with those obtained by ORCA with the methods explained hereunder.
ORCA uses two distinct algorithms for tracks and cascades. The track reconstruction is directly adapted from the main reconstruction of ANTARES [104] and focuses on the muon direction using the combined information of the PMT spatial positions and the Cherenkov photon arrival times. The neutrino energy estimation is mainly given by the reconstructed muon track length, which is complemented by the number of hits used in the track reconstruction algorithm. Muon tracks produced in neutrino interactions at E ν ≥15 GeV are not always fully contained, which turns the estimate into a lower limit above these energies, as shown in Fig. 20. The time residuals under a spherical emission profile (shower-like) or according to a Cherenkov cone (track-like) are used to obtain sensitivity to the inelasticity in the track channel.
The cascade reconstruction in ORCA takes advantage of the long scattering length in sea water, which preserves the structure of the Cherenkov light cone, and tries to identify the leading lepton in the cascade. An example of the distribution of the expected number of photons as a function of emission angle for different inelasticity intervals is shown in Fig. 22. A peak is always visible at the Cherenkov angle (42 • ), whose height with respect to the off-peak region depends on y. Cascades are reconstructed in two separate steps using maximum likelihood fits. First the interaction vertex is obtained with a resolution of about 0.5-1 m by an algorithm based on hit time residuals. It is then followed by a fit of the direction, energy and inelasticity of the event. The performances of the cascade reconstruction are summarized in the Figures 20 and 21.
In ORCA the inelasticity of about 60% of the tracks with true y ≤ 0.25 or y ≥ 0.75 is reconstructed correctly; the accuracy of the inelasticity estimator of cascades is slightly worse. The inelasticity could be used for potential statistical separation between neutrinos and anti- For ORCA individual resolutions for neutrinos and antineutrinos are shown, while a mixture of both is given for PINGU.
Resolutions are better for anti-neutrinos than for neutrinos due to the smaller average inelasticity, leading to a smaller intrinsic scattering angle between the neutrino and the leading lepton. Values taken from [94,96].
neutrinos, which can be exploited for the mass ordering measurement [105]. It can also be tested to separate charged current interactions from neutral current interactions. While both PINGU and ORCA are studying this possibility, inelasticity estimates are not yet part of the current analyses that are discussed in the following sections.

Particle identification and background rejection
VLVNTs measuring atmospheric neutrinos should be able to identify and reject atmospheric muons, the largest source of background, and differentiate between events with track-like and cascade-like topologies. PINGU plans to tag atmospheric muons following the strategy developed in DeepCore, i.e. using the outer detector strings to identify particles that enter the fiducial volume, and restricting the analysis to starting and up-going events (see [9] and Fig. 14). The cosmic muon background is expected to be on the level of a few percent, similar to DeepCore. Event reconstruction and selection in PINGU do not rely on direct hits, the single largest impact on signal efficiency in the latest DeepCore results. Signal efficiency in PINGU, therefore, is expected to be minimally affected by background rejection and reconstruction methods, and largely defined by the number of photons observed from an interaction.
The ORCA detector does not rely on an outer detector to tag muons. Current analyses reduce the impact of these muons by selecting only up-going events and rejecting the mis-reconstructed ones using variables such as their reconstruction quality and the position of their reconstructed interaction vertex. The topology of neutrino interactions, track-like or cascade-like, is identified using the distribution of hit time residuals, distances between reconstructed vertices at various reconstruction steps, the quality of the reconstructions, and topological variables, among others. A single multivariate method which incorporates the parameters listed above is applied to the data, and classifies events as tracks, showers or atmospheric muons. The procedure achieves ∼1% muon contamination in the final sample without a severe signal loss.
The approach followed in PINGU to separate tracks from cascades also uses a multivariate method with vari- ables describing the reconstruction quality of the event under the track vs cascade hypothesis, as well as the reconstructed muon track length as input. Figure 23 shows a comparison of the performance of these classification methods for neutrino interactions around the energies relevant for mass ordering measurements for PINGU and ORCA. Current methods differ at low energies, with PINGU showing a bias towards classifying low energy tracks as cascades and ORCA exhibiting the opposite behavior. Above 10 GeV both classification schemes result in a similar outcome. The results suggest that the behavior of the particle identification algorithms at low energy can be tuned for optimizing sensitivity to the NMO measurement. In both cases, the final performances are subject to further optimization.

C. Physics potential and systematics
The preliminary performances described above are used by the PINGU and ORCA collaborations as inputs to estimate the confidence level with which the projected experiments will be able to reject a given NMO. This is done by drawing several thousands of pseudoexperiments generated under each mass ordering hypothesis, as outlined in [106]. The analysis is conducted by comparing the two-dimensional histograms of pseudodata and simulation as a function of the reconstructed energy and zenith. The pseudo-data sets are generated varying different parameters, such as the values of the mixing angles, in order to study the impact of degeneracies in the measurement.
A full log-likelihood ratio (LLR) method is used by both collaborations to report their expected sensitivity. In this method each pseudo-experiment is analyzed by performing a log-likelihood fit with the oscillation parameters as free parameters (mostly θ 23 , ∆m 2 32 and θ 13 ) and assuming both hierarchies in turn. Sources of systematic uncertainty are incorporated as additional parameters in the fit (see Table II). As such methods can be quite CPU expensive, in particular when studying various sources of systematics, the PINGU collaboration also implemented a simplified ∆χ 2based approach. This method is a parametric analysis based on the Fisher information matrix, which relies on the partial derivatives of the event counts in each bin with respect to all parameters under study. Inverting the Fisher matrix yields the full covariance matrix between the parameters. The covariance matrix of the mixing angle θ 23 is calculated at several values to overcome the limitations of the method. The results obtained with the Fisher matrix are in agreement with the LLR method, and are also used to report the projected sensitivity of PINGU.
The parameters of the fits performed by ORCA and PINGU, presented in Table II, are the oscillation parameters of interest plus a set of parameters related to uncertainties on the detection process, neutrino fluxes, cross sections and the remaining oscillation parameters. The oscillation parameters, in particular θ 23 , have the largest impact on the achievable precision. The overall normalization has the second largest impact on the precision. This absorbs uncertainties on the efficiency of the detector, the absolute atmospheric neutrino flux and interaction cross sections. PINGU has recently studied uncertainties on the neutrino flux by using a more refined de-scription, which involves a set of 18 parameters [23]. The impact found was a reduction of the three year sensitivity by 0.2σ [107] (not yet included in Fig. 24). Cross sections have been also studied in more detail by modifying the six most relevant parameters of the model implemented in GENIE. The reduction in sensitivity was found to be negligible. Studies within ORCA and PINGU have tested the impact of δ CP and found an additional reduction of up to 0.5σ at the three-year benchmark [79,94]. Note that all figures in this review do not include this effect.
The LLR (and ∆χ 2 for PINGU) resulting from fits to the pseudo-experiments are used to calculate the separability of the two possible mass orderings. The median (i.e with 50% statical power) sensitivities to the NMO are shown in Fig 24 (top) after 3 years of data taking. The results are obtained by fixing δ CP to zero, and are shown as a function of θ 23 . Both collaborations observe that constraining θ 23 to either octant while doing a fit artificially increases the sensitivity to the NMO, thus the parameter is left unconstrained in these studies.
Though ORCA and PINGU sensitivities should be compared with caution, as the various inputs are slightly different, both studies find a better sensitivity to the NMO for a true value of θ 23 in the second octant in the case of normal mass ordering. For the case of inverted ordering, the sensitivity has a much weaker dependence on the value of θ 23 . The consistency of the two results is encouraging, as they have been obtained with completely independent analysis chains.
The expected improvement in sensitivities with running time, which does not yet include the effects of δ CP nor the reconstructed inelasticity, are shown in Fig 24  (bottom). Once more, the discrimination power of both detectors is comparable.
The identification of the mass ordering devised by both collaborations also produces a measurement of θ 23 and the absolute value of the atmospheric mass splitting. Projections of the sensitivity to sin 2 θ 23 have a strong dependence on the assumed true values. For sin 2 θ 23 = 0.45 both PINGU and ORCA expect to achieve errors of the order of 0.05 after three years of operation. The precision achievable on the absolute value of the mass splitting is roughly independent of the true value, and the expected error on the measurement for both projects is about 0.05 × 10 −3 . Both experiments are expected to produce measurements with better precision to those projected for NOvA and T2K by the year 2020.
The results shown in the present paper are a compilation of the most recent, publicly shown projections of both collaborations, and include most leading systematics effects [79,94,103,108]. Recently a thorough study of the interplay between the oscillations parameters has been reported in [109], consistent with the recent results from ORCA and PINGU. The authors also introduced uncertainties in the estimated energy and zenith resolutions, as well as additional (conservative) uncorrelated uncertainties. Their results show that after 5 years of data taking, the loss in sensitivity ranges from 24% to Median significances as a function of time for the benchmark detectors described in the text. The oscillation parameters injected are close to those found in [3] (θ23 = 42 • for a NMO, θ23 = 49 • for an IMO). From [94,103].
40% under pessimistic assumptions (such as fully uncorrelated errors in each analysis bin), thus leaving room for a measurement of the NMO by ORCA and PINGU on a reasonable time scale.

VI. SUMMARY
Atmospheric neutrinos are a versatile tool to study neutrino oscillations. This naturally occurring beam covers baselines as large as Earth's diameter, and has an energy range which spans over the regimes of oscillations in vacuum, and with resonant and saturated matter effects. Current VLVNTs, ANTARES and IceCube, can detect neutrinos in the latter regime, and have already produced measurements of the atmospheric oscillation parameters, θ 23 and |∆m 2 32 |. Constant improvements in the understanding and modeling of the detector and media, as well as more sophisticated data analysis techniques have lead to promising results, which have started to become comparable with those of other, more mature experimental set-ups.
Proposed VLVNTs, ORCA and PINGU, aim to lower the energy threshold and access the resonant regime, with the goal of measuring the sign of ∆m 2 31 and completely determining the neutrino mass ordering. While both projects are on the way of optimizing their detector geometries and/or analysis techniques, current studies are nevertheless mature and indicate that they could provide a significant measurement (≥ 3σ, depending on the true value of θ 23 ) of the neutrino mass ordering after 3 − 4 years of operation.