Intrinsic limits on resolutions in muon- and electron-neutrino charged-current events in the KM3NeT/ORCA detector

Studying atmospheric neutrino oscillations in the few-GeV range with a multimegaton detector promises to determine the neutrino mass hierarchy. This is the main science goal pursued by the future KM3NeT/ORCA water Cherenkov detector in the Mediterranean Sea. In this paper, the processes that limit the obtainable resolution in both energy and direction in charged-current neutrino events in the ORCA detector are investigated. These processes include the composition of the hadronic fragmentation products, the subsequent particle propagation and the photon-sampling fraction of the detector. GEANT simulations of neutrino interactions in seawater produced by GENIE are used to study the effects in the 1 - 20 GeV range. It is found that fluctuations in the hadronic cascade in conjunction with the variation of the inelasticity y are most detrimental to the resolutions. The effect of limited photon sampling in the detector is of significantly less importance. These results will therefore also be applicable to similar detectors/media, such as those in ice.


Introduction
It has recently been suggested that the neutrino mass hierarchy can be determined with multi-megaton neutrino detectors in water and ice by measuring the flux of atmospheric neutrinos that have passed through the Earth [1]. Due to the influence of matter on neutrino oscillation during propagation through the Earth, the expected interaction rate of neutrinos in the energy regime of ∼ 1-20 GeV differs between the normal and inverted hierarchies. This is the main science goal pursued by KM3NeT/ORCA, the densely instrumented component of KM3NeT 2.0 being built off the southern coast of France [2]. In its current design, ORCA will consist of a 3D array of 2070 optical modules, each housing 31 3-inch photomultiplier tubes (PMTs), with about 20 m horizontal and 9 m vertical spacing between optical modules, instrumenting in total a volume of about 5.5 × 10 6 m 3 of seawater. The PINGU infill array of the IceCube detector at the South Pole has also been proposed to perform this measurement [3]. Additionally, such detectors are expected to further constrain neutrino-oscillation parameters, and have been proposed to study sterile neutrinos, non-standard interactions, dark matter, supernovae, and the internal structure of the Earth.
In order to resolve the mass hierarchy, detectors must determine the neutrino energy and arrival direction 1 -the latter defining the neutrino path through the Earth -and identify the type of interaction, as electron-neutrino charged current (CC), muon-neutrino CC, or neutral current (NC). 2 With ORCA, this will be done by studying the Cherenkov light signature using an array of PMTs. A first estimate of the accuracy of reconstruction methods targeting the primary outgoing lepton (e or µ) in CC interactions is given in the Letter of Intent for KM3NeT 2.0 [2].
Improving reconstruction methods will require using the hadronic component of CC interactions, in particular since it has been suggested to resolve the inelasticity y of the interaction to further add to the mass hierarchy sensitivity [4]. However, the instrumentation density of the ORCA detector will likely be insufficient to identify the Cherenkov cone of each individual particle (as with Super-Kamiokande [5]), or measure the paths of particles as a tracking detector. Effects causing fluctuations in the Cherenkov light signature will therefore limit the achievable reconstruction accuracy.
In this paper, the limits imposed by these intrinsic physical fluctuations on the energy and direction resolution of hadronic cascades are derived (as per section 2) by analysing the output of simulations for seawater in section 3. Limits on resolutions of muon tracks and electromagnetic cascades are much smaller and hence less significant, and are derived in appendices B and C. Combining these fluctuations produces limiting resolutions for the neutrino in ν e and ν µ CC interactions, as presented in section 4. 3 The effects of these limits on the reconstruction of the inelasticity in ORCA are investigated in section 5. Section 6 compares the derived limits to the reconstruction performance of ORCA using a full detector simulation, discusses their applicability to reconstruction in ice-based detectors to the neutrino mass hierarchy. The outcomes of this investigation are summarised in section 7.

Definitions
With the exception of scattering on atomic electrons, all neutrino interactions in the GeV range produce hadronic cascades, although the 'cascade' may only consist of a few particles emerging from the neutrino interaction. There are several possible definitions of the energy and momentum of such a hadronic cascade, p h = (E h , p h ), which may or may not include the target momentum, any sub-relativistic nuclear remnant, and may also weight outgoing particles by their expected Cherenkov light yield. The convention chosen here is to define p h as the momentum transfer q, i.e. through four-momentum conservation applied to the neutrino interaction vertex in the laboratory frame, as shown in figure 1. Thus, where p ν is the original neutrino four-momentum, and p the four-momentum of the outgoing lepton. Similarly, the inelasticity y (nominally, Bjorken y) of the interaction is defined as These definitions are chosen because they are most relevant for reconstructing the neutrino interaction properties. By not including the (unobservable) target four-momentum, variations in this contribution between interactions with the same p h will be included in the fluctuations due to different hadronic final states. In ORCA and like detectors, events are reconstructed by measuring the number of emitted photons, their arrival times and their spatial distribution. Unlike detectors using a dense array of PMTs on a 2D plane, such as Super-Kamiokande [5], the sparser 3D layout of ORCA makes it difficult to distinguish Cherenkov cones from individual particles. Therefore, it is assumed that signatures from the outgoing e, µ in ν e,µ CC interactions can be distinguished from the accompanying hadronic cascade as demonstrated in ref. [2], but that individual particles within the hadronic cascade cannot be identified.
Moreover, while the elongation of the cascade is important for a dense detector configuration, effects due to this elongation are not studied here, and only fluctuations in the emission direction and number of photons are considered. Furthermore, it is assumed that only unscattered photons retain directional information. Therefore, the total number of detected photons can be used to estimate the cascade energy E h , while the emitted directions of photons that are detected before being scattered can be used to estimate the cascade direction, Fluctuations in reconstructed energy and direction are considered according to the following scheme: • Hadronic-state fluctuations. Hadronic final states with identical p h can differ in the number and type of particles, and in their energies and directions. This will result in different photon signatures, and hence fluctuations in reconstructed quantities.
• Propagation fluctuations. The stochastic nature of energetic particle propagation implies that even identical initial particles will produce different photon signatures.
• All-photon limit. Even in the case that every emitted photon is detected, the above two effects will combine to give the all-photon limit to reconstruction accuracy.
• Photon-sampling fluctuations. Sampling only a fraction of the emitted photons will further limit the ability of reconstruction methods to determine event properties.
• Overall limit. The combination of hadronic-state, propagation, and photon-sampling effects results in the overall limit to reconstruction accuracy.
Values of the relative energy reconstruction error, ∆E/E, and the direction reconstruction error, ∆θ, will be derived for each of these cases.
The total and unscattered sampled photon fractions are characterised only by the wavelength-dependent probability of detecting each photon, given in figure 2 for the ORCA benchmark detector studied in ref. [2]. Note that this detector has a 6 m vertical spacing between optical modules, while a 9 m spacing was found to be optimal for the neutrino mass hierarchy determination.
As all sources of fluctuations are independent, the all-photon and overall limits are calculated by adding their component contributions in quadrature.

Simulations
The simulation methods used in this paper to calculate limiting resolutions differ significantly from reconstruction methods based on full detector simulations. As the goal is to derive limiting resolutions in a robust manner, the following simplifications are used when compared to full simulations of ORCA [2]. In all cases, these are intended to be optimistic, to ensure that the resulting limits are indeed limiting.
• PMT response ignored. The precise number and arrival times of photons are assumed to be measured. This ignores the specifics of PMT response and readout, timingand charge-calibration uncertainties, and imperfections in reconstruction methods.
[nm] λ KM3NeT Figure 2. Detection probability, P , for all unabsorbed photons (red) and for those that are additionally not scattered (blue), calculated from in-situ measurements of the inherent optical properties of deep-sea water [6], and the optical module effective area and density for the KM3NeT/ORCA benchmark detector design [2]. The peak near 440 nm arises from the absorption minimum.
• Generalised detector geometry. Detected photons are randomly sampled according to wavelength-dependent probabilities only. This ignores the specifics of the detection geometry, such as partially contained events, and that PMTs are located in 'clumps' (i.e. optical modules) at specific positions. This simplified geometry is equivalent to a detector composed of an infinite number of infinitely small detection elements, and thus ignores the non-uniform and correlated photon-detection induced by 'clumpiness'.
• Known photon origin. Photons from the hadronic cascades and the outgoing lepton in a CC interactions are assumed to be distinguishable, and background photons from atmospheric muons and natural light sources are ignored.
• No simulation uncertainties. It is assumed that neutrino interactions, particle propagation, and Cherenkov light production can be modelled perfectly, and reconstructions are based on perfect simulations.
Initial hadronic states are taken from simulations of neutrinos interacting in seawater using gSeaGen [7], which implements GENIE 2.8.4 [8]. Samples of 1000 such states for each of several discrete energies between 1 GeV and 20 GeV were selected in a sufficiently unbiased manner (see appendix A.3 and table 1 for further details). Variations about the mean properties of each sample characterise the hadronic-state fluctuations.
Each state was simulated 1000 times with an implementation of GEANT 3.21 [9] in the standard ANTARES simulation package [10]. Cherenkov emission was generated over the 300-710 nm range via the Frank-Tamm formula [11] using the wavelength-dependent phase velocity in seawater [12].
Variations between simulations of the same cascade are used to characterise fluctuations due to propagation. Implementations of both GHEISHA ('G-GHEISHA') [13] and FLUKA ('G-FLUKA') [14,15] were used for hadron tracking, with differences (discussed in appendix A.2) serving as a measure of systematic uncertainties. In most cases, this difference was found to be small. Neutrino resolutions presented in section 4 use G-GHEISHA results.
For each simulation of each hadronic final state, the properties of all emitted photons, and a sample randomly selected according to the detection probabilities of figure 2, are recorded, and used to estimate corresponding reconstruction errors as described below. Distributions of parameters related to reconstructed energy are approximately Gaussian, and fluctuations are characterised by their root-mean-square (RMS). Those related to direction reconstruction tend to have large tails due to individually significant scattering events, and are characterised using 68.27% quantiles in space angle.
In case of neutrino resolutions, the Monte Carlo truth values of p l and p h are taken from gSeaGen. The errors in reconstructing them are calculated by considering the leptonic and hadronic components separately, i.e. it is assumed that the photons from the lepton can be identified. Correlations between energy and direction errors for each component are accounted for. A brief investigation has shown correlations with other properties of neutrino interactions to be very small (appendix A.3), so they are ignored in this work.
A further approximation is required when E h or E e,µ do not equal one of the discrete simulated energies. In such cases, the expected properties of ∆E/E and ∆θ are interpolated for E > 1 GeV, and extrapolated for E < 1 GeV. Thus, neutrino resolution limits depending on components -especially hadronic cascades -with E < 1 GeV should be interpreted with care.

Energy resolution
An example of the fluctuations in the photon yield for a hadronic cascade with E h = 5 GeV is given in figure 3. The effect of different hadronic states is asymmetric, with relatively many events producing a large number of photons. A brief investigation suggests that these are hadronic cascades with a large electromagnetic component (γ, e ± , π 0 ). Particle propagation effects however are more symmetric and centrally-peaked. Their combined effect in the all-photon limit is compared to a Gaussian fit.
The hadronic cascade energy E h can be estimated from the total number of detected photons, with fluctuations therein producing fluctuations in the estimated energy. Since the total number of emitted photons is not linearly proportional to hadronic cascade energy (as in the case of electrons), the hadronic cascade energy is estimated using a fit to the total number of photons, N γ,h . This can be expressed as a fraction f h relative to that from electromagnetic cascades of the same energy, N γ,e (appendix C), i.e.
Hadronic state Propagation All-photon limit KM3NeT Figure 3. Probability density distributions of the relative deviations from the mean emitted photon number, ∆N γ /N γ , in hadronic cascades with E h = 5 GeV. Shown are the contributions from hadronic state (purple), particle propagation (green), and their combined effect in the allphoton limit (black). A Gaussian fit (black line) to the all-photon limit is shown to illustrate its appropriateness in characterising these fluctuations.
uses the same functional form as ref. [16], and is shown in figure 4 (black solid line). The energy dependence of f h means that a deviation in the number of photons will have a non-linear effect on the reconstructed energy. The relative error in reconstructed energy can be calculated via where the coefficient of ∆N γ /N γ is the energy error reduction factor and is also shown in figure 4. As f h is an increasing function, the resulting relative errors in reconstructed energy tend to be smaller than the intrinsic relative N γ fluctuations. For hadronic state and propagation fluctuations, variation in the Monte Carlo values of all N γ about the mean value are used for ∆N γ , while for photon sampling, the Poisson variation in the number of detected photons, N det γ , is used: The resulting relative energy errors, calculated by applying eq. 3.2 to relative variations in photon number from simulations (e.g. figure 3), are shown in figure 5. In the all-photon limit, fluctuations in hadronic state dominate over the whole energy range. The effect due to photon sampling is relatively small, so that the overall limiting resolution is only slightly larger than the all-photon limit.
The light yield of hadronic cascades (figure 4) and its intrinsic fluctuations (figure 5) were also investigated by ref. [16] down to 10 GeV, with broadly similar results.

Direction resolution
A robust method to estimate the mean hadronic cascade direction is to use the mean direction of detected unscattered photons -the validity of this method is discussed in appendix C.2. In the all-photon limit, this estimate will be made using all emitted photons, while for photon sampling, only the detected unscattered fraction is used. An illustration of the resulting effects on the direction resolution is given in figure 6 for a single hadronic final state. The distributions of direction errors for E h = 5 GeV are illustrated in figure 7. Limiting direction resolutions for hadronic cascades between 1 GeV and 20 GeV are presented in figure 8. Below ∼ 2 GeV, variation in the initial hadronic state dominates the direction resolution in the all-photon limit, while above this energy particle propagation effects dominate. Over the entire energy range, their combined effect in the all-photon limit is less significant than the effects of photon sampling, so that the overall limiting resolution is worse. The relative difference in results between the G-FLUKA and G-GHEISHA cases is less than ∼ 10%. Hadronic state Propagation All-photon limit Photon sampling Overall limit KM3NeT Figure 5. Relative errors (RMS) in reconstructed energy for hadronic cascades, due to hadronic state (purple), particle propagation (green), their combined effect in the all-photon limit (black), and the additional variation introduced due to photon sampling (blue) in the overall limit (red). The mean is calculated from the average between the G-GHEISHA (larger) and G-FLUKA (smaller) results, and the error bars cover the range between them.
The neutrino energy E ν and normalised direction, u ν , can be reconstructed using the interaction kinematics of eq. 2.1, i.e. by combining the reconstructed lepton properties with that of the hadronic cascade: Naively, the weights w and w h could be set to unity, at which point eq. 4.1 closely resembles eq. 2.1. However, there are several reasons for having w and w h differ from unity. Firstly, while the hadronic cascade momentum magnitude | p h | is assumed not to be reconstructable, it can be estimated using the reconstructed value of E h and the expected ratio which is plotted in figure 9. Fluctuations in this ratio can be considered as an additional hadronic-state contribution to neutrino direction resolution, which is not included in either the hadronic or lepton components.
A second factor influencing the weight is the relative accuracy of u and u h . In general, the best estimator of u ν will be found by reducing w h and increasing w to reduce the influence of the much larger fluctuations in u h , even if this biases the resulting estimate in the direction of the lepton. rarer, ν τ CC events can be approximated by combining electron, muon, and hadronic cascade resolutions according to τ decay products.
[deg]  Figure 6. Photon distributions from a hadronic cascade with E h = 10 GeV. The true cascade direction is in the u z direction, at (θ x , θ y ) = (0, 0) (green circle). The mean photon direction, averaged over 1000 simulation iterations with the same hadronic final state, is given by the purple star -the 3 • offset from (0, 0) reflects the intrinsic bias for this particular state. The distribution of all emitted photons for one particular propagation iteration is shown by the grey shading. The mean direction of all these photons is given by the black plus sign -the 8 • offset of this from (0, 0) measures the total intrinsic variation in the case of the all-photon limit, while the 5 • offset from the purple star measures the variation due to random particle propagation only. A random sample of unscattered photons, chosen according to the probabilities of figure 2, are shown by blue triangles, along with their mean direction (red diamond). The 25 • offset of the red diamond from the black cross measures the variation due to photon sampling only, while the offset of 18 • from (0, 0) gives the overall variation.
A third and final factor is the use of eq. 2.1 to constrain the measured values of p h and p , both by requiring they add to a physical value of p ν , and by imposing prior likelihoods on the implied interaction kinematics.
The general situation of choosing arbitrary weights between u and u h is illustrated in figure 10. Ultimately, the optimum method would define weights based on the reconstructed energies of both components, and the angle between them, to account for the three aforementioned factors. However, such an optimisation is beyond the scope of this paper.
In the following, three methods are investigated. The first method (section 4.1) is the naive 'four-momentum conservation' approach using w h = w = 1. This is motivated by the momentum-magnitude and statistical down-weighting factors mentioned above tending to cancel each other, and because it is simple and robust. The second method (section 4.2) [deg] θ ∆    Figure 9. Ratio E h /| p h | (calculated according to eq. 2.1) as a function of E h , for hadronic cascades extracted from both ν CC (red) and ν CC (blue) events. Error bars represent the variation characterised by the RMS. Note that the definition of p h ≡ q (eq. 2.1, and figure 1) requires the electron and hadronic cascades in ν e CC events cannot be distinguished. In all cases, the 'track length' method for energy reconstruction of muons is used (appendix B). ∆E/E ν e,µ CC: all-photon limit ν e CC: overall limit ν µ CC: overall limit Figure 11. Limitations on relative neutrino energy resolution (RMS) as a function of neutrino energy E ν and inelasticity y. Resolutions are shown as contour lines for ν e,µ CC in the all-photon limit (black), for ν e CC in the overall limit (red), and for ν µ CC in the overall limit (blue). The region E h < 1 GeV (grey shading) has been calculated using extrapolated hadronic cascade results, and should be interpreted with care.

Using naive four-momentum conservation
are very similar for ν e and ν µ CC events. The exception is the overall limits on energy resolution, which are shown separately. Figures 13 and 14 show the neutrino energy and direction resolutions for ν e , ν e , ν µ and ν µ CC events as a function of neutrino energy integrated over the corresponding y distributions. Both energy and direction resolutions are better for ν CC than for ν CC events, since the former have on-average smaller y. Because the muon track length is a more reliable energy estimator than the total light yield, the energy resolution for ν µ is slightly better than for ν e (see appendices B and C). For comparison, in case of measuring only the lepton energy -and ignoring the hadronic cascade energy -the resolution is nearly energy independent with ∆E/E ≈ 0.5 (0.3) for ν e,µ (ν e,µ ) CC events. Also, note that while the direction resolution for ν CC events is significantly better than the scattering angle φ ν, between neutrino and lepton, the direction resolution for ν CC events is slightly worse.

Using the expected neutrino-lepton scattering angle
Using the hadronic cascade for neutrino direction reconstruction is limited predominantly by the large error ∆θ h on the direction of the hadronic cascade. Two secondary factors are: the resolution of E h , which leads to errors in u ν through eq. 4.1; and the inability to reconstruct | p h |, which necessitates some assumption about E h /| p h |. An improvement therefore might be made by only using u h to define the plane of the interaction, and ν e,µ CC: all-photon limit ν e,µ CC: overall limit ν e,µ CC: scattering angle Figure 12. Limitations on neutrino direction resolution (68% quantiles) in ν e,µ CC events as a function of neutrino energy E ν and inelasticity y, in both the all-photon limit (black) and the overall limit (red). This is compared to the intrinsic scattering angle between the outgoing lepton and the neutrino (purple). The region E h < 1 GeV (grey shading) has been calculated using extrapolated hadronic cascade results, and should be interpreted with care. The 'kinks' in the bottom-left corner are due to a change in recoil mass (different interaction channels) and therefore different values of E h /| p h |.
using the expected median scattering angle med[φ ν, ] between the neutrino and lepton to define the direction in this plane. Such a method will tend to improve the neutrino direction resolution when the relative scatter in the angle φ ν, about med[φ ν, ] is small. Here, med[φ ν, ] is calculated from the reconstructed energies E reco and E reco h by setting Figure 15 shows the neutrino direction resolution limits for ν e CC events as a function of E ν and y, relative to those shown in figure 12. 5 For a large region of parameter space in the E ν -y plane, this method provides a better neutrino direction reconstruction (ratio < 1), with up to ∼ 25% improvement. It does not perform well at large y due to an intrinsically large scatter in φ ν, .

Treating electron neutrino charged-current events as a single cascade
If the detector is not able to distinguish between photons from the electromagnetic and the hadronic cascade in ν e CC events, the reconstruction must treat all photons as coming from a single cascade. This corresponds to applying the method for hadronic cascades of section 3 to entire ν e CC events, using the approximate 2.5:1 ratio of atmospheric ν e and ν e CC events in this energy range [18,19].   Figure 14. Direction resolution (68% quantiles) for ν µ (black solid), ν µ (black dashed), ν e (red solid) and ν e CC events (red dashed). For comparison, the 68% quantiles of the scattering angle φ ν, between the (anti)neutrino and lepton is shown as purple lines.
In figure 16, the resulting limits on the neutrino energy resolution for ν e and ν e CC events (ratio 2.5:1) are compared to those obtained when the photons from the outgoing  Figure 15. Ratio of neutrino direction resolution limits for applying the expected neutrino-lepton scattering angle method over the limits for the naive four-momentum conservation method (cf. figure 12) for ν e CC events in the overall limit as a function of E ν and inelasticity y. A ratio < 1 indicates the parameter space where the expected scattering angle method is more accurate.
e ± and hadronic cascade can be resolved. Due to the lack of knowledge about the source of the photons, at E ν = 10 GeV the energy resolution worsens from 16.8% to 19.3%. This deterioration emphasises the importance of the capability to identify the electron in ν e CC events in order to allow for more advanced reconstruction procedures.
In the unresolved photon source case, contributions from different sources of variation are separately shown in figure 16. The variation due to the unknown inelasticity of the interaction is calculated by fixing the number of photons from the electron and hadronic cascades to their mean values and varying y using the true distribution. The variation due to photon sampling is given by the average Poisson error on N γ (eq. 3.3), summed over both components. The remaining variation is attributed to the contribution from the hadronic state and propagation of the hadronic cascade, and is calculated by subtracting both previously mentioned contributions in quadrature from the overall energy resolution. In addition, when the photon source is unresolved, the energy resolution for a mixed composition of ν e and ν e CC events becomes worse than a simple linear combination of their individual results would suggest, since fluctuations of N γ about the joint mean are necessarily larger. The light yield of ν e and ν e CC events differs by ∼ 10% in the considered energy range.
Variations in inelasticity affect the e ± unresolved case by smearing the relation between the neutrino energy and the expected N γ . However, the effect is still non-zero in the e ± resolved case, where fluctuations in y imply that hadronic and electromagnetic energies must be estimated independently. This is why the resolutions for both cases are identical KM3NeT Figure 16. Relative neutrino energy resolution (RMS) for ν e and ν e CC events (ratio 2.5:1) in the overall limit, assuming that photons from the outgoing e ± and hadronic cascade can (red solid) and cannot (blue solid) be resolved. For the unresolved case, contributions from variation due to unknown inelasticity y (black dashed), hadronic cascade (purple dotted) and photon sampling (green dashed-dotted) are shown separately.
at low energies: when E ν is small, almost all light will come from the outgoing e ± ; E h will not be directly observable; and both cases will reduce to measuring E e only, and relating that to E ν based purely on the expected value of y.
In the case of direction resolution when the photon source is unresolved, mainly neutrino interactions with small inelasticity (y 0.25) have a worse resolution. Without the additional knowledge about the source of the photons, the neutrino direction resolution degrades from 11.6 • (9.2 • ) to 11.8 • (10.3 • ) for ν e (ν e ) CC at E ν = 10 GeV.

Resolution of interaction inelasticity
The different differential cross sections of ν and ν interactions allow them to be statistically separated, which will add significance to a neutrino mass hierarchy measurement [4]. The energy resolutions obtained in section 3 and appendices B and C can be readily converted into limiting resolutions on the interaction inelasticity y (eq. 2.2), using Figure 17 shows the y true distributions for ν µ and ν µ CC events, and the y reco distributions in the overall limit for E ν = 5 GeV. For ν e CC events, the distributions are very similar. Despite the relatively poor resolution on E h , the differences between the y reco and y true distributions are small. Partly, this is due to the lack of structure in the y distribution,  Figure 17. Probability density distributions, P , of inelasticity y true (solid) and y reco (dashed) for ν µ (red) and ν µ (blue) CC events with E ν = 5 GeV. The dip in the last bin (as y → 1) is caused by the finite muon mass. Due to a fraction of low-energy hadronic cascades producing very few photons, the y reco distribution is shifted to the left, causing the excess in the first bin (as y reco → 0).
particularly for ν e,µ interactions, and partly because the error in y reco tends to zero with both E h and E e , as can be seen from eq. 5.1. The statistical resolution between the y distributions from neutrinos P ν (y) and antineutrinos P ν (y) can be characterised using the correlation coefficient, c, by defining the separation power s as Perfect separation is indicated by s = 1, while no separation corresponds to s = 0. The degradation due to limiting accuracies on a detector's ability to distinguish between ν and ν is given by comparing s calculated using y true and y reco . Following the same procedure as in section 4.1, figure 18 shows the separation power as a function of neutrino energy both with and without reconstruction errors. Above E ν = 5 GeV, where the intrinsic separation power is highest, the relative reduction in separation power for both all-photon and overall limits is 5% or less. the confounding effects of the optical backgrounds, detector layout, PMT response, event triggering and selection, and uses events that are not necessarily fully contained in the detector.
Both ν µ and ν e CC event reconstruction methods described in ref. [2] aim to identify the direction of the outgoing lepton, with a limiting resolution for the neutrino direction of the scattering angle θ ν, (figure 14). Accounting for different presentation methods (space angle vs. zenith angle, median vs. 68% quantile), the resulting direction resolutions (figures 68, 81 and 82 in ref. [2]) are very close to this limitation. Thus, all prospects for improving the resolution using information from the accompanying hadronic cascade discussed in section 4 apply, and the resolution is not downgraded significantly in the case of a detector with a larger spacing between photon sensors. Furthermore, it indicates that the direction resolution of ORCA for µ and e is not significantly degraded by the aforementioned confounding effects. Indeed, the electron direction resolution of ORCA (figure 82 of ref. [2]) is close to the '1D' limit (figure 28 in appendix C) presented in this paper.
A comparison to the performance of the energy reconstruction for ν e CC events is shown in figure 19. Besides the resolution for the nominal event sample from ref. [2], the influence of poorly contained events is eliminated by applying a strict containment criterion based on Monte Carlo truth information. These resolutions are compared to both the e ± resolved and e ± unresolved case of ν e CC energy reconstruction limits described in section 4. The event selection using strict containment cuts shows an accuracy comparable to the limiting e ± unresolved accuracy -this is not unexpected, since event reconstruction in a real detector cannot perfectly differentiate between photons from the outgoing lepton and from the hadronic cascade. Relative energy resolution (RMS) for KM3NeT/ORCA, calculated using a full detector simulation on atmospheric ν e and ν e CC events. Results are shown using both the event sample from ref. [2] ('nominal LoI'), and one based on strict containment criteria on Monte Carlo truth neutrino interaction vertex and direction ('MC truth contained') [17]. This is compared to the overall limiting resolutions for ν e and ν e CC events (ratio 2.5:1) calculated here, both for the 'e ± resolved' and 'e ± unresolved' case, as discussed in section 4.3.
In the case of ν µ CC events, the length of the muon track makes it difficult to define a fully contained event for a significant range of the parameter space, and the results are currently too strongly affected by selection effects to make comparisons useful.

Applicability to in-ice detectors
The simulations in this paper have been performed in seawater, and assumptions on reconstruction methods are most applicable to the techniques proposed for ORCA. However, since ice and seawater are predominantly composed of H 2 O, both primary neutrino interactions and subsequent particle propagation are similar. The refractive index of ice is also similar to that of seawater at optical and near-UV wavelengths (between 1.3 and 1.4) [6,20]. Hence, limits due to different hadronic final states and particle propagation (i.e. all-photon limits) will therefore be almost directly applicable to in-ice experiments. Since ice has a lower absorption coefficient than seawater, the photon sampling effects on the energy reconstruction of electrons and hadronic cascades will be significantly lower than those presented here, though perhaps could be calculated with a simple scaling factor. Limits on the neutrino energy resolution (figure 11) however are dominated by fluctuations due to different hadronic states (cf. section 3.1) and not by limited photon statistics. For direction reconstruction, the photon sampling effect in ice is difficult to estimate: a much higher scattering coefficient reduces the number of unscattered photons, but as discussed above, scattered photons do contain some directional information.
Therefore, it is expected that both all-photon and overall limits on the energy reconstruction accuracy will apply equally well in ice, while direction reconstruction limits are strictly applicable only in the all-photon case.

Comparison to estimates in the literature
The limiting resolutions derived here can be used to assess the validity of several estimates of the sensitivity of PINGU, and hence approximately ORCA, to the neutrino mass hierarchy. Some, such as the calculations of refs. [21][22][23], use the published resolutions from the experiments, for which the comparisons above apply. Others, in particular refs. [1,24,25], use their own resolution estimates, or approximations of the published results. None of these authors take into account the correlations between energy and direction resolutions, or the differences between ν and ν interactions.
The neutrino direction resolutions from refs. [1,24,25] are compatible with the limits presented in this work. However, direction errors are often [1,25] assumed to be Gaussian in zenith angle and not in space angle (the limits presented in this work are functions of space angle).
The energy resolutions assumed by these authors [1,25] range from 15% to 35%, and tend to become incompatible with the limits (cf. figure 16) at low energies, i.e. they are too optimistic in this range.
Effects governing neutrino resolution have also been investigated in ref. [24]. However, the authors do not consider the fluctuations in hadronic cascades due to different hadronic final states or particle propagation. Additionally, they overestimate the light yield of hadronic cascades, assuming f h = 80% (eq. 3.1), independent of energy. This is significantly too large for the considered energies (cf. figure 4). Consequently, they underestimate the contribution due to variation in inelasticity y. Furthermore, the neglected effects contribute significantly to the limiting resolutions on neutrino energy. As a consequence, the energy resolutions assumed in ref. [24] are incompatible with the limits derived here.

Sensitivity to physical uncertainties
The logic of this paper is to assume that all physical processes are perfectly modelled, and fluctuations within these models are studied. In particular, this applies to the use of GENIE for neutrino interaction as well as hadronisation, and GEANT for particle propagation. Deviations in this modelling from the truth will have two effects. Firstly, it may lead to systematic shifts between reconstructed and true parameters. These can only act to worse the resolution, but might be accounted for through the use of nuisance parameters in global fits (see e.g. section 4.6 in ref. [2]). Secondly, intrinsic fluctuations about the mean behaviour could either increase or decrease, causing a corresponding worsening or improvement in resolutions. The sensitivity of ORCA to this second effect can be gauged from the relative effects of hadronic state, particle propagation, and photon sampling fluctuations.
As shown in section 3, the energy resolution on the hadronic cascade is dominated by the intrinsic fluctuations in the emitted Cherenkov light due to effects of hadronic state variations. Therefore, experiments such as ORCA might be sensitive to systematic uncertainties in the neutrino interaction and hadronisation process in the range of a few GeV. These processes are very difficult to model. GENIE [8] uses scaling methods [26] to extend PYTHIA [27,28] to this region. Nonetheless, there are several remaining discrepancies between experimental data and theory, and work is ongoing to improve the performance of both GENIE and other neutrino interaction simulation tools (see e.g. ref. [29], and references contained therein). The impact of variation in PYTHIA parameters has been found to be small for inelasticity measurements [29]. As discussed in section 5 however, inelasticity is also relatively insensitive to intrinsic fluctuations, and so is not a good indicator of model sensitivity. If experiments such as ORCA are affected by such systematic uncertainties, then data collected with ORCA might also help to constrain the models.
Determining this sensitivity however would require dedicated studies, which would be beyond the scope of this paper.

Conclusion
An investigation of the intrinsic limits on resolutions for neutrino events in the 1-20 GeV energy regime for the KM3NeT/ORCA water Cherenkov detector in the Mediterranean Sea has been performed. Limits on the energy and direction resolution of muon tracks, and of hadronic and electromagnetic cascades have been calculated. These have been combined to derive limits on primary ν µ and ν e CC resolutions, and the interaction inelasticity y.
The results indicate the best reconstruction accuracies achievable with the ORCA detector under the assumption that Cherenkov cones from individual hadronic particles will not be reconstructable. The results are presented for both the benchmark ORCA detector, and in the theoretical case that all photons are detected, which also limits resolutions for a detector with a significantly increased photocathode density.
It is found that the limits imposed by the methods described are close to the ORCA resolutions for ν e CC events obtained using full simulations, indicating that the influence of effects such as natural optical background light and the time-and charge-resolution of KM3NeT PMTs are at most small.
The main result is that the energy resolution of few-GeV neutrinos is primarily dictated by intrinsic fluctuations in the number of emitted Cherenkov photons. The uncertainty due to detecting only a small fraction of them plays only a minor role. Light yield fluctuations are dominated by the fluctuations in the hadronic final state in conjunction with the variation of the interaction inelasticity y, emphasising the importance of accurate neutrino interaction models. The neutrino direction resolution is dominated by the large errors in hadronic cascade direction reconstruction, which is mainly driven by detected photon statistics. Due to the significantly worse resolutions for hadronic cascades than for muons and electrons, the neutrino energy and direction errors, are strongly correlated via the inelasticity y. The methods allowing an optimum reconstruction which have been identified in this paper thus also depend on y. This emphasises the importance of resolving the outgoing lepton in a ν e,µ CC interaction -a capability which has already been demonstrated for ORCA [2]. Such a result is also important for the statistical discrimination between ν and ν using their different y distributions and, hence, increased mass hierarchy sensitivity. Importantly, it has been shown that intrinsic fluctuations do not significantly degrade this discrimination power.
The nature of these limits means that these conclusions will hold for neutrino detection in ice, and for different detector configurations, up to the point that a Super-Kamiokandestyle reconstruction of the Cherenkov cone from each individual particle becomes possible.
The results of this investigations also explain why the resolutions for ORCA presented in the Letter of Intent for KM3NeT 2.0 [2] do not significantly degrade for larger spacings between optical modules. This motivates an increase in the nominal vertical spacing between optical modules of the ORCA benchmark design from 6 m to 9 m. Moreover, the results allow the broader community to make improved estimates of the potential sensitivity of ORCA to the neutrino mass hierarchy, including for alternative detector layouts. Digitised results are available as supplementary material in the online version of this article.

A.1 GEANT settings
GEANT was set to track electrons down to their Cherenkov threshold of 0.25 MeV, and photons down to 0.4 MeV, at which the maximum energy of Compton knock-on electrons is 0.25 MeV. The simulation also allows the explicit production and tracking of δ-electrons above 0.25 MeV. To ensure the necessary accuracy of electron tracking at low energies, the standard GEANT 3.21 routine 'gtelec.f' had to be modified (setting IABAN to 0) to avoid the premature stopping of electrons. Additionally, slow neutrons with kinetic energies below ∼ 100 MeV were artificially removed from the simulation to avoid an unphysical number of neutron captures and subsequent nuclei decays. The GEANT implementations of both GHEISHA [13] and FLUKA [14,15] were used for hadronic tracking.

A.2 G-GHEISHA vs. G-FLUKA
Hadron tracking in GEANT 3.21 can be performed by either GHEISHA [13] or a preliminary version of FLUKA [14,15], often respectively termed 'G-GHEISHA' and 'G-FLUKA' to distinguish their implementations in GEANT from independent distributions. A comparison of these packages, both with each other, an independent version of FLUKA, and experimental data was performed in ref. [30] in the relevant range around 1-100 GeV. The most relevant observed deviations from expected behaviour were the non-conservation of energy and baryon number by G-GHEISHA in pion-nucleon interactions; and an underestimation of up  3, 3.5, 4, 4.5, 5, 6, 7, 8, 9, 10, 12, 14, 16, 18, 20  to 30% by G-FLUKA in the particle multiplicity from pion-carbon interactions below 7.5 GeV (G-GHEISHA accuracies were of order 10%).
In order to characterise the effects of the choice of G-GHEISHA or G-FLUKA, the resulting distributions of the mean number of detected photons and its standard deviation have been studied. For hadronic cascades with E h = 10 GeV, G-FLUKA shows on average a ∼ 6% higher light yield (cf. figure 4) and ∼ 10% less fluctuations than G-GHEISHA.

A.3 Selection of hadronic cascades
The chosen hadronic events were sampled from GENIE 2.8.4 simulations of ν µ CC events with original neutrino energies between 1 and 100 GeV. The first thousand events having a required hadronic energy E h (defined according to eq. 2.1) within 1% of the nominal values were chosen for each of the energies specified in table 1. Preliminary studies indicated that for a fixed E h , differences in event properties (mean light yield, particle content) due to a different choice of primary neutrino (flavour, energy) were negligible, in accordance with ref. [31].
In section 2.2, it is noted that more-detailed kinematical properties of neutrino events (four-momentum transfer squared q 2 for instance) may be correlated with energy and direction reconstruction errors, and that this association is lost when selecting hadronic final states using E h only. It was found that the correlation between direction error and q 2 is small, with correlation coefficient |c| < 0.05. The correlation between N γ (cascade energy, section 3.1) and q 2 is slightly larger, at c ≈ 0.15. This is still much smaller than the correlation between N γ and the direction error (c ≈ −0.5), which is accounted for in the results of section 4. Figure 20 compares the neutrino direction resolution obtained using the method of section 4 with that using the exact hadronic final state from each interaction, i.e. accounting for any/all correlations, including that with q 2 . A very good agreement is achieved.
Similarly, investigations have shown that the intrinsic energy and direction resolutions for hadronic cascades induced by Z 0 bosons are very similar to that induced by W ± bosons [17], with relative differences in all resolutions below 5% for E h > 3 GeV.
A further approximation is required when E h or E e,µ do not equal one of the simulated energies of table 1, in which case the expected properties ∆E/E, ∆θ must be interpolated between simulated values. When E < 1 GeV, the values of ∆θ and ∆E/E were set to their values at 1 GeV. In the case of energy resolution in the overall limit (including photon sampling), uncertainties due to photon sampling are extrapolated according to the photon statistics ( N γ ), and uncertainties due to muon track length determination are assumed to be independent of initial muon energy. The effects of these approximations are included in figure 20 for E ν = 3 GeV, y = 0.2 and y = 0.8. The agreement with full simulations is only slightly worse.

B.1 Energy resolution
In the considered energy range of 1-20 GeV, muons behave very similarly to minimum ionising particles: they experience a roughly constant energy loss per unit track length, and travel in approximately straight lines. Therefore, their signature in water and ice Cherenkov detectors is a straight track with a nearly uniform luminance.
This behaviour leads to two obvious methods of muon energy reconstruction: the track length, and the total Cherenkov light yield.
Light yield. The same principles as for hadronic cascade energy reconstruction discussed in section 3.1 are used, simplified by the number of emitted photons N γ being almost directly proportional to E µ . Thus, relative all-photon energy reconstruction errors are calculated as The sampling errors are given by eq. emission of muons in seawater. These two sources of variation (all-photon and photon sampling) are independent, and can be added in quadrature to find the overall uncertainty.
Track length. Assuming in all cases that the starting point of a muon track can be resolved perfectly using information from the accompanying hadronic cascade, the general strategy of a muon track length measurement is to detect the endpoint of a muon track using information from unscattered photons, i.e. by back-projecting onto the track assuming the Cherenkov angle. In the all-photon limit, it is assumed that the track endpoint can be determined precisely. Fluctuations in the muon track length L µ itself will then correspond directly to an error in the energy reconstruction. If ∆L µ (E µ ) is the fluctuation of L µ , the relative energy resolution in the all-photon limit is given by The distribution of the relative track length for 3 GeV and 10 GeV muons is shown in figure 21.
In the overall limit, the relatively small number of sampled unscattered photons will limit the resolution on the muon path. The muon track length measurement is then determined by the distance from the starting point to the position where the last detected, unscattered photon is emitted by the muon. The corresponding uncertainty, ∆L sampling µ , is estimated from the variation about the mean offset between the position of last emission and the true muon endpoint. For the assumed photon detection probability (figure 2), this results in roughly 0.1 GeV uncertainty on the muon energy. The overall limit on the relative energy resolution is calculated by adding these components in quadrature. 0.14 All-photon limit (track length) Overall limit (track length) All-photon limit (light yield) Overall limit (light yield) KM3NeT Figure 22. Limits on the relative muon energy resolution (RMS) as a function of muon energy E µ for the muon track length method (red circles) and the total light yield method (blue squares).
Resolutions are shown for the all-photon limit (hollow markers) and overall limit (filled markers).
Results. All-photon and overall limits on the relative muon energy resolution for both discussed energy estimates are shown in figure 22.
In the case of the all-photon limit, fluctuations in total photon emission are almost negligible, whereas those in the muon track length are not, so the total light yield method is more accurate.
In the case of the overall limit, the additional fluctuations due to photon sampling in the total light yield are significantly larger than the uncertainty due to detecting the track endpoint, so that even when combined with intrinsic track length fluctuations, the muon track length method is more accurate for this detector below 13 GeV. Therefore, the track length method is used when calculating results for ν µ CC throughout sections 4 and 5.

B.2 Direction resolution
In the case of the all-photon limit, under the assumption that the precise emission points of each photon can be reconstructed, in theory it would be possible to measure the initial direction of the muon track immediately after the neutrino interaction, giving perfect direction resolution. Here, a single constraint is introduced to obtain intrinsic limits: that the muon track fit assumes a linear muon path over its entire length, so that the linearity of the track itself limits the accuracy of the direction fit.
The effects of photon sampling in the overall limit to the resolution is unclear, since only a few unscattered photons (with perfect timing and position resolution) need to be detected in order to fully constrain a track fit. It might be possible to obtain such limits by considering the timing uncertainty of the PMTs, or the dispersion characteristics of seawater. However, such limits will be very small indeed, as verified by the very good  Figure 23. Distributions of the slopes dx/dz (dashed) and dy/dz (dotted) from fits to muons which began traveling along the z-axis with an initial energy of 3 GeV (red) and 10 GeV (blue). Gaussian fits (solid) are also shown.
direction resolution of muon tracks in ANTARES [32]. Hence, only the intrinsic fluctuations from the straightness of a muon track itself will be considered, and applied equally to both the all-photon and overall limits. The slopes dx/dz and dy/dz from fits to the true muon path (which always begins travelling along the z-axis) are used as the errors of such a muon track fit. Their distributions are fitted with Gaussians to derive their standard deviations, σ x,y . This is shown for 3 GeV and 10 GeV muons in figure 23.
The resulting direction error, ∆θ = √ 2σ x,y , is both the all-photon and overall limit, and is shown as a function of muon energy in figure 24. Note that these direction errors are dominated by multiple scattering, although the fitting procedure results in a different direction error than more common measures of multiple scattering angles [33].

C.1 Energy resolution
The energy of an electron can be best estimated by counting all detected photons produced by its subsequent cascade, which scales almost linearly with the electron energy E e . The resulting fluctuations can be calculated analogously to the total light yield method for muons discussed in appendix B.1. For the all-photon limit, distributions of the total number of detected photons N γ for 3 GeV and 10 GeV electrons are given in figure 25. For the overall limit, Poisson variation is calculated from the total number of detected photons N det e = 20.8 × E e /GeV. KM3NeT Figure 24. Limiting resolutions on muon direction reconstruction (68% quantiles) as a function of muon energy E µ .  All-photon limit Photon sampling Overall limit KM3NeT Figure 26. Limiting relative energy resolutions (RMS) as a function of electron energy E e , for the all-photon limit (black), the error due to photon sampling (blue, overplotted by red), and for the overall limit (red).

C.2 Direction resolution
The method of estimating errors in the reconstructed direction of hadronic cascades presented in section 3.2 ('2D' method) is motivated by the intrinsic asymmetry of the light pattern about the hadronic cascade axis. As this is not the case for a single electron, an alternative ('1D') method is presented below, which is deemed most appropriate for electron cascades. Electron cascades tend to produce an azimuthally symmetric distribution of Cherenkov light around the cascade axis, which can be described by a one-dimensional function of the emission angle θ with respect to the cascade axis. A direction reconstruction can then fit this function to trial cascade axes. Here, the direction error is characterised via fluctuations in the mean values of the θ distributions. Distributions of θ, and of the mean values, are shown in figure 27 for 5 GeV electrons.
Results for the full energy range are shown in figure 28. For comparison, the direction error resulting from the 2D method, as used for the hadronic cascade resolution in section 3.2, is also shown. Errors due to photon sampling, and hence the overall limits, are dominant over the entire energy range.
The 1D method underestimates fluctuations by using Monte Carlo truth information to define an axis of symmetry, while the 2D method introduces artificial fluctuations in the reconstruction of an azimuthally symmetric photon distribution by allowing a non-uniform distribution in the azimuthal angle to influence the reconstructed direction.
For both the 1D and 2D cases, the errors in the all-photon limit are nonetheless small compared to that due to photon sampling, so that the detected photon distribution will be approximately azimuthally symmetric. This is in contrast to the case for hadronic cascades All-photon limit Photon sampling Overall limit 5) × Mean distribution ( 5) × Single event (   Figure 27. Angular probability distributions of unscattered light from 5 GeV electrons, shown for the range of emission angles θ near the Cherenkov peak at θ C ≈ 42 • . The angular distribution averaged over 1000 cascades (purple circles) is compared to that from a single cascade (green squares); both are multiplied by 5 for clarity. The means of the angular distributions (calculated between 0 • and 180 • ) are in the approximate range 50 • -55 • , due to there being more photons above the Cherenkov angle than below. The distribution of these means from 1000 cascades is also shown, calculated using all photons (black), and a sample of photons to give the overall limit (red). The effect due to photon sampling only (blue, overplotted by red) is calculated relative to the mean θ value for all photons in the cascade.
(figure 8), where the large direction error in the all-photon limit is comparable to that due to photon sampling. For electron cascades therefore, azimuthal symmetry (i.e. the 1D method) could be used to reduce the error due to the small number of sampled unscattered photons (∼ 10/GeV), while for hadronic cascades in this energy range it could not be, and the 2D method must be used. The appropriateness of this choice of methods for both electron and hadronic cascades has been verified by full (and computationally expensive) likelihood reconstruction studies on subsets of data, producing results within O(10%) of those shown here. All-photon limit (1D) Photon sampling (1D) Overall limit (1D) All-photon limit (2D) Photon sampling (2D) Overall limit (2D) KM3NeT Figure 28. Limiting resolutions on electron direction reconstruction (68% quantiles), using both the 1D method described in appendix C.2, and the 2D method described in section 3.2, in the case of the all-photon (black) and overall (red) limits. The effect of photon sampling only is shown in blue (overplotted by red).