Neutrino Physics with Non-Standard Interactions at INO

Non-standard neutrino interactions (NSI) involved in neutrino propagation inside Earth matter could potentially alter atmospheric neutrino fluxes. In this work, we look at the impact of these NSI on the signal at the ICAL detector to be built at the India-based Neutrino Observatory (INO). We show how the sensitivity to the neutrino mass hierarchy of ICAL changes in the presence of NSI. The mass hierarchy sensitivity is shown to be rather sensitive to the NSI parameters $\epsilon_{e\mu}$ and $\epsilon_{e\tau}$, while the dependence on $\epsilon_{\mu\tau}$ and $\epsilon_{\tau\tau}$ is seen to be very mild, once the $\chi^2$ is marginalised over oscillation and NSI parameters. If the NSI are large enough, the event spectrum at ICAL is expected to be altered and this can be used to discover new physics. We calculate the lower limit on NSI parameters above which ICAL could discover NSI at a given C.L. from 10 years of data. If NSI were too small, the null signal at ICAL can constrain the NSI parameters. We give upper limits on the NSI parameters at any given C.L. that one is expected to put from 10 years of running of ICAL. Finally, we give C.L. contours in the NSI parameter space that is expected to be still allowed from 10 years of running of the experiment.


INTRODUCTION
The 50 kton magnetised iron detector (ICAL) to be built at the India-based neutrino observatory (INO) will be mainly observing muon neutrinos coming from the Earth's atmosphere [1]. Amongst the most important goals of this experiment is the determination of the neutrino mass hierarchy (MH) through the observation of Earth matter effects in the expected data sample. For the mass-squared difference ∆m 2 31 > 0 (normal hierarchy) 1 , one expects matter enhanced oscillations in the neutrino channel in the energy range around (5)(6)(7)(8)(9)(10) GeV, while the antineutrino channel does not experience any such matter induced enhancement. On the other hand for ∆m 2 31 < 0 (inverted hierarchy), matter enhanced oscillations are expected in the antineutrino channel, while the neutrino channel does not obtain any such enhancement. ICAL being magnetised will be able to tell its neutrino signal from its antineutrino signal, giving the detector an added handle on measuring these Earth matter effects, and hence, the neutrino mass hierarchy. The sensitivity reach of this experiment for measuring standard neutrino oscillation parameters have been studied extensively in Refs. [2][3][4][5][6].
So far, there has been no signal of physics beyond the Standard Model in any of the accelerator-based experiments including LHC. However, we have unambiguous evidence from complimentary experiments that the Glashow-Weinberg-Salam model of elementary particles is at best a low-energy effective theory and that there exists physics beyond the Standard Model. The Standard Model of particle physics is unable to convincingly explain data from neutrino oscillation experiments. In addition, it also fails to provide explanation of the existence of dark matter and dark energy in the Universe, as well as baryogenesis. There are also theoretical issues with the Standard Model which demand its extension. Any extension of the model could include addition to its gauge or particle sector, or both. It is therefore pertinent to envisage that such an extended theory would also have new (effective) interactions between the particles, beyond what is included in the Standard Model. Such interactions are expected to change the predicted outcome of experiments and existing data can be used to put limits on the strength of these interactions. In this work, we are primarily interested in non-standard interactions (NSI) affecting neutrino oscillations as neutrinos propagate inside 1 We define ∆m 2 ij = m 2 i − m 2 j .
Earth matter. These NSI, if present, would modify the transport of atmospheric neutrinos inside Earth matter, and hence alter the signal at the ICAL detector.
The currently running Super-Kamiokande (SK) atmospheric neutrino experiment has looked for possible presence of these NSI in its event sample, and has found the data to be completely consistent with the Standard Model. Through a statistical analysis, the SK collaboration converts this into an upper bound on the relevant NSI parameters [7]. Expected constraints from other (future) atmospheric neutrino experiments have been studied previously in the literature, see e.g. Refs. [8][9][10][11][12][13] (see Ref. [14][15][16][17][18][19] for earlier works).
In this paper, we will study in detail the impact of NSI on the atmospheric neutrino signal in the ICAL detector at INO. We analyse the prospective data at ICAL in terms of the measured muon energy and muon angle through a binned χ 2 analysis. ICAL is expected to also measure the energy deposited in the associated hadron shower. Inclusion of the hadron energy information improves the energy reconstruction of the events and hence in general improves the sensitivity of ICAL [5]. We have not included the hadron energy information in this work. This is being studied in a follow-up work by the INO collaboration [20]. We use the Nuance event generator with the ICAL detector geometry for generating muons from atmospheric neutrinos. The ICAL energy and angle resolutions and reconstruction and charge identification efficiencies are obtained from the Geant-based detector simulation code developed for ICAL. We generate muon events in the range (1-100) GeV and show the increase in the sensitivity to NSI parameters with the increase of the neutrino energy, and hence the muon energy, as was pointed out in Ref. [10]. We will quantify the extent of this modification in the expected muon signal at ICAL. Using that we will study the expected limits that ICAL could impose on NSI parameters if there is no evidence of NSI in the data.
If on the other hand the NSI parameters are large enough, we would see a signal of new physics in the ICAL data. We give the lower limit on the NSI parameters which is needed in order to allow their discovery in ICAL at any given C.L. Likewise, the presence of NSI could change the sensitivity of ICAL to other neutrino oscillation parameters. In particular, we will show how the NSI parameters alter the mass hierarchy sensitivity in this class of experiments and present the revised sensitivity limits.
The paper is organised as follows. In Section II, we discuss the neutrino oscillation probabilities in the presence of NSI. In particular, we study the impact of NSI on the difference in the probabilities between the NH and IH cases. In Section III, we describe the ICAL detector, our simulation techniques, and the statistical analysis procedure. We present our main results in Section IV. All results are shown for 500 kton-year of data in ICAL.
Subsection IV A is devoted to study the impact of NSI on the mass hierarchy sensitivity of ICAL. In Subsections IV B, IV C, and IV D, we discuss the expected constraints on NSI parameters, the expected lower limit allowing for discovery of NSI, and the allowed areas in NSI parameters space, respectively. We end in Section V with our conclusions.

II. IMPACT OF NSI ON OSCILLATION PROBABILITIES
As outlined in the Introduction, an extension of the Standard Model of particle physics in its gauge sector and/or particle sector is likely to give rise to additional (effective) interactions between Standard Model particles. In particular, in this work, we are concerned with such interactions experienced by the neutrinos when they propagate inside Earth matter.
This effective term in the Lagrangian is of the form [21][22][23][24] where f is a fermion, P C = (1 ± γ 5 )/2 (C = R, L) are the chiral projection operators, G F is the Fermi constant, and ǫ f C αβ are the corresponding NSI parameters. Since Earth matter is made up of the first generation fermions only, the NSI parameters corresponding to e, u, and d are the only ones which contribute towards modifying the neutrino propagation inside the Earth. For the neutral-current NSI what is relevant is the sum ǫ f αβ = ǫ f L αβ +ǫ f R αβ . Furthermore, since only the incoherent sum of the NSI contributions is important, we combine the NSI effects coming from ǫ e αβ , ǫ u αβ , and ǫ d αβ as where n f is the number density of the fermion f and we have normalised the effective contribution to the number density of electrons in Earth. 2 While, in principle, the NSI parameters are complex, we consider only values which have phases either 0 and π. Throughout this 2 Note that the convention followed in defining the NSI parameters is crucial to interpret the actual constraints on them from a given experiment and there are some places in the literature where a different convention has been followed (see Ref. [25] for a discussion).
work, we will use this assumption. Note that we sometimes refer to the NSI parameters as ǫ αβ for simplicity, even though we work with only the real values of these parameters.
Each of the NSI parameters has been constrained from existing data. The corresponding model-independent upper bounds at 90 % C.L. are [26] |ǫ αβ | < where we have arranged the parameters in the form of a matrix with the rows and columns corresponding to {e, µ, τ }. Note that only the NSI parameter ǫ µµ is well constrained in this phenomenological approach, while constraints on all other NSI parameters are rather loose.
In particular, large values of ǫ ee , ǫ eτ , and ǫ τ τ are still allowed. These bounds are generally referred to in the literature as indirect bounds as these bounds on parameters affecting neutrino oscillations come from non-neutrino experiments. The only neutrino experiment that has provided bounds on (some) of these parameters that are better than these indirect bounds is the SK experiment, which puts direct bounds on the NSI parameters |ǫ µτ | and |ǫ τ τ − ǫ µτ | in the framework of the so-called two-flavor hybrid model. From the combined SK I and SK II data sets, the 90 % C.L. upper bounds are given by [7,25] |ǫ µτ | < 0.033 , |ǫ τ τ − ǫ µµ | < 0.147 .
The MINOS experiment has also used its data to set the following bound −0.2 < ǫ µτ < 0.07 at 90 % C.L. [27]. However, this bound is weaker than the one set by SK.
In what follows, we will work with the exact three-flavor neutrino oscillation probabilities and consider the following ranges for the relevant NSI parameters: The NSI parameter(s), if present, will alter the neutrino oscillation probabilities. Oscillograms showing the impact of NSI on the neutrino oscillation probabilities have appeared vastly in the literature. The muon neutrino survival probability P µµ is affected most by the parameters |ǫ µτ | and |ǫ τ τ − ǫ µµ |, while the transition probability P eµ depends on |ǫ eµ | and |ǫ µτ |. This dependence can be used to discover NSI parameters using neutrino oscillation data or constrain them.
It is also known that the dependence of the neutrino oscillation probabilities on NSI parameters is different for the NH and IH cases. Indeed, since measurement of the neutrino mass hierarchy is one of the prime goals of the INO experiment, it is pertinent to ask how the mass hierarchy sensitivity of the experiment alters in the presence of NSI. For the mass hierarchy determination, what matters is the difference in the oscillation probabilities between NH and IH. Therefore, it is appropriate to ask how the difference in the probabilities between NH and IH changes in presence of NSI. In order to show the impact of NSI on the mass hierarchy sensitivity, we present in Figs. 1 and 2 the contours of the relative difference A MH αβ between the neutrino oscillation probabilities P αβ (including NSI) corresponding to NH and IH. We define the relative probability difference A MH αβ as follows (cf. the definition of the total CP-asymmetry in Ref. [28]) where P NH αβ and P IH αβ being the ν α → ν β oscillation probability for NH and IH, respectively. In each case, we calculate A MH αβ for a particular zenith angle θ z , while the energy dependence is integrated out in the range (1-100) GeV.
In Fig. 1, we show the relative probability difference A MH µµ in the ǫ µτ -ǫ τ τ plane, keeping ǫ eµ and ǫ eτ fixed at zero. 4 The probabilities are calculated numerically within the full three-generation oscillation paradigm, using the PREM [29] density profile for the Earth matter.
We compute this for two benchmark zenith angles of cos θ z = −1 and −0.55 corresponding to neutrino baseline lengths of L = 12742 km and L = 7000 km, respectively. The colors represent the contours corresponding to the values of A MH µµ shown in the colorbar. The (0,0) point in the ǫ µτ -ǫ τ τ plane corresponds to neutrino oscillations without NSI (i.e. standard neutrino oscillations). At all other points, NSI are included in the model, and this can be observed to alter the mass hierarchy sensitivity parameter A MH µµ . Note that for standard oscillations we have A MH µµ ∼ −5%. This small relative difference is what the atmospheric neutrino experiments observing ν µ are exploiting to determine the neutrino mass hierarchy.
When the NSI parameters are switched on, A MH µµ changes. The relative difference A MH µµ is seen to increase for ǫ µτ < 0 and decreases further to larger negative values for ǫ µτ > 0.
However, since for the hierarchy measurement what is relevant is the absolute difference |A MH µµ |, we will later see that for all |ǫ µτ | > 0, the hierarchy sensitivity increases as long as all the parameters are kept fixed between the NH and IH cases. The change in A MH µµ with |ǫ µτ | is seen to be in the same direction for both the zenith angles shown, though its magnitude is seen to be larger for cos θ z = −1 case. The change in A MH µµ with ǫ τ τ is less pronounced. In particular, for the cos θ z = −0.55 case, the dependence on ǫ τ τ is marginal. Even for the cos θ z = −1 case, the dependence of A MH µµ on ǫ τ τ for ǫ τ τ = 0 is negligible. For larger values of |ǫ τ τ |, the role of |ǫ τ τ | is to reduce the overall change in A MH µµ due to NSI, and this happens for both positive and negative ǫ τ τ .
In order to understand the change of the probability difference as a function of the NSI parameters, we can series expand the neutrino oscillation probabilities in orders of the NSI parameters and keep only the first-order terms. The expression for the difference in the muon neutrino survival probability between NH and IH, keeping only leading-order terms in NSI parameters and neglecting the standard matter effects, is given by [30][31][32] where A = 2 √ 2G F n e E. Note that Eq. (8) depends only on the parameters Re(ǫ µτ ) and |ǫ µµ | − |ǫ τ τ | to leading order. Dependence on ǫ eµ and ǫ eτ appear only at higher orders, which can be neglected unless these parameters are taken to be large. Therefore, in Fig. 1, we show the relative probability difference for the survival channel in the ǫ µτ -ǫ τ τ plane keeping the other NSI parameters at zero. For ǫ τ τ = 0, the expression clearly shows that |∆P MH µµ |, and hence |A MH µµ |, grows with |ǫ µτ | and flips sign when the sign of ǫ µτ changes. The quantity A MH µµ is positive for ǫ µτ < 0 and negative for ǫ µτ > 0. This agrees fairly well with the exact results shown in Fig. 1. The impact of ǫ τ τ on the other hand is more involved. Using Eq. (8), we note that for any given large value of ǫ µτ , we should have the highest possible |A MH µµ | for ǫ τ τ = 0, and since the dependence on this parameter comes in the form of |ǫ τ τ |, we should have lower |A MH µµ | on both sides of ǫ τ τ = 0. On the other hand, for ǫ µτ = 0, the ∆P MH µµ obtains contribution only from the second term, and there is a relative sign between the two terms in the parentheses. As a result for ǫ µτ = 0 we do not expect large contribution to ∆P MH µµ from NSI. These features can be observed in the exact result in Fig. 1.
In Fig. 2, we present the A MH eµ contours in the ǫ eµ -ǫ eτ plane with ǫ µτ and ǫ τ τ fixed at zero. The probability P eµ depends crucially on the NSI parameters ǫ eµ and ǫ eτ at leading order, and hence, NSI bring significant change to |A MH eµ |. In this case, the corresponding analytic expression is complicated and we refer the reader to Ref. [30] for a related expression for the approximate formula. However, the exact results shown in the figure tell us that the presence of the NSI parameters ǫ eµ and ǫ eτ could bring substantial change to the mass hierarchy sensitivity of atmospheric neutrino experiments. In fact, ǫ eµ and ǫ eτ could either increase or decrease the mass hierarchy sensitivity compared to what we expect from standard oscillations.

III. EVENT SPECTRUM AT INO
The ICAL (Iron CALorimeter) detector at INO will be a 50 kton detector with layers of magnetised iron interleaved with glass Resistive Plate Chambers (RPC), which will serve as the active detector element. The atmospheric neutrinos in ν µ , ν µ , ν e , and ν e species come from decay of pions and kaons produced from cosmic ray interactions with the Earth's atmosphere. These neutrinos can interact with the detector nucleons producing the corresponding charged lepton. The dense iron material of ICAL helps to detect muons through their long tracks 5 , while the magnetic field allows the identification of their charge. Since Earth matter effects develop only in either the neutrino or the antineutrino channel for a 5 The electrons give rise to an electromagnetic shower in the detector, which cannot travel far and is therefore more difficult to observe in this class of detectors.
given mass hierarchy, this charge identification capability gives ICAL an edge to better observe the Earth matter effects, and hence, the neutrino mass hierarchy. The capability of this experiment to help discover the mass hierarchy has been studied before by the INO collaboration using information on muon energy and angle in Ref. [2] and using both the muon energy and angle information as well as hadron energy information in Ref. [5]. In this work, we only use the muon energy and angle information and quantify the change in the mass hierarchy sensitivity of ICAL in presence of NSI. We also study the prospects of constraining or discovering the NSI parameters with the muon event sample.
For calculating the predicted number of µ − and µ + events in ICAL, we use the same prescription as in the earlier INO collaboration papers. The unoscillated events are calculated using the Nuance event generator modified for ICAL. The oscillation probabilities, with and without NSI effects, are implemented through a re-weighting algorithm. Finally, the muon reconstruction efficiency, charge identification efficiency, and muon energy and angular resolutions are folded in as described in Refs. [2][3][4][5][6]. The new ingredient in the simulations performed for this work is that while all the earlier papers used muon sample in the energy range (1-11) GeV, we consider in this work a much larger energy range of (1-100) GeV. In order to do that, we extend the earlier study [33] for muon detector response to 100 GeV from detector simulations done with the GEANT-based code developed for ICAL. The muon energy and zenith angle resolutions, as well as the charge identification efficiency and reconstruction efficiency are obtained as a function of muon energy and zenith angle. This is then folded with the oscillated events to obtain the final muon spectrum expected in ICAL. We generate raw events corresponding to 1000 years of running of ICAL in order to reduce the Monte Carlo fluctuations and normalise the final events to 10 years of running. This event sample is then binned in energy and zenith angle bins as follows. For the zenith angle we have 20 equal size bins in the cos θ z range (−1,1). For the energy, we take variable bin sizes to ensure that there are reasonable number of events in each bin. Between muon energies (1-10) GeV, we take 9 energy bins of size 1 GeV, and between (10-100) GeV, we take 3 energy bins of size 30 GeV.
The predicted events are then analysed by a statistical procedure identical to the one used in the earlier papers. A χ 2 function is defined as where N ′ th ij (µ ± ) and N ex ij (µ ± ) being the predicted and 'observed' number of µ ± events in ICAL, respectively, π k ij the correction factors due to the k th systematic uncertainty, and ξ ± k the corresponding pull parameters. In this analysis, we include 5 systematic uncertainties. These This resultant χ 2 is then marginalised over the oscillation parameters, and when applicable, over the NSI parameters. We assume for the oscillation parameters the true values given in Eq. (6) and marginalise our χ 2 over their current 3σ ranges. We include priors defined as where p true is the assumed true value of the parameter p and σ p the 1σ error on it. We include priors on |∆m 2 31 |, sin 2 θ 23 , and sin 2 2θ 13 with 1σ errors of 1 %, 2 %, and 0.005, respectively [34,35]. For the NSI parameters, the χ 2 is marginalised over their range given in Eq. (5).
For all results presented in this work, we use 500 kton-year of statistics in ICAL. We next define the different χ 2 that we compute in this work for the different physics studies we perform.
Sensitivity to Neutrino Mass Hierarchy: To find the sensitivity of ICAL to the neutrino mass hierarchy, we compute the χ 2 according to Eqs. (9) and (10), where N ex ij (µ ± ) correspond to the right hierarchy and N ′ th ij (µ ± ) correspond to the wrong hierarchy. The ∆χ 2 MH thus obtained is then marginalised over the oscillation and NSI parameters. We do this for different assumed true values of the NSI parameters.
Bounds on NSI parameters: In the event that there is no signal for NSI in the ICAL data, one will be able to give upper bounds on the NSI parameters at a given C.L. In order to find the expected sensitivity of ICAL to constrain the NSI parameters, we compute the χ 2 by generating N ex ij (µ ± ) for standard oscillations (with all NSI parameters set to zero) and fitting this with N ′ th ij (µ ± ) computed with non-zero NSI parameters. The corresponding ∆χ 2 S , where S stands for sensitivity, obtained after marginalising over oscillation gives a measure of the sensitivity reach of ICAL to NSI.

Discovery of NSI parameters:
If on the other hand, one finds a signal of NSI in the ICAL data, this would be a discovery of NSI, and hence, physics beyond the Standard Model. Of course, the NSI parameters in this case have to be above a certain value to be able to produce a discoverable signal at ICAL. We find this lower limit on the NSI parameters needed to be discovered at ICAL for a given C.L. by generating N ex ij (µ ± ) with NSI and fitting this with N ′ th ij (µ ± ) corresponding to standard oscillations. The corresponding ∆χ 2 D obtained after marginalising over oscillation gives a measure of the discovery reach of ICAL to NSI.
Precision on NSI parameters: Finally, for a given set of NSI parameters, one can use the ICAL data to produce C.L. contours in the NSI parameter space. We will show these contours at the 68 %, 95 %, and 99 % C.L. in the ǫ µτ -ǫ τ τ and ǫ eµ -ǫ eτ planes. For this, we will generate N ex ij (µ ± ) for a certain set of NSI parameters and fit it with all values of the NSI parameters in the given plane while marginalising over the oscillation parameters.

IV. EXPECTED RESULTS FROM INO
We now present our main results. We first show the impact of NSI on the mass hierarchy sensitivity of ICAL, which is the main goal of the experiment. We next give the sensitivity reach of this experiment in constraining NSI parameters. Subsequently, we quantify the NSI discovery potential at INO. Finally, we briefly discuss with what precision the NSI parameters could be measured at INO if they were indeed above the discovery limit. MH , giving the expected mass hierarchy sensitivity from 10 years of running of ICAL, as a function of the true value of NSI parameters. We keep only one ǫ αβ (true) to be non-zero at a time, while others are set to zero. The ∆χ 2 is obtained as explained in the text. However, the resultant ∆χ 2 is not marginalised over the oscillation parameters as well as NSI parameters.

A. Impact of NSI on Mass Hierarchy Sensitivity
We noted in Figs. 1 and 2 that the difference in the oscillation probabilities between NH and IH changes in the presence of NSI. Since the mass hierarchy sensitivity is defined in terms of the difference in the event distribution for NH and IH, it is therefore expected that the mass hierarchy sensitivity of the experiment would change in the presence of NSI. The mass hierarchy sensitivity for standard oscillations using only the muon events in ICAL is given in Ref. [2]. We revisit the mass hierarchy sensitivity in ICAL in the presence of NSI parameters and show our results in Figs. 3 and 4. The ∆χ 2 MH corresponds to the difference in χ 2 of the fit with the wrong and the right hierarchy as a function of the true value of the parameters are set to zero. This is then fitted with a theory prediction corresponding to IH.
In Fig. 3, we present the ∆χ 2 MH obtained when all oscillation and NSI parameters in the fit are fixed at their assumed true values. In Fig. 4, we marginalise the ∆χ 2 over the oscillation parameters |∆m 2 31 |, sin 2 θ 23 , and sin 2 2θ 13 with priors. The ∆χ 2 MH is also marginalised over the NSI parameter which is taken as non-zero in the data, while the other NSI parameters are kept fixed at zero. For instance, in the top-left panel, the ∆χ 2 is marginalised over ǫ eµ , while the other NSI parameters are kept fixed at zero. In all cases, the ∆χ 2 MH is marginalised over the oscillation parameters |∆m 2 31 |, sin 2 θ 23 , and sin 2 θ 13 with priors included as described in the previous section. The other panels are also obtained in a similar way.
The horizontal black dashed lines in the four panels of Figs ǫ eµ 0 |A MH eµ | decreases, while for ǫ eµ 0 it increases. This is less clear in the core-crossing bin, however, since the largest mass hierarchy sensitivity at ICAL comes from zenith angle bins close to cos θ z = −0.55, this feature stays in the final ∆χ 2 MH . The effect of ǫ τ τ on the mass hierarchy sensitivity is seen to be less severe from the bottomright panels of Figs. 3 and 4. This NSI parameter affects the muon neutrino survival channel the most. Figure 1 reveals that the impact of ǫ τ τ on |A MH µµ | (when ǫ µτ = 0) is very small for both the core-crossing and the cos θ z = −0.55 bin. The impact of ǫ µτ on the mass hierarchy sensitivity is more interesting and has been discussed in Refs. [10][11][12]. For ǫ τ τ = 0, we see that |A MH µµ | could change up to 20 % for the core-crossing bin and a few percent for the cos θ z = −0.55 bin, as we change ǫ µτ . Note that for standard oscillation |A MH µµ | is already a very small number, and hence, the relative change of |A MH µµ | due to ǫ µτ is significant. This is reflected in the bottom-left panel of Fig. 3, where we see a large increase in ∆χ 2 MH with ǫ µτ . However, once we marginalise over oscillation and ǫ µτ in the fit, this increase is washed out and we obtain no significant impact of ǫ µτ on ∆χ 2 MH in the bottom-left panel of Fig. 4. figure is obtained as follows. We use as data the event distribution at ICAL corresponding to standard oscillations by setting all NSI parameters to zero. This is then fitted with the predicted event distribution which includes one non-zero NSI parameter at a time, and the corresponding ∆χ 2 S calculated. In Fig 5, we show this ∆χ 2 S as a function of the NSI parameter that is allowed to be non-zero in the fit. The ∆χ 2 S is marginalised over the oscillation parameters |∆m 2 31 |, sin 2 θ 23 , and sin 2 2θ 13 . Priors on the three oscillation parameters were included as described in the previous section. The resultant ∆χ 2 S shows little change as a result of marginalisation over them. The black solid curves are obtained when the data are considered corresponding to NH, while the red dashed curves are for data corresponding to IH. We keep the hierarchy fixed to its assumed true value in the fit.
The expected sensitivity for IH is only marginally worse than that for NH. At the 90 % (3σ) C.L., the expected bounds on the NSI parameters from 500 kton-years of statistics in ICAL for NH can be read from the figure as For the IH case, the bounds are comparable and can be read from the figure.
In Fig. 6, we show the improvement that we obtain in the sensitivity reach of ICAL to the NSI parameters when we increase the muon energy range considered in the analysis from 20 GeV (red dashed curve) to 100 GeV (black solid curve). The 3σ bound on ǫ µτ improves from −0.033 < ǫ µτ < 0.033 to −0.027 < ǫ µτ < 0.027, when we increase the muon energy from 20 GeV to 100 GeV in the data.

C. Discovery Reach for NSI Parameters
In the previous section, we looked at how well ICAL will be able to constrain NSI parameters if its data were consistent with just standard oscillations. Next, we take the complementary view and ask ourselves that if NSI parameters were indeed non-zero, at what C.L. would ICAL be able to tell them apart from standard oscillations. In other words, we are looking for the limiting true values of the NSI parameters above which the data at ICAL would be a signal for NSI at a certain C.L. For that, we now consider data for various (assumed) true values of the NSI parameters and fit them with a predicted event spectrum corresponding to standard oscillations and compute the resultant ∆χ 2 D . We present this in Fig. 7. For simplicity, we consider only one non-zero NSI parameter at a time in the data.
We marginalise over the oscillation parameters |∆m 2 31 |, sin 2 θ 23 , and sin 2 θ 13 with priors imposed on each one of them as discussed before. The black solid curves correspond to the case for NH, while the red dashed curves are for IH. We keep the hierarchy to be the same in the theory as in the data.   left, middle, and right panels, respectively. These assumed true points are shown by black dots in the plots. The values ǫ eµ and ǫ eτ are assumed to be zero in both the data as well as in the fit, and we do not show the contours in planes involving these parameters. The lower panels are similar to the upper panels except that now we show the C.L. in the ǫ eτ -ǫ eµ plane when the true values of (ǫ eτ , ǫ eµ ) are taken as (0,0.1), (0.15,0), and (0.15,0.1) for the left, middle, and right panels, respectively. For these panels, the values ǫ µτ and ǫ τ τ are assumed to be zero in both the data as well as in the fit. Again, the figures are for NH, however, the ones for IH are similar and we do not present them for brevity.

V. CONCLUSIONS
The study of the physics potential of the ICAL detector at the proposed India-based Neutrino Observatory is underway. As a part of this on-going effort, we probe in this work the impact of NSI parameters on the expected signal at ICAL and the physics conclusions that one can draw out of it. The neutral-current NSI if present, could alter the propagation of atmospheric neutrinos inside the Earth matter changing the signal at ICAL. This change due to NSI can be used to study the NSI parameters. On the other hand, one needs to estimate how much the potential of ICAL to standard physics gets modified in the presence of NSI. In this work, we have taken both these considerations into account and studied the physics potential of ICAL in presence of NSI.
Measurement of the neutrino mass hierarchy is the primary goal of the ICAL atmospheric neutrino experiment. We showed how the difference in the neutrino oscillation probabilities between NH and IH change in presence of NSI. We defined the relative probability difference A MH αβ for the oscillation channel ν α → ν β and showed the oscillograms for A MH µµ and A MH eµ , the two oscillation channels relevant for the atmospheric neutrinos in ICAL. These oscillograms (and all other results shown in this paper) were obtained from an exact numerical calculation of the three-generation neutrino oscillation probabilities using the PREM profile for the Earth matter density [29]. The oscillograms show that the relative difference A MH µµ changes significantly with ǫ µτ compared to its Standard Model value, while A MH eµ is seen to vary sharply with the values of ǫ eµ and ǫ eτ . The impact of the NSI parameter ǫ τ τ is seen to be less important.
We next simulated µ − and µ + events in the ICAL detector in presence of NSI and defined a χ 2 function, including energy and zenith angle correlated as well as uncorrelated systematic uncertainties, to give C.L. predictions for the estimated sensitivity of ICAL. 7 The χ 2 is marginalised over the NSI parameters and the oscillation parameters |∆m 2 31 |, sin 2 θ 23 , and sin 2 2θ 13 after putting priors on them. Using this we presented the change in ∆χ 2 MH if NSI was assumed to be a certain true value in Nature. We showed that the ∆χ 2 MH increases rapidly for ǫ eµ (true)> 0 and ǫ eτ (true)> 0, while it decreases for ǫ eµ (true)< 0 and ǫ eτ (true)< 0 compared to what we expect for standard oscillations. This behavior can be understood from the oscillograms we showed. The impact of the NSI parameter ǫ τ τ is small, however, the ∆χ 2 MH could vary significantly with ǫ µτ . However, if we allow for marginalisation over the oscillation (especially |∆m 2 31 |) and NSI parameters, the ∆χ 2 MH comes to be around the value predicted by the Standard Model.
We next showed the potential of ICAL in discovering or constraining NSI. If the case that ICAL was consistent with no NSI in the data, we presented the expected upper limit on the NSI parameters. At the 90 % (3σ) C.L. we have for the NH the limits The limits for IH are similar. Compared to the current 90% C.L. bounds given in Eqs. (3) and (4) the expected bounds from ICAL are promising. We next considered the case where the data at ICAL is consistent with NSI and we gave the expected statistical significance with which ICAL will rule out the theory with no NSI. We calculated the range of the NSI parameters that would lead to 90% (3σ) C.L. discovery of NSI at ICAL. Finally, we presented the C.L. contours in the two-parameter NSI planes, for different choices of true values of the NSI parameters. 7 It has been shown that the inclusion of hadron energy information in the analysis of ICAL data improves the mass hierarchy sensitivity of ICAL. The impact of the hadron energy information on the sensitivity of ICAL to NSI is being studied in an independent work [20].