Active-sterile neutrino oscillations at INO-ICAL over a wide mass-squared range

We perform a detailed analysis for the prospects of detecting active-sterile oscillations involving a light sterile neutrino, over a large $\Delta m^2_{41}$ range of $10^{-5}$ eV$^2$ to $10^2$ eV$^2$, using 10 years of atmospheric neutrino data expected from the proposed 50 kt magnetized ICAL detector at the INO. This detector can observe the atmospheric $\nu_{\mu}$ and $\bar\nu_{\mu}$ separately over a wide range of energies and baselines, making it sensitive to the magnitude and sign of $\Delta m^2_{41}$ over a large range. If there is no light sterile neutrino, ICAL can place competitive upper limit on $|U_{\mu 4}|^2 \lesssim 0.02$ at 90\% C.L. for $\Delta m^2_{41}$ in the range $(0.5 - 5) \times 10^{-3}$ eV$^2$. For the same $|\Delta m^2_{41}|$ range, ICAL would be able to determine its sign, exploiting the Earth's matter effect in $\mu^{-}$ and $\mu^{+}$ events separately if there is indeed a light sterile neutrino in Nature. This would help identify the neutrino mass ordering in the four-neutrino mixing scenario.

Most of the data from the above mentioned experiments fit well into the standard three-flavor oscillation picture of neutrinos [30][31][32]. The three-neutrino (3ν) oscillation scheme is characterized by six fundamental parameters: (i) three leptonic mixing angles (θ 12 , θ 13 , θ 23 ), (ii) one Dirac CP phase (δ CP ), and (iii) two independent mass-squared differences 1 (∆m 2 21 and ∆m 2 32 ). However, not all the oscillation data are compatible with this three-flavor oscillation picture [33]. There are anomalous experiments (for recent reviews see [34,35]), which point towards oscillations with substantially large mass-squared difference (∆m 2 ∼ 1 eV 2 ) as compared to the well-known solar and atmospheric mass splittings. For example, the LSND experiment observed the appearance ofν e events in aν µ beam, which can be explained in a two-flavor oscillation framework with ∆m 2 ∼ 1 eV 2 [36,37]. Another indication came from the short-baseline (SBL) Gallium radioactive source experiments GALLEX [6] and SAGE [7], where the disappearance of ν e was observed. Several SBL reactor antineutrino experiments noticed a deficit in the observedν e event rates [38] in comparison with that expected from the calculation of reactor antineutrino fluxes [39,40]. This so-called reactor neutrino anomaly also strengthened the case in favor of neutrino oscillations driven by ∆m 2 ∼ 1 eV 2 . Recently, new model-independent hints in favor of SBLν e oscillations have emerged from the reactor experiments NEOS [41] and DANSS [42]. By performing the combined analysis of the spectral ratios of these two experiments, the authors of Ref. [43] have obtained an indication (∼ 3.7σ) in favor of SBL ν e oscillations 2 with ∆m 2 ≈ 1.3 eV 2 . A recent reanalysis of the ILL reactor data [54] also claims a ∼ 3σ evidence for a light sterile neutrino with ∆m 2 ∼ 1 eV 2 .
All these anomalies recorded at the SBL experiments cannot be explained with the help of three sub-eV massive neutrinos ν 1 , ν 2 , ν 3 and require the existence of a fourth mass eigenstate ν 4 at the eV scale. This cannot couple to W and Z bosons due to the LEP bounds [55] on the number of the weakly interacting light neutrino state. Hence it has to be a gauge singlet, and can reveal its presence only through its mixing with the active neutrinos. This non-interacting neutrino is popularly known as the "sterile" neutrino. New more sensitive SBL experiments have been planned (see the reviews in Refs. [56][57][58][59]) to test the existence of sterile neutrino by observing the typical L/E pattern in the oscillations driven by the new mass-squared splitting ∆m 2 41 ∼ 1 eV 2 . Apart from the SBL oscillation experiments, the light sterile neutrino can also have visible effects in ongoing [60][61][62] or upcoming [63][64][65][66][67] long-baseline experiments. Signals of a sterile neutrino may be observed at the existing or planned multi-purpose water-or ice-based large detectors like Super-Kamiokande [49], IceCube [50,68,69], or DeepCore [70]. The β-decay and neutrinoless double β-decay experiments can also feel the presence of a light sterile neutrino (see the review [57] and the references therein). The possible existence of a light sterile neutrino can also have profound implications in cosmology [56,71,72].
The sterile neutrino is an elementary particle beyond the Standard Model (SM) and the possible discovery of a light sterile neutrino would prove that there is new physics beyond 1 We define ∆m 2 ij ≡ m 2 i − m 2 j , where mi's are the masses of the neutrino mass eigenstates νi's, arranged in the decreasing order of electron flavour content. To address the solar neutrino anomaly, we need ∆m 2 21 ≈ 7.5 × 10 −5 eV 2 and to resolve the atmospheric neutrino anomaly, we require |∆m 2 32 | ≈ 2.5 × 10 −3 eV 2 . 2 However, there are SBL experiments like CCFR [44], CDHSW [45], and SciBooNE-MiniBooNE [46,47] that have observed null results while searching for muon flavor disappearance associated with this high oscillation frequency. The same is also true for the long-baseline experiment MINOS [48] and the atmospheric neutrino experiments Super-Kamiokande [49] and IceCube [50]. Thus, there are tensions among various data sets, while trying to fit all of them in a four-neutrino framework [34,35,[51][52][53]. Therefore, the existence of a light sterile neutrino is still inconclusive. the SM at low-energies, which is completely orthogonal to new physics searches at highenergies at the Large Hadron Collider (LHC). There are several interesting motivations to search for a light sterile neutrino at mass scales different than the eV scale. For instance, the possible existence of a super-light sterile neutrino [73][74][75][76], which weakly mixes with the active neutrinos and has a mass very close to the active ones (∆m 2 41 ≈ 10 −5 eV 2 ), may explain the suppression of the upturn in the energy spectrum of solar neutrino events below 8 MeV [11,12,14]. A very light sterile neutrino at a mass scale smaller than 0.1 eV could affect the oscillations of reactor antineutrinos [77][78][79][80]. A light sterile neutrino having a mass of a few eV can help in nucleosynthesis of heavy elements inside the Supernova [81][82][83]. A relatively heavy sterile neutrino with mass around keV scale can act as warm dark matter in the context of the Neutrino Minimal Standard Model (νMSM) [84][85][86][87]. Therefore, it makes perfect sense to investigate the presence of a light sterile neutrino over a wide range of mass scale. Recently, the Daya Bay Collaboration also has performed such a search and has obtained bounds on the active-sterile mixing parameters in the mass range of 2 × 10 −4 eV 2 |∆m 2 41 | 0.3 eV 2 [80]. In this paper, we perform for the first time a detailed analysis of the search for a light sterile neutrino over a wide ∆m 2 41 range, starting from 10 −5 eV 2 and going all the way up to 10 2 eV 2 , using 10 years of atmospheric neutrino data expected from the proposed 50 kt magnetized Iron Calorimeter (ICAL) detector [88] under the India-based Neutrino Observatory (INO) project [89]. The magnetized INO-ICAL detector would detect the atmospheric ν µ andν µ separately over a wide range of energies and baselines. These features would help the ICAL detector to probe the presence of ∆m 2 41 over a wide range and also to determine its sign. The charged-current (CC) interactions of ν µ andν µ inside the ICAL detector will produce µ − and µ + particles, respectively. The ICAL detector would be able to measure the energies and directions of these muons to a very good precision [90]. Moreover, the ICAL detector would be sensitive to neutrinos in the multi-GeV range, and the energies of hadron showers produced by the interactions of these multi-GeV neutrinos can also be well-estimated [91,92]. While the main aim of this detector is to resolve the issue of neutrino mass ordering [88,[93][94][95] and to improve the precision on atmospheric neutrino mixing parameters [88,94,[96][97][98], the above features can also be instrumental in the search for a light sterile neutrino.
A sensitivity study of ICAL to a sterile neutrino in the range 0.1 eV 2 ∆m 2 41 10 eV 2 has recently been performed [99], where exclusion limits on the mixing parameters with 500 kt-yr exposure have been shown. In this paper, we not only extend this analysis to much lower values of ∆m 2 41 (all the way down to 10 −5 eV 2 ), but also point out for the first time that in the ∆m 2 41 range of (0.5 − 5) × 10 −3 eV 2 , ICAL has a far greater sensitivity to the magnitude as well as sign of ∆m 2 41 . We also explore the impact of hadron energy calibration, the actual value of |U e4 |, up-going vs. down-going events, spectral information, and the signs of ∆m 2 31 as well as ∆m 2 41 , on the sensitivity of ICAL to a sterile neutrino. In addition, in the scenario where a light sterile neutrino exists in Nature, we study the role of ICAL, in particular its muon charge identification capability, in identifying the mass ordering configuration of the four-neutrino spectrum.
We plan this article in the following fashion. We start Sec. 2 with a discussion on  different mass ordering schemes and configurations which are possible in the four-neutrino (4ν) framework. Then we describe the lepton mixing matrix in the 4ν scheme and mention the present constraints that we have on |U e4 | 2 , |U µ4 | 2 , and |U τ 4 | 2 from various oscillation experiments. We further elaborate on the Earth matter effects in the presence of three active neutrinos and one light sterile neutrino. We end this section by drawing the neutrino oscillograms in E ν -cos θ ν plane in the 4ν paradigm. In Sec. 3, we mention our analysis procedure and present the expected event spectra at the ICAL detector as a function of E µ and cos θ µ for several benchmark choices of active-sterile oscillation parameters. We identify the regions in E µ and cos θ µ plane which give significant contributions in constraining active-sterile oscillations. In Sec. 4, we determine the parameter space in the |∆m 2 41 |-|U µ4 | 2 plane that can be excluded by the ICAL data if there is no light sterile neutrino in Nature. We also study the dependence of these constraints on |U e4 |, mass ordering schemes, as well as energy and direction measurements. In Sec. 5, we assume that there is a light sterile neutrino in Nature and then address several interesting issues, like the precision in the determination of ∆m 2 41 and the chances of measuring its sign at the INO-ICAL atmospheric neutrino experiment. Finally in Sec. 6, we summarize and draw our conclusions.

Mass ordering schemes and configurations
The three active neutrinos ν e,µ,τ and the sterile neutrino ν s give rise to four mass eigenstates, ν 1,2,3,4 . Here ν 1,2,3 are dominated by active flavors, in the decreasing order of electron flavor fraction, while ν 4 is dominated by the sterile neutrino. The mass eigenstates may be arranged, according to their relative mass squared values, in eight possible configurations that are allowed by the current data. We label these configurations initially according to the signs of ∆m 2 31 and ∆m 2 41 (N-N, N-I, I-N and I-I, as shown in the top row of Table 1), and further by the relative position of ν 4 on the m 2 i -scale, which will depend on the signs of ∆m 2 42 and ∆m 2 43 , as indicated in the lower rows of Table 1. The configurations have also been shown graphically in Fig. 1.
In this paper, we shall focus on exploring these configurations by observing atmospheric neutrinos over a wide range of energies and baselines at INO-ICAL. Some of these   Table 1. configurations, viz. N-N-1 and I-N-1, correspond to the sterile neutrino solution required to satisfy the SBL anomalies as discussed in the previous section, while still being consistent with the cosmological bounds on the sum of neutrino masses [100]. Some others, viz. N-N-3 and I-N-2, may be relevant for a super-light sterile neutrino that has been proposed [73][74][75][76] to explain the non-observation of the upturn in the low energy solar data [11,12,14]. However, in this work, we shall not be looking for a solution to any of these anomalies, rather we shall perform an agnostic search for a sterile neutrino that may belong to any of these configurations. Note that the configurations N-I and I-I-2 will be disfavoured by the cosmological bounds [100], unless the magnitude of ∆m 2 41 is very small.
The Super-Kamiokande Collaboration limits |U µ4 | 2 to less than 0.04 for ∆m 2 41 > 0.1 eV 2 at 90% C.L. in N-N-1 configuration using 4438 live-days of atmospheric neutrino data and assuming |U e4 | 2 = 0 [49]. The IceCube neutrino telescope at the South Pole searches for active-sterile oscillation in their measured atmospheric muon neutrino spectrum as a function of zenith angle and energy in the 320 GeV -20 TeV range [50]. The muon antineutrinos in their data sample do not experience a strong MSW resonance effect, which allows them to obtain a bound on the active-sterile mixing. As an example, for ∆m 2 41 ≈ 1.75 eV 2 (N-N-1 configuration), which is the best-fit value of ∆m 2 41 according to Ref. [34], the IceCube data give a bound of sin 2 2θ 24 0.06 at the 90% C.L. 3 . In the same ∆m 2 41 range, the MINOS Collaboration places a new upper limit of sin 2 θ 24 < 0.03 at 90% C.L. using their long-baseline data, suggesting that |U µ4 | 2 0.03 [48,102].
The combined analysis of the available data from the reactor experiments Daya Bay and Bugey-3 provides a new upper limit of sin 2 2θ 14 0.06 at 90% C.L. 4 around ∆m 2 41 ≈ 1.75 eV 2 [102]. Recently, the reactor antineutrino experiments NEOS [41] and DANSS [42] have provided new hints in favor of short-baselineν e oscillations, which are model independent in the sense that these indications are not dependent on the precise estimate of the reactor ν e fluxes. In N-N-1 configuration, the combined analysis of the DANSS and NEOS spectral ratios predict a narrow-∆m 2 41 island at ∆m 2 41 ≈ 1.75 eV 2 with sin 2 2θ ee = 0.049 ± 0.023 at 2σ [43], which means that the best-fit value of |U e4 | 2 is around 0.012. Another recent study on the same topic can be found in Ref. [52].
As far as the constraint on |U τ 4 | 2 is concerned, the Super-Kamiokande experiment places an upper limit of 0.18 on this parameter at 90% C.L., for ∆m 2 41 > 1 eV 2 [49]. The IceCube DeepCore Collaboration sets upper bounds of |U µ4 | 2 < 0.11 and |U τ 4 | 2 < 0.15 at 90% C.L. for ∆m 2 41 ∼ 1 eV 2 , using their three years of atmospheric neutrino data in the range 10 to 60 GeV [103]. One can test the mixing between sterile and tau neutrino by looking for a depletion in the neutral-current event rates at the far detector of a longbaseline setup. Both the MINOS [104] and NOνA [105] experiments place competitive constraints on mixing between sterile and tau neutrinos using this approach. According to the global fit study performed in Ref. [35], |U τ 4 | 2 0.014 at 90% C.L. for any value of ∆m 2 41 in N-N-1 configuration. The authors of Ref. [34] obtain similar bounds on |U τ 4 | 2 around ∆m 2 41 ≈ 6 eV 2 . The INO-ICAL experiment will probe a combination of the ν e → ν µ and ν µ → ν µ oscillation channels, mainly via the charged-current (CC) interactions that produce muons. The active-sterile mixing elements of relevance here are therefore U e4 and U µ4 (and to some extent U τ 4 , due to Earth matter effects). In order to illustrate the impact of activesterile oscillations for different choices of ∆m 2 41 , we shall be using the benchmark values |U e4 | 2 = 0.025, |U µ4 | 2 = 0.05, and |U τ 4 | 2 = 0 while showing our oscillograms and event plots. For most of our sensitivity plots, we shall take |U e4 | 2 = 0, unless otherwise mentioned. We shall also demonstrate later that |U e4 | 2 = 0 yields the most conservative bounds in the (∆m 2 41 , |U µ4 | 2 ) parameter space. As far as the three-neutrino oscillation parameters are concerned, Table 2 shows the benchmark values that we consider in this work.
leptonic mixing matrix elements as [106,107]: For antineutrinos, the oscillation probability follows Eq. 2.2 with the replacement of U with its complex-conjugate matrix. While propagating through the Earth, the elements of the leptonic mixing matrix change, since the forward scattering of neutrinos on the Earth matter gives rise to effective matter potentials between ν e and ν s , in the flavor basis. Here N e is the electron number density and N n is the neutron number density inside the Earth. These modified elements need to be used while computing the neutrino oscillation probabilities through layers of different densities inside the Earth. As neutrinos start crossing longer distances through the Earth, the Earth matter effects start becoming important, as the relative matter potentials V es and V µs are ∼ ∆m 2 31 /(2E) for atmospheric neutrinos (e.g. for E = 5 GeV and L = 5000 km). As a consequence, for ∆m 2 41 ∼ ∆m 2 31 , the Earth matter affects active-sterile conversion probabilities significantly. In this ∆m 2 41 region, the similar magnitudes ∆m 2 41 and ∆m 2 31 also lead to an interference between the oscillation frequencies governed by them, which will manifest itself in the sensitivity of the ICAL experiment in the ∆m 2 41 -|U µ4 | 2 plane, as we shall see later in our results. For ∆m 2 41 ∆m 2 31 , which is the parameter space usually focused on (since it is relevant for explaining the LSND anomaly), matter effects are not very important. However for ∆m 2 41 ∆m 2 31 , a region of the parameter space that we specially explore, the Earth matter effects would play a major role.
In order to compute the neutrino conversion probabilities in the presence of Earth matter numerically, in both the three-flavor and four-flavor scenarios, we use the GLoBES software [108,109] along with its new physics tools. To take into account the Earth matter effects, we take the Preliminary Reference Earth Model (PREM) profile for the density of the Earth [110], with five density steps (to keep the computation time manageable). For the down-going neutrinos (L 450 km), we take all the neutrinos to be produced at an uniform height of 15 km, and consider them to be travelling through vacuum. For the upward-going neutrinos (L 450 km), we neglect the height of the atmosphere, which is small compared to the total distance travelled. We have checked that these approximations do not affect the oscillation probabilities, and hence the event distribution in energy and zenith angle, to any appreciable extent.

Oscillograms in E ν -cos θ ν plane
In order to get an idea of the impact of active-sterile mixing on the ν µ → ν µ survival probability P µµ and the ν e → ν µ conversion probability P eµ , we define the quantities where "(4f)" and "(3f)" denote the quantities calculated in the 4-flavour and 3-flavour mixing scenario, respectively. In Fig. 2, we plot these quantities in the E ν -cos θ ν plane, where θ ν is the zenith angle and E ν is the energy of the neutrino. We present our results for two different values of ∆m 2 41 , viz. ∆m 2 41 = 1 eV 2 (corresponding to the N-N-1 configuration), and ∆m 2 41 = 10 −3 eV 2 (corresponding to the N-N-2 configuration). In the former scenario, we average over the fast oscillations when neutrinos travel large distances through the Earth. From the figure, the following observations may be made: • Let us first see the impact of the active-sterile mixing in the ν µ → ν µ survival channel, which plays an important role in atmospheric neutrino oscillations. The effect of sterile neutrino mixing on P µµ is clearly significant. For ∆m 2 41 = 1 eV 2 , the effect is observed over a wide region in E ν and cos θ ν . For ∆m 2 41 = 10 −3 eV 2 , the impact is mostly confined to the upward going neutrinos, while its magnitude is significant around E ν ∈ [3, 10] GeV and cos θ ν ∈ [−1, −0.7].
• While fast oscillations for the high ∆m 2 41 ∼ 1 eV 2 will get averaged out due to the finite energy and angular resolutions of the detector, at ∆m 2 41 ∼ 10 −3 eV 2 the oscillation probability dependence on energy and direction of neutrinos would be more clearly resolvable. Hence we expect that the spectral and angular information should play a role in identifying signatures of sterile neutrinos at low ∆m 2 41 values .
• The impact of 4ν mixing on P eµ is confined to a narrow region to the upwardgoing neutrinos, and it is more prominent at the lower ∆m 2 41 . This would indicate a significant contribution of the matter effects due to the Earth, which influence the upward-going neutrinos.
The much higher effects of active-sterile mixing on P µµ , combined with the much higher value of P µµ as compared to P eµ , and the higher fraction of ν µ in the atmospheric neutrino flux, indicates that the ν µ → ν µ conversions will play a dominant role in the sensitivity of atmospheric neutrino detectors to sterile neutrinos. The corresponding antineutrino oscillation probability plots show the same features. In fact, the oscillograms at the higher ∆m 2 41 values are virtually identical with those for neutrinos, shown in the top panels of Fig. 2. Therefore 5 , we should be able to interpret most of our results based on P µµ .
The oscillograms of P µµ at low ∆m 2 41 indicate the presence of significant Earth matter effects. This points to the possible sensitivity of the data to the sign of ∆m 2 41 , a point we shall explore further in Sec. 5.2.

Event distributions in energy and zenith angle
In order to perform the simulations of the events in INO-ICAL, we generate the events using the same procedure as described in [88,93,94]. The NUANCE event generator [111] is used to generate unoscillated atmospheric neutrino events for 1000 years exposure of the 50 kt ICAL, which are later rescaled to the relevant exposure. The oscillations are implemented with the four-flavor probabilities in matter, using the reweighting algorithm, as described in [93]. The energy and direction resolutions of muons and hadrons are used, as per the latest results of the ICAL Collaboration [90,91]. We estimate the physics reach of ICAL in probing active-sterile oscillation parameters by combining the muon momentum information (E µ , cos θ µ ) and the hadron energy information (E had ≡ E ν − E µ ) on an event-by-event basis. For this work, it is assumed that the muon and hadron hits can be separated with 100% efficiency, and that the background and noise are negligible. The 3f  4736  2070  3108  1361  1628  709  4f (1 Table 3: Number of events for 500 kt-yr exposure of ICAL, for Eµ ∈ [1,11] GeV and NH in the active sector. We have taken |Ue4| 2 = 0.025, |Uµ4| 2 = 0.05, and |Uτ4| 2 = 0. The information on hadron energy is not used. For sterile neutrinos, the results have been shown for three different ∆m 2 41 values, corresponding to the N-N-1, N-N-2 and N-N-3 configurations, respectively. events are binned into 10 uniform E µ bins in the range [1,11] GeV, and 20 uniform cos θ µ bins in the range [−1, 1]. The hadrons are binned in 4 E had bins: 1-2 GeV, 2-4 GeV, 4-7 GeV and 7-11 GeV, as in [94]. Table 3 gives the number of µ − and µ + events, in 3ν case and with three benchmark values of ∆m 2 41 for 4ν configurations. The hierarchy in the active sector has been taken to be normal (∆m 2 31 > 0). It may be observed that for ∆m 2 41 = 1 eV 2 (N-N-1 configuration), the number of events with |U µ4 | 2 = 0.05 is about 10% less than the number with vanishing |U µ4 |, as expected after averaging over fast ∆ 41 oscillations. For ∆m 2 41 = 10 −3 eV 2 (N-N-2 configuration), the fractional difference in the total number of events with finite vs. vanishing sterile mixing is rather small. However Up/Down ratios of the µ − and µ + events are different in these two scenarios. Thus for low ∆m 2 41 , the information on angular distribution of these events may be useful.
The last row in Table 3 shows the scenario where ∆m 2 41 is extremely small, ∆m 2 41 = 10 −5 eV 2 (N-N-3 configuration). Note that even at this extremely low ∆m 2 41 value, the effect of oscillations is still visible through the loss in number of events as well as the difference in the ratio of Up/Down events. This may be attributed to the fact that even when ∆m 2 41 → 0, active-sterile oscillations may still occur through the frequencies corresponding to |∆m 2 42 | ∼ 10 −4 eV 2 and |∆m 2 43 | ∼ 2.4 × 10 −3 eV 2 . As a result, the ∆m 2 41 → 0 limit in the 4ν mixing case does not correspond to the decoupling of sterile neutrinos, and limits on active-sterile mixing may still be obtained in this case. In our analysis, we go down to values of ∆m 2 41 as low as 10 −5 eV 2 . Figure 3 shows the distribution of µ − events as a function of E µ , the reconstructed muon energy in two sample cos θ µ -bins, where the efficiencies and resolutions of the detector are folded in. Clearly, the number of events should be smaller in the 4ν mixing scenario as compared to the 3ν mixing scenario, due to the nonzero |U µ4 |. The fractional difference is almost uniform over energy for ∆m 2 41 = 1 eV 2 (N-N-1 configuration). The absolute difference in the number of events is the largest in the lowest energy bin and decreases at larger energies. This indicates that the low-energy data will contribute more significantly to the identification of sterile mixing. For ∆m 2 41 = 10 −3 eV 2 (N-N-2 configuration), on the other hand, the difference in the number of events as a function of energy is small and the energy dependence is mild. We have explicitly checked that the features for ∆m 2 41 = 10 −5     N-N-1 configuration), the fractional loss in the number of events is observed to be almost uniform throughout the cos θ µ range. On the other hand, ∆m 2 41 = 10 −3 eV 2 (N-N-2 configuration) yields almost no loss of downward-going events (cos θ µ 0) as compared to the three-flavour scenario, while for the upward-going events (cos θ µ 0), a clear loss of events is observed. This is due to the fact that at such low ∆m 2 41 values, the downward-going neutrinos do not have enough time to oscillate to the sterile flavor. We have explicitly checked that the features for ∆m 2 41 = 10 −5 eV 2 (N-N-3 configuration) are the same as those for N-N-2 configuration. This non-trivial dependence of the event spectrum on the direction of the muon is instrumental in providing sensitivity to the value of |U µ4 | at low ∆m 2 41 .
3.2 Useful regions in the E µ -cos θ µ plane We constrain the active-sterile neutrino mixing by using where the χ 2 are defined in a manner identical to that in [93] while performing the analysis using only the muon energy and direction (the so-called "2D" analysis). While using the additional information on hadron energy (the so-called "3D" analysis), the χ 2 definition as in [94] is used. Following [93], the systematic errors on flux normalization (20%), zenith angle distribution (5%), tilt (5%), and cross sections (10%), as well as an overall systematics (5%) have also been included. In our analysis, the values of the six oscillation parameters of active neutrinos are taken to be fixed, in the simulated data as well as in the fit (see Table 2). In order to judge the impact of θ 23 and ∆m 2 31 , the two parameters that are expected to influence the atmospheric neutrino measurements strongly, we have also performed the analysis by marginalizing over these two parameters in their current 3σ allowed ranges. The other three parameters that appear in the oscillation probabilities, viz. ∆m 2 21 , θ 12 , and θ 13 , are known to a very good precision, and anyway appear in the atmospheric neutrino oscillation probabilities as sub-leading terms [112]. Also, the impact of δ CP in the ICAL experiment is very mild [88]. We observe that, due to the marginalization over θ 23 and ∆m 2 31 , the sensitivity of ICAL to sterile neutrino decreases marginally at low ∆m 2 41 , though it remains unchanged at ∆m 2 41 1 eV 2 . Even at lower ∆m 2 41 values, with the present relative 1σ precision of ∼ 5% on sin 2 2θ 13 [30][31][32], the impact of these marginalizations is very small and the results presented in this paper stay valid.
Before we present our final sensitivity results on active sterile oscillation parameters, we try to locate the regions in (E µ -cos θ µ ) plane which contribute significantly towards ∆χ 2 . In Fig. 5, we show the distribution of ∆ χ 2 , i.e. ∆χ 2 per energy-cos θ interval (as in [94]), in the reconstructed (E µ -cos θ µ ) plane. Note that for this particular figure, we have not included the effects of the five systematic uncertainties mentioned above. The following insights are obtained from the figure.
• For ∆m 2 41 = 1 eV 2 (N-N-1 configuration), we expect there to be a overall suppression, with significant contributions coming from low energy downward going events. This is because at such high values of ∆m 2 41 , active-sterile oscillations develop even from downward neutrinos that travel ∼ 10 km in the atmosphere, while the active ν oscillations would not have developed. This feature would be prominent in the µ − events, but not so much in the µ + events due to lack of statistics (See Table 3). Overall, we expect there to be an averaged overall suppression of events. • For ∆m 2 41 = 10 −3 eV 2 (N-N-2 configuration), the information in ∆χ 2 is concentrated in the low-energy upward-going events. This is because at low ∆m 2 41 , it takes ∼ 1000 km propagation for the oscillations to develop. Also it may be noticed that the ∆ χ 2 in this region is more significant for µ + events. This follows from the effects of the MSW resonance of active-sterile mixing in matter, which appears in the antineutrino channel for the N-N-2 configuration.
• The events around cos θ µ ∼ 0 do not have much information, since the efficiency of the detector is less for horizontal events [90].
In the next section, we discuss the constraints on the active-sterile mixing parameters in the ∆m 2 41 -|U µ4 | 2 space, and gain insights into the dependence of these constraints on other mixing parameters like |U e4 | and the mass ordering configuration.

Constraining active-sterile mixing
In this section, we present the results on the reach of ICAL for excluding the mixing parameter |U µ4 | 2 , as a function of ∆m 2 41 . We generate the data in the absence of sterile neutrino, and try to fit it with the hypothesis of the presence of a sterile neutrino, with ∆m 2 41 ranging from 10 2 eV 2 all the way down to 10 −5 eV 2 . This range covers, and goes well beyond, the ∆m 2 41 range relevant for the LSND anomaly on the higher side and the solar ∆m 2 41 on the lower side. We use an exposure of 500 kt-yr, and present the results in terms of 90% C.L. exclusion contours in the ∆m 2 41 -|U µ4 | 2 plane. We take the values of active neutrino mixing parameters from Table 2, and take U τ 4 to be zero throughout the analysis. We have checked that, when the accuracy in the measurement of active mixing parameters expected after 10 years is taken into account (using priors), the variation in the values of these parameters does not affect the results to any significant extent.
In Sec. 4.1, to begin with we restrict ourselves to |U e4 | = 0, and the N-N ordering scheme (∆m 2 31 > 0 and ∆m 2 41 > 0) for the sake of clarity. In Sec. 4.1.1, we study the dependence of the exclusion contours on nonzero |U e4 | and find that the bounds in the (∆m 2 41 ,|U µ4 | 2 ) parameter space are the most conservative with |U e4 | = 0. We also check the dependence of the results on the mass ordering schemes by computing the exclusion contours for the remaining ordering schemes (N-I, I-N, and I-I), and find it to be very mild. We therefore focus on the N-N ordering scheme and |U e4 | =0 in the rest of Sec. 4. In Sec. 4.1.2 we study the relative importance of the total rate of events and the shape of their energy as well as angular distributions, and in Sec. 4.2 we compare the INO-ICAL sensitivity in the ∆m 2 41 -|U µ4 | 2 plane with that of the other ongoing experiments. 2 41 -|U µ4 | 2 plane Figure 6 shows the exclusion contours obtained using only the muon momentum information (2D analysis), and those obtained by adding the hadron energy information (3D analysis), in the N-N mass ordering scheme. The following features of the reach of ICAL, expressed in terms of the exclusion contours, may be observed and interpreted in terms of broad characteristics of neutrino oscillations.

Exclusion contours in the ∆m
• For ∆m 2 41 1 eV 2 , the ICAL can exclude active-sterile mixing for |U µ4 | 2 > 0.08. The reach of ICAL in this region is independent of the exact value of ∆m 2 41 . This is the high-∆m 2 41 region where neutrinos from all directions have undergone many oscillations, such that only the averaged oscillation probability is observable, and there is no ∆m 2 41 -dependence.
• Our results for 0.1 eV 2 ∆m 2 41 10 eV 2 are consistent with those obtained in [99]. of magnitude. As a result, interference effects between these frequencies leads to a further improvement in the reach of ICAL in this parameter range: around ∆m 2 41 ∼ 10 −3 eV 2 , the ICAL reach can be up to |U µ4 | 2 > 0.03.
• When ∆m 2 41 10 −4 eV 2 , the active-sterile oscillations due to the frequency ∆ 41 will be suppressed. However in this case ∆ 43 ≈ − ∆ 31 , and the oscillations due to this frequency will continue to be present. These oscillations will be independent of the value of ∆m 2 41 , and hence the reach of ICAL in this lowest ∆m 2 41 range, |U µ4 | 2 0.07, would be virtually independent of the exact ∆m 2 41 value.
• The addition of hadron information, as described in [94], results in a small improvement in the reach of ICAL for excluding |U µ4 |. The improvement is negligible for ∆m 2 41 1 eV 2 , the region where the ∆χ 2 information is almost uniformly distributed in the whole E µ -and cos θ µ -range, and hence additional details of the hadron information do not make much difference. At lower ∆m 2 41 values, though, the hadron information may contribute significantly. At the lowest ∆m 2 41 range, it increases the ICAL reach from |U µ4 | 2 0.07 to |U µ4 | 2 0.05.

Dependence on |U e4 | and the mass ordering scheme
While calculating the exclusion contours in Sec. 4.1, the value of |U e4 | was taken to be zero. A possible nonzero value of |U e4 | could change the results. This change may be discerned from the left panel of Fig. 7, where we show the change in the exclusion contours when the value of |U e4 | 2 is taken to be 0.025. The reach of ICAL clearly improves with a nonzero |U e4 |. The constraints given in the vanishing-|U e4 | scenario are thus the most conservative, and in this paper, we continue to give our constraints in this conservative limit. In case the sterile neutrinos are discovered through the oscillations ofν e at reactor experiments, the measured value of |U e4 | would help in improving the ICAL constraints. The right panel of Figure 7 shows the dependence of the exclusion contours on the mass ordering scheme. The figure suggests that the exclusion reach of ICAL does not depend on the mass ordering scheme for ∆m 2 41 10 −2 eV 2 , while it would depend mildly on the scheme for ∆m 2 41 10 −2 eV 2 .

Contribution of rate, energy and direction measurements
The atmospheric data on charged-current muon neutrinos consist of the event rates, as well as the energy and direction of the muon in these events. It is an instructive exercise, and an overall check of our analysis, to compare the relative importance of these quantities, and interpret them in terms of analytical approximations of oscillation probabilities. In Fig. 8, we show the exclusion contours if certain information were omitted, which would give us an idea of how important that information is for the analysis.
The left panel of Fig. 8 indicates the relative importance of upward-going and downwardgoing events. The right panel of the figure indicates the relative importance of the event rate and the spectral shape (this includes the energy as well as direction information for muons.) The following observations may be made from these figures.
• For ∆m 2 41 1 eV 2 , the information about the active-sterile mixing is slightly more  2 41 plane using the information on both the rate and shape of the muon spectrum, and using the information only on the rate of muon events (integrated over the range [1][2][3][4][5][6][7][8][9][10][11] GeV and over all the directions). The results have been presented for the N-N mass ordering scheme, and using only muon information.
in the downward-going events than that in the upward-going events. Since we know that in this parameter range the direction of neutrinos does not matter, this difference should be simply due to the larger number of downward-going muon neutrinos (since about half the atmospheric ν µ are converted to ν τ on their way through the Earth). The reach of ICAL for |U µ4 | 2 values with the downward-going and upward-going data sets shown in the left panel is indeed observed to follow the approximate number of events in this range. This is also consistent with the observation from the right panel, that in this ∆m 2 41 range all the information is essentially in the event rates, and almost no additional information comes from the energies or directions of the muons.
• For 10 −2 eV 2 ∆m 2 41 1 eV 2 , the left panel shows that the information from upward-going neutrinos is independent of ∆m 2 41 , which may be understood by observing that the upward-going neutrinos in this ∆m 2 41 range still undergo many oscillations inside the Earth, and only their averaged effect is observed at the detectors. On the other hand, the downward-going neutrinos undergo a small number of oscillations ( 1) in the atmosphere of the Earth, and hence the exact value of ∆m 2 41 is relevant. At the higher end of this ∆m 2 41 range, most of the information about sterile mixing is in the downward-going neutrinos, while at the lower end, the information on sterile oscillations in down-going neutrinos is very small since these neutrinos do not get enough time to oscillate.
• For 10 −2 eV 2 ∆m 2 41 1 eV 2 , the right hand panel shows results that support the above observation. At the lower end of this range, the loss of neutrinos due to oscillations to the sterile species is small, as we have also observed in Sec. 3.1. So   [49], IceCube [50], MINOS and MINOS+ [113], SciBooNE and MiniBooNE [47]. The Super-Kamiokande contour, shown with a dotted line, is the result of the analysis using a one-parameter fit, while all the other contours are with analyses using two-parameter fits.
the information in only event rates is insignificant. On the other hand, most of the information is now in the neutrino direction (and hence the muon direction).
• At extremely low ∆m 2 41 values, only the upward-going neutrinos oscillate, and even their number of oscillations is small, so that almost all the information is in the energy and direction distribution of the events.
The analysis in this section shows that the information on muon direction is crucial for ∆m 2 41 1 eV 2 . The accurate muon direction measurements at ICAL [88,90] thus make it a suitable detector for probing low ∆m 2 41 values. Figure 9 shows the projected 90% C.L. reach of ICAL for |U µ4 | 2 with 500 kt-yr of data, and its comparison with the current bounds from the atmospheric neutrino experiments Super-Kamiokande [49] and Icecube [50], the short-baseline experiments CCFR [44], MiniBooNE and SciBooNE [47], as well as the long-baseline experiments MINOS and MINOS+ [113].

Comparison with other experiments
Super-Kamiokande (SK) is the experiment most similar to the INO-ICAL, in the sense that it is sensitive to the atmospheric neutrinos in a similar energy range. In addition to the ν µ → ν µ and ν e → ν µ oscillation channels that ICAL is sensitive to, SK can also detect the charged-current interactions of ν e , and hence can analyze the ν e → ν e and ν µ → ν e channels. (However, in the SK exclusion limits shown in the figure, the information from ν µ → ν µ channel dominates.) It also has a lower energy threshold than ICAL, so it has the advantage of being able to detect the low energy events where the depletion in the number of events is the most prominent for higher ∆m 2 41 (for example, at ∆m 2 41 1 eV 2 , see Fig. 3). On the other hand, it cannot distinguish neutrinos from antineutrinos, and does not have as good muon direction resolution as ICAL, which would play a crucial role in identifying sterile neutrinos for lower ∆m 2 41 values (see Fig. 4). An important difference between the sterile neutrino analysis of SK and the other experiments needs to be noted. Since the SK analysis has been performed only for ∆m 2 41 > 0.1 eV 2 , where the oscillations from all directions are averaged out, the data is not sensitive to the value of ∆m 2 41 . The fit in the analysis is therefore carried out for only one parameter (|U µ4 | 2 ), and therefore the 90% exclusion contour corresponds to ∆χ 2 > 2.71. On the other hand, for all the other experiments including ICAL, the 90% exclusion contours have been evaluated with two-parameter analyses, and hence correspond to ∆χ 2 > 4.61. If the sensitivity analysis of ICAL were to be performed only in the ∆m 2 41 > 1 eV 2 region with a single-parameter fit, the 90% exclusion reach of ICAL would improve from |U µ4 | 2 > 0.08 to |U µ4 | 2 > 0.06.
Some of the most stringent bounds on the sterile mixing parameters are given by IceCube [50]. Indeed for ∆m 2 41 ≈ 0.3 eV 2 , IceCube puts a bound of |U µ4 | 2 0.004. However since IceCube has a lower energy threshold of ∼ 100 GeV, it is sensitive to only higher values of ∆m 2 41 10 −2 eV 2 , as is evident from Fig. 9. The fluxes at the fixed-baseline experiments are more well-determined than the atmospheric neutrino fluxes at ICAL. However ICAL has the advantage of a large range of L/E, and hence sensitivity to a large range of ∆m 2 41 values. Indeed, this paper shows how qualitatively different features of oscillation probabilities come into play for different ∆m 2 41 values, and how ICAL can address all these features so as to be sensitive to a large range of ∆m 2 41 . The CCFR experiment at Fermilab [44] had large neutrino energies ( 100 GeV) and a small baseline (∼ 1 km), and could only probe ∆m 2 41 1 eV 2 . The SciBOONE and MiniBOONE experiments have nearly the same baselines as CCFR, but neutrino energies of ∼ GeV, so they could be sensitive to ∆m 2 41 0.1 eV 2 , but no lower [46]. MINOS [113], having a longer baseline of ∼ 735 km, is sensitive to even lower ∆m 2 41 values, up to ∆m 2 41 10 −4 eV 2 . At ∆m 2 41 10 −2 eV 2 , ICAL would be one of the few experiments sensitive to active-sterile neutrino mixing.

Exploring features of ∆m 2 41
In the previous section, we derived the constraints in the ∆m 2 41 -|U µ4 | 2 space in the scenario where ICAL observes no signal for active-sterile oscillations. In this section, we assume that there is a sterile neutrino with 10 −5 eV 2 ≤ ∆m 2 41 ≤ 10 2 eV 2 , with the mixing parameters |U e4 | 2 = 0.025, |U µ4 | 2 = 0.05 and |U τ 4 | 2 = 0. Note that while exploring the features of ∆m 2 41 in this section, we keep the values of the mixing parameters fixed while generating the data in the 4ν scenario. Even while fitting the data at test values of ∆m 2 41 , we take the active-sterile mixing parameters to be fixed and known. The values of the 3ν mixing parameters have been kept fixed to the values shown in Table 2, and the sign of ∆m 2 31 has been taken to be positive, while generating the data and performing the fit. In this situation, we attempt to address the following two questions in the context of ICAL for the first time: • With what precision can we measure a given value of ∆m 2 41 ?
• Can we identify the sign of ∆m 2 41 ?
In this process, we will also try to quantify the contributions of two of the unique features of ICAL, viz. charge-identification (CID) capability and the hadron energy measurement on the event-by-event basis.

Precision in the determination of ∆m 2 41
In order to illustrate the capability of ICAL to measure the value of ∆m 2 41 precisely, we generate data for an exposure of 500 kt-yr at ∆m 2 41 = 10 −3 eV 2 , which corresponds to the N-N-2 configuration, The quantity is then calculated for test values of ∆m 2 41 in the whole range 10 −5 eV 2 < ∆m 2 41 < 0.1 eV 2 , −0.1 eV 2 < ∆m 2 41 < −10 −5 eV 2 , which gives us an idea of how well a particular value of ∆m 2 41 may be excluded from the data, if indeed there is a sterile neutrino.
The results have been shown in Fig. 10 for three different cases: (i) 3D CID: information on muon charge identification (CID) as well as hadron energy is used; (ii) 2D CID: information on CID is used, but not on the hadron energy; and (iii) 3D without CID: information on hadron energy is used, but not on CID. The muon momentum information is naturally used in all cases. The comparison of these three cases will help us quantify the importance of hadron energy calibration and muon charge identification, two of the unique features of ICAL. The following observations may be made from this figure.
• The value of ∆χ 2 SP indeed shows a sharp dip near the true value of ∆m 2 41 . For the benchmark parameter values, it is expected that the other candidate mass ordering configurations, N-N-1, N-N-3 and N-I, could be ruled out at 3σ in the 3D CID case.
• Note that there is a narrow dip in ∆χ 2 SP near ∆m 2 41 ≈ (3-4)×10 −3 eV 2 . This is likely to be an effect of the interference between the oscillation frequencies corresponding to ∆m 2 41 and ∆m 2 31 .
• It is observed that the values of ∆m 2 41 ≈ −3 × 10 −3 would be strongly disfavoured. The reason for this lies in the observation that in such a scenario, there would have with 500 kt-yr of ICAL data, when the true value is ∆m 2 41 = 10 −3 eV 2 (corresponding to N-N-2). The mixing parameters are taken to be |U e4 | 2 = 0.025, |U µ4 | 2 = 0.05, and |U τ 4 | 2 = 0. The results have been presented with and without muon charge identification, and with and without using the hadron energy information.
been an active-sterile MSW resonance, that would have given rise to significant differences from the actual ∆m 2 41 ≈ +3 × 10 −3 scenario, which has no such resonance.
• For |∆m 2 41 | 10 −2 eV 2 , the ∆χ 2 SP values are independent of the sign of ∆m 2 41 , since the active-sterile oscillations get averaged out at such a high magnitude of ∆m 2 41 , and the information about the sign of ∆m 2 41 is lost.
• In the absence of hadron energy information (2D CID), the value of ∆χ 2 SP reduces by more than 25% at almost all test values of ∆m 2 41 . This illustrates the importance of using the hadron energy information, even though the measurement of hadron energy is not very accurate in ICAL.
• Similarly, the absence of muon charge identification (3D without CID) would also reduce ∆χ 2 SP by more than 25% at almost all test values of ∆m 2 41 . In particular, it would be difficult to rule out the N-N-3 configuration, or the N-N-1 configuration with ∆m 2 41 0.05 eV 2 , to more than 3σ if the CID capability were absent. This illustrates the significant advantage ICAL would have due to its excellent muon charge identification.
• The CID capability would also allow us to analyze the µ − and µ + events (equivalently, ν µ andν µ ) separately, and determine the ∆m 2 41 values with two independent data sets. This would also serve as a test of CPT violation in the 4ν framework. 2)]. The mixing parameters are taken to be |U e4 | 2 = 0.025, |U µ4 | 2 = 0.05, and |U τ 4 | 2 = 0. The results have been presented with and without muon charge identification, and with and without using the hadron energy information.
As has been observed above, ICAL may be able to determine ∆m 2 41 to an accuracy that would enable us to identify the correct m 2 i configuration from N-N-1, N-N-2, N-N-3, and N-I. In particular, for ∆m 2 31 > 0 and ∆m 2 41 ≈ 10 −3 eV 2 , it may be able to identify the sign of ∆m 2 41 . In the next section, we explore this possibility further, focusing on the feasibility of such an identification, with true ∆m 2 41 taken over a wide range.

Determination of the sign of ∆m 2 41
In the last section, we observed that ICAL is expected to be quite sensitive to the mass ordering configurations in the four-neutrino mass spectrum. Earth matter effects would play a significant role in this identification, since whether the active-sterile resonance takes place in the neutrino or antineutrino channel depends crucially on the sign of ∆m 2 41 . Matter effects thus magnify or suppress the effects due to active-sterile mixing, and make the atmospheric data sensitive to the mass ordering in the sterile sector.
In this section, we shall focus on quantifying the sensitivity of ICAL to the sterile sector mass ordering, i.e. the sign of ∆m 2 41 . Fig. 10 already shows that for ∆m 2 41 = 10 −3 eV 2 , the wrong mass ordering would be ruled out to ∆χ 2 10 with 500 kt-yr of exposure. We shall now examine how this result depends on the actual value of ∆m 2 41 . We define the ICAL sensitivity to sterile neutrino hierarchy with the quantity where χ 2 min [∆m 2 41 (wrong sign)] denotes the minimum value of χ 2 when the test value of ∆m 2 41 is varied over all values with the wrong sign. Note that it is not enough to distinguish ∆m 2 41 from −∆m 2 41 , or some other effective value. We would like to ensure that all possible values of ∆m 2 41 with the wrong sign are excluded. The results are shown in Fig. 11 for the three cases 3D CID, 2D CID, and 3D without CID, as earlier. We may observe the following.
• For |∆m 2 41 | 10 −2 eV 2 , there is no sensitivity to the sterile hierarchy. This is expected, since the matter effects are negligible for such large values of ∆m 2 41 • For ∆m 2 41 10 −4 eV 2 , there is no sensitivity to the sterile hierarchy. This may be attributed to the insensitivity of the data to extremely low values of ∆m 2 41 .
• In the intermediate regime 10 −4 eV 2 |∆m 2 41 | 10 −2 eV 2 , ICAL is sensitive to the sign of ∆m 2 41 . The sensitivity grows as the value of |∆m 2 41 | gets closer to 10 −3 eV 2 . This is because at these values of ∆m 2 41 , there is a ν µ -ν s (ν µ -ν s ) resonance inside the Earth due to MSW effects, for ∆m 2 41 < 0 (∆m 2 41 > 0). As a result, the neutrino oscillation probabilities will be widely different with the two signs of ∆m 2 41 .
• In the same intermediate range, three sources contribute to the oscillation probabilities: the oscillation frequencies corresponding to ∆m 2 31 and ∆m 2 41 , as well as the matter potential. The contributions from these sources are of the same order of magnitude in this parameter range. The interference of these three contributions results in the smaller peaks and valleys as observed in the figure.
• In addition, it is observed that for ∆m 2 41 < 0, the ∆χ 2 SH values are much larger than those for ∆m 2 41 > 0. This may be attributed to the facts that the MSW resonance occurs in the ν µ -ν s channel for ∆m 2 41 < 0, and the number of neutrino events at ICAL are almost twice the number of antineutrino events.
• In the absence of hadron energy information (2D CID), the ∆χ 2 SH values reduce by more than 25%, as compared to the 3D CID case, in the intermediate range.
This illustrates the virtue of using the hadron energy information, even though the measurement of hadron energy is not very precise at ICAL.
• Similarly, the absence of CID (3D without CID) would also reduce ∆χ 2 SH by more than 25%, as compared to the 3D CID case, in the intermediate range. This illustrates the advantage ICAL would have due to its excellent CID.
The above analysis reinforces the expectation that in the range |∆m 2 41 | = (0.5 − 5) × 10 −3 eV 2 , the ICAL detector, due to its capacity to estimate the hadron energy and to determine the muon charge, would be highly sensitive to the sterile mass hierarchy.

Summary and Concluding Remarks
The proposed ICAL detector in INO will be able to measure the energy and direction of muons, induced by the interactions of atmospheric ν µ /ν µ in the detector, to a very good precision. It will also be able to determine the muon charge, hence distinguish the incoming muon neutrinos from antineutrinos. Its ability to reconstruct the energy of hadron showers on an event-by-event basis will help it observe neutrinos at multi-GeV energies efficiently. All these features help ICAL to be sensitive to the oscillations of atmospheric neutrinos, and to explore the Earth matter effects on them.
In this paper, we explore the sensitivity of the ICAL detector to the active-sterile neutrino mixing parameters, in the presence of a single light sterile neutrino. If we can confirm the existence of such a light sterile neutrino, it would be a revolution in our understanding of neutrinos, as vital as the discovery of 3ν flavor oscillation. There is no doubt that exploring the properties of such a light sterile neutrino would provide crucial information on the new physics that is being looked for, at the terrestrial neutrino oscillation experiments as well as in astrophysics and cosmology. It is therefore important to explore the sterile neutrino parameter space over a wide range of ∆m 2 41 . The large L/E range scanned by the atmospheric neutrinos makes ICAL sensitive to ∆m 2 41 even as low as 10 −5 eV 2 . An important ingredient of this 4ν scenario is the relative mass squared values, m 2 i 's, of the four neutrino mass eigenstates, i.e. the neutrino mass ordering. While ∆m 2 21 > 0, the signs of ∆m 2 31 and ∆m 2 41 lead to four possible mass ordering schemes (N-N, N-I, I-N, and I-I). Since we would like to explore a wide range of ∆m 2 41 values, we have eight possible mass ordering configurations (that are the subsets of the above four schemes): N-N-1, N-N-2, N-N-3, N-I, I-N-1, I-N-2, I-I-1, I-I-2. One of the main themes in this paper is to quantify the sensitivity of ICAL to the active-sterile mixing parameters, for these mass ordering configurations.
Our analysis is based on realistic energy and angular resolutions and efficiencies for muons induced in CC ν µ interactions, as obtained by the ICAL collaboration through detector simulations. The oscillation channels responsible for the signals at ICAL are ν µ → ν µ andν µ →ν µ , and to a smaller extent, ν e → ν µ andν e →ν µ . The atmospheric neutrino data will therefore be sensitive mainly to the mixing parameter |U µ4 | 2 (which affects all the channels above), with a sub-leading sensitivity to |U e4 | 2 (which only affects the appearance channels) 6 . We explore the sensitivity of ICAL to the parameters ∆m 2 41 and |U µ4 | 2 , while pointing out the changes due to the sub-leading contributions from |U e4 | 2 . In order to achieve this, we solve the neutrino evolution equation in the four-neutrino framework in the presence of Earth matter, using the GLoBES package.
In the scenario where no sterile neutrino is present, we present our sensitivity results in terms of exclusion contours in the ∆m 2 41 -|U µ4 | 2 plane, for an exposure of 500 kt-yr. It is found that the sensitivity reach of ICAL at very high and very low values of ∆m 2 41 is independent of ∆m 2 41 . For ∆m 2 41 1 eV 2 , it is due to the averaging of oscillations; while for ∆m 2 41 10 −4 eV 2 , it is because the oscillations are governed by |∆m 2 43 | ≈ |∆m 2 atm |. In these two ranges, ICAL can exclude the active-sterile neutrino mixing |U µ4 | 2 up to 0.08 and 0.05, respectively, at 90% C.L. with the two-parameter fit in the 3D analysis mode. For intermediate ∆m 2 41 values, the shapes of exclusion contours show more interesting features, which are governed by the matter effects and the interference between atmospheric and sterile mass-squared differences. Here the sensitivity is enhanced: for ∆m 2 41 ∼ 10 −3 eV 2 , values of |U µ4 | 2 > 0.02 can be excluded. In the range ∆m 2 41 0.5 × 10 −3 eV 2 , ICAL can place competitive constraints on |U µ4 | 2 compared to other existing experiments. The exclusion reach is found to be rather insensitive to the value of |U e4 | 2 for ∆m 2 41 10 −3 eV 2 , while at lower ∆m 2 41 , nonzero |U e4 | 2 increases the sensitivity reach of ICAL for |U µ4 | 2 . We also explore the effects of using only the up(down)-going muon event samples to identify which kind of events primarily contribute to the ICAL sensitivity in different ∆m 2 41 regions. As expected, the sensitivity for ∆m 2 41 1 eV 2 is dominated by the down-going events, where neutrinos oscillate even over short distances. For low ∆m 2 41 10 −2 eV 2 , the sensitivity is dominated by the up-going events. We also find that at high ∆m 2 41 , the event rate information is sufficient to achieve the optimal sensitivity, whereas at low ∆m 2 41 , the sensitivity primarily comes from the shape+rate analysis. Including the hadron energy information in the χ 2 analysis also helps to improve the sensitivity marginally. It is observed that the sensitivity reach is more or less independent of the mass ordering scheme assumed (N-N, N-I, I-N or I-I).
If a sterile neutrino exists, then the muon charge identification and hadron energy estimation in ICAL help in the precision measurement of ∆m 2 41 , and even the identification of its sign. We quantify the performance of ICAL in addressing the above important issues for the first time in this paper. For a given true ∆m 2 41 = 10 −3 eV 2 , we explicitly calculate the extent to which test values of ∆m 2 41 can be excluded. We find that the mass ordering configurations N-N-1, N-N-3, and N-I can be ruled out at 3σ for the benchmark parameter values. Since the identification of the mass ordering in the sterile sector would be one of the priorities of any neutrino physics program once a sterile neutrino is discovered, we explore the capability of ICAL to perform this task over a wide range of possible ∆m 2 41 values. We find that for ∆m 2 41 ≈ (0.5 − 5.0) × 10 −3 eV 2 , ICAL has a significant capability for identifying the right mass ordering configuration in the sterile sector.
In summary, ICAL has unique strengths which can play an important role in addressing some of the pressing issues in active-sterile oscillations, involving a light sterile neutrino over a wide mass-squared range. It is one of the few experiments sensitive to ∆m 2 41 10 −3 eV 2 and can put strong limits on |U µ4 | 2 in this range. Its muon charge identification and hadron energy estimation capabilities also help in pinning down the magnitude and sign of ∆m 2 41 in the range (0.5 − 5.0) × 10 −3 eV 2 . We hope that the analysis performed in this paper will take ICAL a step forward in its quest to look for new physics beyond the three-neutrino oscillations.