Search for the sterile neutrino mixing with the ICAL detector at INO

The study has been carried out on the prospects of probing the sterile neutrino mixing with the magnetized iron calorimeter (ICAL) at the India-based Neutrino Observatory (INO), using atmospheric neutrinos as a source. The so-called 3 +\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$+$$\end{document} 1 scenario is considered for active–sterile neutrino mixing and lead to projected exclusion curves in the sterile neutrino mass and mixing angle plane. The analysis is performed using the neutrino event generator NUANCE, modified for ICAL, and folded with the detector resolutions obtained by the INO collaboration from a full GEANT4-based detector simulation. A comparison has been made between the results obtained from the analysis considering only the energy and zenith angle of the muon and combined with the hadron energy due to the neutrino induced event. A small improvement has been observed with the addition of the hadron information to the muon. In the analysis we consider neutrinos coming from all zenith angles and the Earth matter effects are also included. The inclusion of events from all zenith angles improves the sensitivity to sterile neutrino mixing by about 35%\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\%$$\end{document} over the result obtained using only down-going events. The improvement mainly stems from the impact of Earth matter effects on active–sterile mixing. The expected precision of ICAL on the active–sterile mixing is explored and the allowed confidence level (C.L.) contours presented. At the assumed true value of 10∘\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10^\circ $$\end{document} for the sterile mixing angles and marginalization over Δm412\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta m^2_{41}$$\end{document} and the sterile mixing angles, the upper bound at 90% C.L. (from two-parameter plots) is around 20∘\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$20^\circ $$\end{document} for θ14\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{14}$$\end{document} and θ34\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{34}$$\end{document}, and about 12∘\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$12^\circ $$\end{document} for θ24\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{24}$$\end{document}.

While data from all of the above mentioned experiments fit nicely into the standard three-generation picture, there are indications (sometimes referred as anomalies) from other neutrino experiments which provide a motivation to extend the three-generation paradigm to include a fourth neutrino mixed with the three standard neutrinos. The first such indication comes from the LSND experiment [25], which showed an excess of electron anti-neutrino events above the expected background with a 3.8σ significance, giving a hint forν μ →ν e oscillations. However, the m 2 required to explain the data through neutrino oscillations is ∼1 eV 2 , making it impossible to fit LSND data along with the solar and atmospheric neutrino data within the three-generation framework. Therefore, in order to explain the LSND results in terms of neutrino oscillations one has to postulate the existence of a fourth (or more) neutrino state(s) which is referred to as "3 + N" model, where '3' stands for active flavors neutrinos and 'N' for sterile neutrino [26][27][28]. On the other hand, from the e + e − collider searches, it has been measured the number of light neutrinos coupled to be 2.92 ± 0.05 [29]. Hence, the additional light neutrino(s) may not couple to the standard model particles through weak currents. They are therefore referred to as sterile neutrinos.
Attempts have been made at confirming or refuting the LSND hint for sterile neutrino oscillations. The KARMEN experiment has earlier looked forν μ →ν e oscillations [30] and no excess observed in their data sample, conflicting with the LSND claim. However, due to background issues, KAR-MEN was unable to rule out the entire LSND region. The MiniBooNE experiment [31] was subsequently built to check the LSND signal but with baseline L and energy E of the experiment changed so that the L/E matched with that of the LSND. The MiniBooNE experiment looked for signatures of ν μ → ν e andν μ →ν e oscillations. In the ν μ → ν e study, the MiniBooNE found no evidence for an excess of ν e candidate events above 475 MeV; however, a 3σ excess of electron like events was observed below 475 MeV, which so far has not been be explained convincingly by any known physics. On the other hand, in their anti-neutrino run, MiniBooNE [32] did see an excess ofν e indicatingν μ →ν e oscillations with m 2 0.1-1.0 eV 2 range, consistent with the evidence for anti-neutrino oscillations from the LSND experiment [25].
In addition to the long-standing LSND results, recently anomalies have emerged in reactor anti-neutrino [33][34][35][36] and gallium-based solar neutrino experiments the GALLEX and the SAGE [37,38]. The reactor anti-neutrino flux calculations were revisited and resulted in a decrease in the ratio of observed to predicted event rate from 0.976 ± 0.024 to 0.943 ± 0.023 leading to deviation from unity at a 98.6% confidence level (C.L.) [33]. This in turn meant that all the earlier reactor anti-neutrino experiments has seen a deficit compared to expectation. Hence, m 2 ∼ 1 eV 2 -driven active-sterile oscillations were put forth as an explanation. The solar neutrino experiments SAGE and GALLEX, during their calibration process with a 51 Cr source of known strength, observed an event rate which was somewhat lower than expected. This again hinted toward oscillations on short distance scales, active-sterile neutrino oscillation with m 2 ∼ 1 eV 2 has been proposed as an explanation. We refer the reader to (SAGE) [39] for a recent status review of these anomalies (another reference including reactor and gallium is [40]).
Cosmological measurements can also provide information on the existence of sterile neutrinos. The sterile neutrinos having non-zero mass could impact the cosmic microwave background (CMB) power spectrum and modify large scale structure formation. The existence of sterile neutrinos which have been thermalized in the early Universe may contribute to the number of relativistic degrees of freedom (effective number of neutrino species). The combined analysis of data from CMB + lensing + baryon acoustic oscillation (BAO) experiments [41] provides a robust frequentist upper limit m ν ≤ 0.26 eV with 95 % C.L. This paper presents the results of an investigation where the muon neutrino flux is used to constrain a possible mixing of a single sterile neutrino with the 3 known neutrinos, viz. the (3 + 1) model. At a high value of m 2 41 where m 2 41 = m 2 4 − m 2 1 , the oscillation probability is averaged out due to the phase part of oscillation probability. Hence, for simplicity, we have considered only the case of "normal" mass hierarchy both for the active and sterile neutrinos. In particular, the study quantifies the sensitivity of up-coming atmospheric neutrinobased magnetized iron calorimeter (ICAL) detector at the India-based Neutrino Observatory (INO) in constraining the active-sterile mixing parameters. The presence of sterile neutrinos with m 2 ∼ 1 eV 2 leads to fast oscillations causing the suppression of the down-going atmospheric neutrinos, otherwise absent in the standard three-generation paradigm. In Ref. [42] the authors studied the constraints on the activesterile mixing expected from atmospheric neutrinos in liquid argon detectors and magnetized iron calorimeters. They carried out their analysis in terms of neutrino energy and zenith angle, assuming fixed values of detector resolutions and efficiencies. In this work, a similar but a more extensive study is carried out. The NUANCE [43] event generator with the ICAL detector geometry has been used to simulate the event spectrum at ICAL. These are then folded with the detector resolution functions and efficiencies obtained by the INO collaboration using the GEANT4-based detector simulation framework developed for ICAL [44,45]. Here, we carried out the analysis in three different ways. First, we binned the data in muon energy and muon zenith angle. Next we took into account the hadron energy information along with the muon energy and muon zenith angle. We also compared the results from these two set of analyses. Finally, we considered the neutrino events coming from all zenith angles in contrast to only down-going neutrino events as used in Ref. [42] and took the Earth matter effects into account as well.
The outline of the paper is as follows. In Sect. 2 we discuss the ICAL detector and its physics goals. The sterile neutrino oscillation formalism including the Earth matter effect is introduced in Sect. 3. The incorporation of detector resolutions using the Monte Carlo method on neutrino induced raw events is given in Sect. 4. The procedure adopted for estimating the oscillated events and the binning scheme considered for estimating the χ 2 using muon energy, its zenith angle and hadron energy are briefly explained in the same section. The definition of χ 2 with pull is given in Sect. 5. The sensitivity to sterile neutrino mixing at an exposure of 1 Mt-yr is discussed in Sect. 6. The detailed results are summarized in Sect. 7.

INO ICAL detector
The 51 kton magnetized ICAL detector, which will use iron as a target, will be placed in an underground cavern, with a minimum all round rock cover of 1 km, in order to reduce the cosmic ray background. The shape and dimensions of the ICAL detector have been decided keeping in mind the cavern dimensions of 132 m × 26 m × 32 m. The ICAL has a rectangular shape with dimensions of 48 m × 16 m × 14.5 m. It consists of three modules each weighing ∼17 kton. The baseline of the ICAL magnet configuration for each of the three modules consists of 151 layers of low carbon steel. The layers are alternated with gaps of 40 mm wherein will be placed active detectors, the resistive plate chambers (RPC), to detect the charged particles produced in neutrino interactions with the iron nuclei. The ICAL RPC detectors give X, Y hit information at ∼0.96 cm spatial resolution, the layer number gives Z-position and timing information with a resolution (σ t ) ∼1 ns [46]. The calorimeter will be magnetized with a piecewise uniform magnetic field (B = 1−1.5 T) [47] thereby being able to distinguish between the oppositely charged particles, μ − and μ + [produced due to charged current (CC) interaction of ν μ andν μ , respectively], from the curvature of their tracks in the presence of the magnetic field.
The ICAL detector is mainly sensitive to muons and hence to the interaction of the ν μ /ν μ . For the electron type neutrino, the detection capability is limited because of the large iron plate thickness (5.6 cm) compared to the radiation length of iron (∼1.76 cm). The production of the tau lepton due to the tau-neutrino interaction is also small due to high threshold for tau production (about 4 GeV). Due to this, the magnetized ICAL detector is most suited to measure muon neutrinos through the tracking of the associated muons and hadrons and reconstruction of their energy and momentum.
The GEANT4-based simulation has been carried out to study the ICAL detector response for the muons and hadrons. At the relevant energies, since the muon is a minimum ionizing particle, a clear track is produced in the ICAL detector. The ICAL simulation program uses the Kalman filter technique, developed by the INO collaboration, to reconstruct the muon track. The typical efficiency of the detector for a 5 GeV muon traveling vertically is about 80%, while the typical charge identification efficiency is more than 95% [44]. The energy of such a muon can typically be reconstructed with an accuracy of about 10%, while its angular resolution is better than 1 • [44]. The energy of the neutrino is the sum of the hadron (E h ) and muon (E μ ) energies for a CC interaction. The hit multiplicity of charged particles distinct from the muon track can be used to estimate the total energy of hadrons in an event. The difference in energies of the interacting neutrino and the outgoing muon, E h ≡ E ν − E μ , has been calibrated against the number of hits in the detector due to the shower produced by hadrons. The measured number of hits can then be used to reconstruct the fractional energy carried by the hadron from the incoming neutrino. From the simulation study, it has been found that the energy resolution is about 85% and 36% for the hadron energies of 1 and 15 GeV, respectively, in the central region of the ICAL detector [45]. Although the energy resolution of hadrons is much lower than that for muons, it still gives an additional information as regards the particular event, which can be used to improve the physics reach of the detector [48,49].
The physics goal of the ICAL detector is to measure precisely the atmospheric neutrino mixing parameters, viz. m 2 31 and sin 2 2θ 23 , [50] and in particular, determine the neutrino mass hierarchy by measuring the sign of m 2 31 [48,51]. In addition, it can be used for octant sensitivity study i.e. whether θ 23 is maximal or not, and if it is indeed nonmaximal, whether θ 23 is less than 45 • or greater than 45 • [49]. The ICAL experiment would also be able to put severe constraints on new physics scenarios like CPT violation [52].

Oscillation probability using + 1 model
The sterile neutrino oscillation probabilities are based on expansion of the 3 generation Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix to 3 + 1 generation, where "3" stands for active and "1" for sterile neutrino, respectively. The neutrino flavors and mass eigenstates are related through ⎛ where U is the mixing matrix. In this analysis the following parametrization has been considered: where R(θ i j ) are the (complex) rotation matrices and θ i j are the mixing angles with i, j = 1, 2, 3, 4; and the order of rotation angles are considered from Ref. [53]. Using the above definition, neutrino flavor change can be described as a function of the mixing matrix elements and neutrino masses in terms of the neutrino oscillation probability, where α, β = e, μ, τ , s; m 2 i j = m 2 i − m 2 j with i > j, and L is the source to detector distance and E is the energy of neutrinos. The order of rotation angles shows that, if all mixing angles are zero, then there is a correspondence between the The global analysis on the sterile neutrino mixing has been carried out in Ref. [54]. The best fit values of the parameters |U e4 | 2 , |U μ4 | 2 and m 2 41 , characterizing the active-sterile neutrino (anti-neutrino) mixing in the 3 + 1 scheme are where m 2 41 = m 2 4 − m 2 1 .

Matter effect on sterile neutrino oscillation
When atmospheric neutrinos travel through matter, they undergo coherent forward scattering on matter. The scattering is mediated by both charged and neutral current processes. The effective mass matrix in matter then changes to where where and Here, A CC and A NC are the matter induced weak charged current (CC) and neutral current (NC) potentials, respectively, which depend on Fermi's constant, G F , matter density ρ, Avogadro number N A , electron fraction Y e in matter and energy of neutrino E. In Eqs. (7) and (8), the "+" and "−" sign corresponds to neutrinos and anti-neutrinos, respectively. Only ν e has CC interaction with electrons and NC interaction with electrons and nucleons, whereas ν μ and ν τ have only NC interaction and sterile neutrinos have no weak interactions. The analysis including matter effect has been carried out considering the varying density profile of the Earth i.e. the preliminary reference Earth model (PREM) [55]. Due to the matter effect, oscillation probability, as given in Eq. 3 will be modified, and in this work the neutrino oscillation probabilities are numerically calculated in matter for exact 3 + 1 generation framework. Figure 1 shows the survival probabilities for neutrinos (left-hand panels) and anti-neutrinos (right-hand panels) for 3 and 3 + 1 generations, as a function of (anti)neutrino energy. The values of the parameters chosen to obtain the oscillation probabilities are given in the figure caption. The black and red lines show the plot for 3 and 3 + 1 generation neutrinos, respectively considering the Earth matter effect. The top panels are for cos θ ν z = 1.0 (down-going neutrinos) while the bottom panels are for cos θ ν z = −1.0 (up-coming neutrinos), where θ ν z is the neutrino zenith angle. We show here the probabilities for the case when only the θ 24 sterile mixing angle is taken as non-zero while θ 14 = 0 = θ 34 . The values of the other oscillation parameters are taken close to their current best-fit values and are given explicitly in the figure caption. The top panels show that while there is no oscillation of neutrinos for the standard 3 generation case whereas the 3 + 1 generations show a small depletion which is the same for neutrinos and anti-neutrinos. This is due to all CP phases being zero and oscillations happening in vacuum. The depletion which comes from active-sterile mixing can be used to constrain the sterile neutrino framework [42]. The bottom panels show the presence of additional Earth matter effects due to active-sterile mixing. The rapid oscillation for 3 + 1 generation neutrinos is due to the | m 2 41 |-driven phase factor of the oscillation probability. The ICAL detector is not expected to decipher these fast oscillations. For this cos θ ν z the neutrino and anti-neutrino survival probabilities have a difference coming mainly from m 2 31 -driven earth matter effects, while the impact due to sign of m 2 41 is hardly evident for this case where the mixing angles θ 14 and θ 34 are taken as zero. Figure 2 shows the appearance P eμ (top panels) and survival, P μμ (bottom panels), neutrino oscillation probabilities for 3 and 3 + 1 generations, as a function of neutrino energy. The values of the mixing parameters chosen to obtain the oscillation probabilities are given in the figure caption. In particular, one can see that here we allow the sterile mixing angles θ 14 and θ 34 to be non-zero and study the impact of these on the oscillation probabilities. The left-hand panels have been obtained for a baseline of L = 7000 km (corresponding to a zenith angle of cos θ ν z = −0.55) and righthand panels are for a baseline of L = 10,000 km (corresponding to a zenith angle of cos θ ν z = −0.785). The black lines show the plot for the 3 generation case and red, green and blue lines represents for the 3 + 1 generation neutrinos case. The red lines correspond to the case where only θ 24 is taken as non-zero (sin 2 θ 24 = 0.03 here). This is similar to the plots shown in Fig. 1. We next show in green lines the change obtained when we make θ 14 non-zero. We notice a small difference between the red and green curves. Finally, in the blue lines we take all the 3 sterile mixing angles to be non-zero and their values are given in the figure caption. We see that the oscillations probabilities change significantly when θ 34 is switched-on. Note that the oscillation probabilities P eμ and P μμ are independent of the value of θ 34 as per the parametrization of the mixing matrix, and the effect of θ 34 comes entirely due to earth matter effects [56]. Both oscillation probabilities, P eμ and P μμ , also show a shift in the position of the maxima or minima with inclusion of θ 34 . We notice from the figure that while each of the sterile mixing parameters changes the neutrino oscillation probabilities, the impact of θ 34 appears to be most dramatic and changes the shape of the oscillation probabilities. It may be noted that in this analysis both the disappearance (ν μ → ν μ ) and the appearance (ν e → ν μ ) oscillation channels are considered while estimating the oscillated muon events for the atmospheric neutrinos.

Simulation technique
In this section we describe the simulation of atmospheric neutrino events in ICAL. The generation of the event spectrum and binning scheme are presented in Sects. 4.1 and 4.2, respectively.

Event reconstruction
The study of sterile neutrino sensitivity has been carried out using the detailed information of neutrino induced events. The atmospheric neutrinos, interacting with the iron target, produce leptons and hadrons through the charge current interaction. There are mainly three processes which contribute to the CC interactions in the ICAL detector. At the sub-GeV energy range, the quasi-elastic process dominates, in which the final state muon carries most of the available energy. As the energy increases from sub-GeV to multi-GeV, hadrons and their showers are produced in resonance (RS) and deepinelastic scattering (DIS) processes. In the RS process, most of the hadron showers consist of a single pion, though multiple pions may contribute in a small fraction of events. On the other hand, in the DIS process multiple hadrons are produced which carry a large fraction of the incoming neutrino energy. In the CC interaction of neutrino, for every hadron shower, there is a corresponding muon coming from the same CC interaction vertex of neutrino.
Again, as shown in Fig. 1, the probability of active-sterile oscillations is valid only for down-going neutrinos, if one uses vacuum oscillations. However, the inclusion of the Earth matter effects for the upward going neutrino may contribute Table 1 Production height (slant distance in km) of neutrinos and its corresponding sigma for six values of cosθ z at three-neutrino energy ranges [58] cosθ to the active-sterile oscillations. In what follows, we will consider neutrino events covering all zenith angles and present a comparison of the improvement coming from up-coming events.
The raw events without oscillations are generated using the NUANCE neutrino event generator modified for the ICAL detector. To reduce statistical fluctuations, initially data are generated for the 1000 years and further normalized to the required exposure during the statistical analysis. In this analysis, an exposure of 1 Mt-yr has been considered. Neutrino oscillations are then introduced through a reweighting algorithm as discussed in [57], which is based on the acceptance-rejection method. While calculating the neutrino oscillation probabilities, we consider the ν production height distribution in the atmosphere as a function of zenith angle and energy of neutrino which is given in [58] only for the down-going neutrino. The path length (L) distribution was incorporated on MC basis by considering Gaussian smearing of the path length whose mean is L and the corresponding zenith angle and energy dependent sigma as given in Table 1. It is to be noted that we have considered the production altitude distribution for down-going neutrinos only. This is because the sterile neutrino oscillations are relevant on small length scales and the uncertainty in the production point for the down-going neutrinos, therefore, becomes crucial. Figure 3 shows the ratio of events with and without oscillations taken into account as a function of neutrino energy. The results are shown for the neutrinos and a particular zenith angle bin of cos θ ν z = 0.45−0.5. Figure 3a Table 1 and σ θ ν z = 10 • in this illustrative figure. It may be noted that the energy and the zenith angle correspond to true values obtained from the NUANCE output. On the other hand, the reconstructed energy and zenith angle correspond to the neutrino events after incorporating detector resolutions. The resolution functions bring a mild smearing of the shape of the event spectrum. Finally, the blue lines are obtained after smearing for the production point in the atmosphere as well according to Table 1. The impact of this smearing is rather large when the data is sensitive to the phase of the m 2 41driven oscillations (left panel). Since the production point becomes uncertain, the path length becomes uncertain, causing a drop in the effective dip of the event spectrum due to active-sterile oscillations.
The ICAL simulations are performed not in terms of the neutrino energy and angle, but in terms of muon energy and angle and hadron energies. Hence the oscillated events are distributed two dimensionally in terms of the final state muon energy and muon zenith angle bins. The hadron events are binned in hadron energy only. The detector energy resolution, angle resolution, charge and reconstruction efficiencies are next incorporated in the simulations. We have incorporated the actual resolutions of the detector for muons and hadrons using INO look up tables [44,45]. The muon energy and zenith angle resolutions are a function of both energy as well as zenith angle. To incorporate the detector response in the true event distribution numerically, we have used the Monte Carlo (MC) methods. At low energy (E < 0.9 GeV), the energy loss of muons in the detector follow a Landau probability distribution function (p.d.f.) and above this a Gaussian p.d.f. For a particular event, the Landau p.d.f. (P L ) for E μ < 0.9 GeV and Gaussian p.d.f. (P G ) above 0.9 GeV have been used to smear the true energy of muon. To incorporate the energy dependent zenith angle resolution of muon, the P G has been used for all energy. The mean of the function are true E μ , cos θ z for energy and zenith angle of muons, respectively, without incorporating the ICAL detector resolutions. The final energy and cosine of the zenith angles are as follows: cos r θ z = P G (cos θ z , σ cos θ z ) (11) where σ E μ is the standard deviation of energy, σ cos θ z is the standard deviation for cosine of zenith angle taken from [44], and E r μ , cos r θ z are reconstructed energy and cosine of zenith angle of muons, respectively. To incorporate the hadron energy resolutions in the ICAL simulation, the Vavilov probability distribution function has been used as described in [45].
In Fig. 4  The expected rate of neutrinos assuming no sterile neutrinos as a function of cosine of zenith angle (cosθ z ) and energy (E μ ) is shown in Fig. 5a. The relative reduction in expected neutrinos events rate for the mixing parameters around the best fit [54] is shown in Fig. 5b. It may be noted here that the event distributions are obtained after incorporating the detector resolution as well as detector reconstruction efficiencies.

Binning scheme
After incorporating the ICAL detector resolutions for muons and hadrons in neutrino induced events, a variable binning scheme has been adopted to study the sensitivity of sterile neutrinos oscillation. In a simulation the bin content may be Table 3 The binning scheme adopted for the reconstructed parameters E μ , cos θ z and E h for muons and hadrons, respectively, where χ 2 is estimated considering combined muon and hadron information of down-going neutrino induced events less than one, but in an experiment, the observed data in various bins should be greater than or equal to one. The bins of different widths are chosen in order to ensure that there Table 4 The binning scheme adopted for the reconstructed parameters E μ and cos θ z for muons produced due to both up-coming and downgoing neutrinos, for the χ 2 estimations with muons only is at least one event in each bin. The binning schemes for down-going neutrinos induced events are given in Tables 2  and 3. Table 2 gives the binning scheme for the down-going events using only the muon energy and zenith angle. Table 3 on the other hand gives the binning scheme for the analysis with the down-going events where the hadron energy information is used in addition to muon energy and zenith angle.
Since the muon reconstruction efficiency is nearly zero for near-horizontal bin, we set a lower limit of cos θ z > 0.1. The production of atmospheric neutrino flux follows a steep power law in energy (∼E −2.7 ), resulting in a smaller number of events at higher muon and hadron energies. Therefore, finer bins at low energies and wider bins at higher energies are considered for both muons and hadrons, respectively. Moreover, it may be noticed that at low energies (E = 1−1.5 GeV) and high zenith angle (cos θ z = 0.1−0.2), a larger bin width is considered due to the fact that the reconstruction efficiency of ICAL is poor in this region. Similar arguments were adopted in choosing the binning scheme when considering events due to neutrinos coming from all directions. The details of the bin combinations for the cases when only muon energy and its zenith angle are used and when hadron energy information is used as well, are given in Tables 4 and 5, respectively. It is to be noted that wider bin combinations are considered for up-coming neutrinos induced events due to the further depletion of events as a result of oscillation of neutrinos in these cases.

χ 2 estimation
To incorporate the muon energy and angle information in the analysis, we define χ 2 as [59] and R th n 1 ,n 2 = R th n 1 ,n 2 where n 1 and n 2 are the numbers of bins for energy and cosine of the zenith angle for muons, R ex n 1 ,n 2 , R th n 1 ,n 2 are observed and theoretically predicted events, π i n is the strength of the coupling between the pull variable ξ i and R th n 1 ,n 2 , which carries the information as regards systematic uncertainties as given in Eq. (13). Equation (12) is minimized with respect to pull variables. Five systematic uncertainties, viz. an overall flux normalization error of 20%, overall normalization of cross-section of 10%, flux tilt factor of 5%, which takes into account the deviation of the atmospheric fluxes from a power law, zenith angle dependence of the flux of 5% and finally an overall 5% systematic error are considered. More discussion of implementation of systematic uncertainty can be found in Ref. [57]. The information on the hadron energy is incorporated along with the muon energy and angle information by defining the χ 2 as 2(R th n 1 ,n 2 ,n 3 − R ex n 1 ,n 2 ,n 3 ) where the index n 3 runs over the hadron energy bins. The other variables are the same as defined before. It may be noted that the theoretically predicted and observed events correspond to with and without sterile neutrino oscillated events, respectively.
The χ 2 so defined is next marginalized over the sterile neutrino oscillation parameters m 2 41 and the mixing angles θ 14 , θ 24 and θ 34 . The total χ 2 , is estimated by adding the priors of m 2 41 and θ 14 , θ 24 and θ 34 , The parameters sin 2 2θ i j and m 2 41 are varied within the range of 0.047-0.22 and 0.82-2.19 eV 2 , respectively. The error on sin 2 2θ i j and m 2 41 are considered as 10% of their assumed true values.

Exclusion limits
The oscillation probabilities of atmospheric neutrinos depend on the active-sterile neutrino mixing. The data from ICAL will therefore be sensitive to this mixing and should be able to constrain them. The upper limit for the sterile neutrino mixing angle θ 24 for an exposure of 1 Mt-yr is shown in Fig. 6. For simplicity we take the mixing angles θ 14 and θ 34 to be zero in this plot. The black lines correspond to the sensitivity obtained in the analysis when only muon energy and angle information is used, while the red lines correspond to the expected sensitivity when the hadron energy together with muon energy and angle information are included as well. A comparison of the black and red lines shows that the addition of hadron energy information in the analysis does not bring any significant improvement in the sensitivity. We show the expected sensitivity for the case where only down-going neutrinos are considered as well as when data from all zenith angles are analyzed. The lines labeled as 'D' in the figure show the expected exclusion plots where only down-going events are considered, while the lines labeled as 'UD' in the figure show the expected sensitivity when neutrino events from all zenith angles are included in the analysis. The inclusion of up-going neutrinos is seen to improve the sensitivity to the sterile mixing angle θ 24 . The improvement comes from the increase of statistics as well as from the inclusion of earth matter effects which change with the inclusion of activesterile mixing. With 1 Mt-yr data, ICAL is expected to limit sin 2 2θ 24 < 0.19 at the 90% C.L. for m 2 41 ∼ 0.1 eV 2 using the down-going events only. This limit is expected to improve to sin 2 2θ 24 < 0.12 when events from all zenith angles are included. This constitutes an improvement of about 35%.
We also show in Fig. 6a, the exclusion plots from other experiments for comparison. The green dashed-dotted line in Fig. 6a shows the exclusion plot from the Sci-BooNE/MiniBooNE [60] experiment which looked for ν μ disappearance. It is found that at lower values of m 2 41 , the ICAL detector has better sensitivity compared to SciBooNE/ MiniBooNE due to the longer path length and lower energies of the atmospheric neutrinos. In addition the violet line shows the result from IceCube [62]. Figure 6b shows the impact of individual systematics on the active-sterile mix- Fig. 7 Comparisons of exclusion limits in the m 2 41 − sin 2 2θ 24 plane for rate-only, shape-only and including all systematics with consideration of down-going neutrinos only ing sensitivity. It is found that the uncertainty due to flux normalization has the maximum effect compared to the others. In addition the exclusion limits are obtained performing the rate-only analysis considering single bin in energy and keeping the zenith angle bins the same as given in Table 2 and also shape-only analysis by removing the pull term on neutrino flux from the chi-square function where we consider the down-going neutrinos only. Figure 7 shows the comparison of exclusion limits obtained from shape-only (red dashed line) and the rate-only (violet line) analysis at 90% C.L. The lower limit on sin 2 2θ 24 ∼ 0.57 at m 2 41 = 0.1 eV 2 obtained from the rate-only analysis. However, while doing the rate only analysis by merging all the angle and energy bins, the limit on sin 2 2θ 24 ∼ 0.70 at m 2 41 = 0.1eV 2 at 68 % C.L. It has been observed that the exclusion limits from shapeonly analysis (sin 2 2θ 24 ∼ 0.59) is less sensitive compare to the results obtained from the analysis of full binned data and including the pull corresponding to all systematic errors (sin 2 2θ 24 ∼ 0.47) at higher values of m 2 41 = 10.0 eV 2 . However, at low value of m 2 41 = 0.1eV 2 the difference is ∼3 %. Figure 8a shows the expected sensitivity of ICAL to the effective mixing angle θ eμ , which for the active-sterile mixing case is defined as sin 2 2θ eμ = 4U 2 e4 U 2 μ4 . We have taken a statistics corresponding to 1 Mt-yr data at ICAL and considered the information on just the muon energy and zenith angles. We show the 90% C.L. expected exclusion limits in the m 2 41 -sin 2 2θ eμ plane for fixed choices of θ 14 [64,65]. The sensitivity of the ICAL experiment is seen to increase significantly depending on the value of θ 14 and could in principle rule out significant parts of the LSND allowed region. The blue dashed line shows results obtained from combined analysis of MINOS and Daya Bay/Bugey-3 experimental data at 90% C.L. [66]. Figure 8b shows the marginalization effect of θ 14 on the sensitivity sin 2 2θ eμ . The 1σ range of θ 14 was considered from Ref. [67]. It is observed that the impact on sin 2 2θ eμ is negligible. Further exclusion plots are generated for various mixing angle combinations considering only muon energies and zenith angles and shown in Fig. 9. The data for these figures are generated for assumed true values of sterile mixing angles of θ 14 = θ 24 = θ 34 = 10 • and at m 2 41 = 1 eV 2 . The assumed true values of the standard oscillation parameters are given in the figure caption. Figure 9 shows the expected constraints from 1 Mt-yr of ICAL data  The point at which the data is generated is shown by the black dots. The different lines show the different marginalizing methods, with added priors and at different C.L., and the details are given in the figure legends. For the cases which include marginalization, we have marginalized the χ 2 over m 2 41 and the remaining sterile mixing angle, viz., θ 34 in the left-hand panel, θ 24 in the middle panel and θ 14 in the right-hand panel, with priors, as was discussed in the previous section. For the case with marginalization, we expect an upper bound at 90% C.L. (from two-parameter plots) of around 20 • for θ 14 and θ 34 , and about 12 • for θ 24 . On the lower side only θ 24 is seen to be bounded at 90% C.L. without priors and at 68% C.L. once priors are included. While for both other mixing angles, the cases θ 14 = 0 and θ 34 = 0 are expected to be compatible with the data at even the 1σ C.L. for all analyses that we considered.

Summary
The sensitivity of ICAL to active-sterile neutrino mixing is studied. The effect of inclusion of hadron energy into the analysis has been probed. The present analysis shows that the inclusion of up-coming neutrinos, which are affected by the intervening matter, improves the sensitivity to sterile neutrino mixing. The inclusion of hadron energy information in the analysis, on the other hand, has almost no effect on the sensitivity to sterile neutrino mixing. From the down-going events alone, one expects an upper bound of sin 2 2θ 24 < 0.16 at 90% C.L. from 1 Mt-yr of data. The sensitivity to the sterile neutrino mixing angle further improves by considering neutrinos coming from all directions at the detector. There is enhancement in sensitivity by about 35% for the whole range of m 2 41 . At lower values of m 2 41 , the ICAL detector has better sensitivity compared to the short baseline experiments like SciBooNE/MiniBooNE. Also the expected bounds on the sterile mixing angles are obtained. For an illustrative case of assumed true value of 10 • for all the three sterile mixing angles and after marginalization, an upper bound at 90% C.L. (from two-parameter plots) of around 20 • for θ 14 and θ 34 , and about 12 • for θ 24 is obtained. The impact of the inclusion of priors on the sterile neutrino parameters is studied. Only for θ 24 could one rule out the zero-mixing possibility at 90% C.L. for this illustrative case, while for both other mixing angles θ 14 = 0 and θ 34 = 0 are compatible with the data at even the 1σ C.L.