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 ϵeμ and ϵeτ , while the dependence on ϵμτ and ϵτ τ is seen to be very mild, once the χ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 [42]. 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-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. [5][6][7][8][9].
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. 1 We define ∆m 2 ij = m 2 i − m 2 j .

JHEP12(2015)126
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 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 [19]. Expected constraints from other (future) atmospheric neutrino experiments have been studied previously in the literature, see e.g. refs. [2-4, 11, 18, 25] (see refs. [31][32][33][34][35][36] 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 also expected to 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 [8]. We have not included the hadron energy information in this work. This is being studied in a follow-up work by the INO collaboration. We use the Nuance event generator [43] 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 Geant4 [44] 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. [11]. 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 2, 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 3, we describe the ICAL detector, our simulation techniques, and the statistical analysis procedure. We present our main results in section 4. All results are shown for 500 kton-year of data in ICAL. Subsection 4.1 is devoted to study the impact of NSI on the mass hierarchy -2 -JHEP12(2015)126 sensitivity of ICAL. In subsections 4.2, 4.3, and 4.4, 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 5 with our conclusions.

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 [26,30,37,38] 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 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 [22] | αβ | < 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 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. [20] for a discussion).

JHEP12(2015)126
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 [19,20] | µτ | < 0.033 , The MINOS experiment has also used its data to set the following bound −0.2 < µτ < 0.07 at 90% C.L. [21]. 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: while for the oscillation parameters we assume the following true values: 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 | µτ |. 4 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 figures 1 and 2 the contours of the relative difference A MH αβ between the neutrino oscillation probabilities P αβ (including NSI) 3 We choose a smaller range than the current 90% C.L. allowed range in eq. (2.3), since outside these ranges the χ 2 corresponding to the ICAL data is already large. 4 In fact, note that in the three-flavor neutrino oscillation probability Pµµ the actual parameters are | µτ |, | µµ|, and | τ τ | and appear independently (see eq. (35) in ref. [16]), whereas in the corresponding two-flavor probability they are | µτ | and | τ τ − µµ|. SK has put upper bounds on the two latter parameters using an effective two-flavor analysis. Since we work with the full three-generation case, a way to have the two sets of NSI parameters to coincide is to set µµ = 0.  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. [39])

JHEP12(2015)126
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 figure 1, we show the relative probability difference A MH µµ in the µττ τ plane, keeping eµ and eτ fixed at zero. 5 The probabilities are calculated numerically within the full three-generation oscillation paradigm, using the PREM [24] 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 [14][15][16] Note that eq. (2.8) depends only on the parameters Re( µτ ) and | µµ | − | τ τ | to leading order. 6 Dependence on eµ and eτ appear only at higher 5 Here and throughout the rest of this work, we keep µµ = 0. In fact, using µµ = 0, we effectively reduce the number of NSI parameters by one, i.e., from µτ , µµ, τ τ to µτ , τ τ , and therefore, τ τ and τ τ − µµ are basically the same parameter. 6 The fact that eq. (2.8) depends only on the difference between the two NSI parameters | µµ| and | τ τ | and not on them individually justifies the choice of setting | µµ| = 0. orders, which can be neglected unless these parameters are taken to be large. Therefore, in figure 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 figure 1. The impact of τ τ on the other hand is more involved. Using eq. (2.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 figure 1.
In figure 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. [16] 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.

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, 7 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 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. [5] and using both the muon energy and angle information as well as hadron energy information in ref. [8]. 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.

JHEP12(2015)126
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. [5][6][7][8][9]. 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 [1] 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 are 20% error on flux normalisation, 10% error on cross-section, 5% uncorrelated error on the zenith angle distribution of atmospheric neutrino fluxes, 5% tilt error, and a 5% overall error to account for detector systematics. 8 The individual contributions from µ − and µ + data samples are calculated by minimising eq. (3.2) over the pull parameters. These are then added [cf. eq. (3.2)] to obtain the χ 2 for a given set of oscillation (and JHEP12(2015)126 NSI) parameters. 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. (2.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 2θ 23 , and sin 2 2θ 13 with 1σ errors of 2%, 0.65%, and 0.005, respectively [45] (see also [40,41]). For all results presented in this work, we use 500 ktonyear 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. (3.1) and (3.2), where N ex ij (µ ± ) correspond to the right hierarchy and N th ij (µ ± ) correspond to the wrong hierarchy. We allow only one NSI parameter to be non-zero in both the data and the fit, while all other NSIs are set to zero. The ∆χ 2 MH thus obtained is then marginalised over the oscillation and the NSI parameter that was allowed to be non-zero. 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. At a time we allow only one NSI parameter to be non-zero and marginalise the corresponding ∆χ 2 S obtained over the oscillation parameters and the NSI parameter under consideration. This 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 two-dimensional NSI plane under consideration, while marginalising over the oscillation parameters.  Figure 3. The ∆χ 2 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.

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.

Impact of NSI on mass hierarchy sensitivity
We noted in figures 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. [5]. We revisit the mass hierarchy sensitivity in ICAL in the presence of NSI parameters and show our results in figures 3 and 4. To find the sensitivity of ICAL to the neutrino mass hierarchy, we compute the χ 2 according to eqs.  Figure 4. The ∆χ 2 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 after marginalisation over the oscillation parameters as well as NSI parameters as explained in the text.
hierarchy. We allow one NSI parameter to be non-zero in both the data as well as in the fit. We compute this ∆χ 2 MH as a function of the assumed (true) values of the NSI parameters and show them in the figures 1 and 2. For the sake of simplicity, we take only one nonzero NSI parameter in the analysis at a time, keeping all other NSIs fixed at zero. For instance, the black curve in the top-left panels of figures 3 is obtained as follows. The data are generated for NH and a given true value of eµ (shown as the x-axis). The oscillation parameters in data are taken from eq. (2.6) and all other NSI parameters are set to zero. This is then fitted with a theory prediction corresponding to IH. In figure 3, we present the ∆χ 2 MH obtained when all oscillation and NSI parameters in the fit are fixed at their assumed true values. In figure 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 figures 3 and 4 show the mass hierarchy sensitivity expected in ICAL for the case when there are no NSI considered in JHEP12(2015)126 either the data or the fit. A comparison of this with the black solid and red dashed curves in the figure reveals that presence of NSI in the data could change the mass hierarchy sensitivity of ICAL. In particular, we see that the ∆χ 2 MH changes sharply with the true value of eµ and eτ . In presence of NSI, we note that the ∆χ 2 MH increases for eµ (true) 0 and eτ (true) 0, while it decreases for eµ (true) 0 and eτ (true) 0, compared to what is expected for standard oscillations.
These features can be understood from figures 1 and 2. Since the NSI parameters eµ and eτ mainly affect the appearance channel P eµ , we refer to figure 2 to understand the upper panels of figures 3 and 4. The (0,0) point of figure 2 refers to standard oscillations and gives the mass hierarchy sensitivity shown by the black dashed lines in figures 3 and 4. If we stay on eτ = 0 and change eµ , we note from the left panel (cos θ z = −0.55) of figure 2 that for 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 bottom-right panels of figures 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. [3,4,11]. 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 figure 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 figure 4.

Expected bounds on NSI
In figure 5, we present the expected sensitivity reach of ICAL in constraining NSI parameters. The figure shows the expected sensitivity for each of the NSI parameters eµ (top-left panel), eτ (top-right panel), µτ (bottom-left panel), and τ τ (bottom-right panel). This 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 figure 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  Figure 5. The ∆χ 2 S , giving the sensitivity reach of 10 years of ICAL data in constraining the NSI parameters in the event that the data show no signal of any new physics, as a function of the fit value of NSI parameters. We keep only one αβ (true) to be non-zero at a time, while others are set to zero. The ∆χ 2 S has been marginalised over the oscillation parameters and NSI parameters as explained in the text.
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 figure 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.034 to −0.027 < µτ < 0.027, when we increase the muon energy from 20 GeV to 100 GeV in the data. This corresponds to 18% improvement in the sensitivity. Likewise, the increasing of the energy range improves the constraint on eµ and eτ by 11% and 6% respectively, while for the τ τ there is no significant improvement. The reason is the following. The difference to the muon neutrino survival probability due to NSI (∆P µµ = P NSI µµ − P SM µµ ) keeps changing significantly due to µτ up to 100 GeV, while for the impact of τ τ is negligible beyond 20 GeV. The NSIs parameters eµ and eτ do change ∆P eµ between 20-100 GeV, albeit less. In figures 7 and 9, we show the projected C.L. area in the NSI parameter space allowed after 10 years of running of the ICAL experiment. In figure 7, we show these for the case where we assume that there are no NSI, or in other words, when the assumed true values of the NSI parameters are taken as zero, shown by the black dots in the figure. The black dotted, blue dashed, and red solid contours show the 68%, 95% and 99% C.L. in twodimensional NSI planes. The contours are marginalised over the oscillation parameters after including priors that are described above. The NSI parameters other than the ones appearing in the two-dimensional plane are kept fixed at zero. NH is assumed for all plots. The corresponding contours for the IH case are very similar and we do not repeat them for brevity.  Figure 7. The expected C.L. contours in the given NSI parameter plane. The other NSI parameters are set to zero. The NH has been assumed to be true. The black dots show the points where the data were generated, which are for no NSI in this case.

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 figure 8. 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.   Figure 8. The ∆χ 2 D , giving the discovery potential of 10 years of ICAL data in observing NSI, as a function of the true value of NSI parameters. We keep only one NSI parameter to be non-zero at a time, while others are set to zero. The ∆χ 2 D has been marginalised over the oscillation parameters as explained in the text.

Precision on NSI parameters
In figure 9, we show the C.L. contours for the case where we assume the true values of the NSI parameters to be non-zero. The upper panels show the C.L. in the µττ τ plane when the true values of ( µτ , τ τ ) are taken as (0.02,0), (0,0.05), and (0.02,0.05) for the 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.

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  Figure 9. The expected C.L. contours in the given NSI parameter plane. The other NSI parameters are set to zero. The upper panels are drawn in the µττ τ plane, while the lower panels are drawn in the eτeµ plane. The black dots show the points where the data were generated.
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 [24]. 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. 9 The χ 2 is marginalised over the NSI parameters and the oscillation parameters |∆m 2 31 |, 9 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.

JHEP12(2015)126
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)(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. In 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. (2.3) and (2.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.