Search for light sterile neutrinos with the T2K far detector Super-Kamiokande at a baseline of 295 km

We perform a search for light sterile neutrinos using the data from the T2K far detector at a baseline of 295 km, with an exposure of 14.7 (7.6)$\times 10^{20}$ protons on target in neutrino (antineutrino) mode. A selection of neutral current interaction samples are also used to enhance the sensitivity to sterile mixing. No evidence of sterile neutrino mixing in the 3+1 model was found from a simultaneous fit to the charged-current muon, electron and neutral current neutrino samples. We set the most stringent limit on the sterile oscillation amplitude $\sin^2\theta_{24}$ for the sterile neutrino mass splitting $\Delta m^2_{41}<3\times 10^{-3}$ eV$^2/c^4$.

We perform a search for light sterile neutrinos using the data from the T2K far detector at a baseline of 295 km, with an exposure of 14.7 (7.6)×10 20 protons on target in neutrino (antineutrino) mode.A selection of neutral current interaction samples are also used to enhance the sensitivity to sterile mixing.No evidence of sterile neutrino mixing in the 3+1 model was found from a simultaneous fit to the charged-current muon, electron and neutral current neutrino samples.We set the most stringent limit on the sterile oscillation amplitude sin 2 θ24 for the sterile neutrino mass splitting ∆m 2  41 < 3 × 10 −3 eV 2 /c 4 .

I. INTRODUCTION
Over the last few decades, the theory of neutrino oscillations has been well established through a series of experiments with neutrinos produced by the Sun [1][2][3][4][5][6], nuclear reactors [7][8][9][10], accelerators [11][12][13][14][15] and in the atmosphere [16,17].Most data from these studies are consistent with the three flavor paradigm where the three weakly interacting neutrino flavors are related to three neutrino mass states by the PMNS mixing matrix [18][19][20].However, deviations from the three flavor scheme have been reported.At LSND [21] and MiniBooNE [22], there were excesses of νe found in short-baseline νµ beams; MiniBooNE also reported an excess in ν e appearance [23]; radioactive calibration sources in gallium experiments [24,25] showed a deficit of ν e flux; and reactor experiments [26] observed less νe than expected.These results could be explained by a fourth neutrino state with a mass difference ∆m 2 ∼ 1 eV 2 /c 4 with respect to the three PMNS states [27][28][29][30].From the measurements of the invisible decay width of the Z 0 boson at the LEP collider, the number of weakly interacting neutrino species with mass below 45 GeV/c 2 is limited to three [31], so the new neutrino state must not couple to the weak interaction and is often referred to as sterile.We can incorporate this additional neutrino state in the simple "3+1" model [28], which involves the three active neutrinos and one sterile neutrino, and study its effect on the oscillation signatures.Currently, the null results, especially in the ν ( -) µ disappearance channels, from short-baseline accelerator experiments like CCFR [32], MiniBooNE/SciBooNE [33] and T2K [34]; long-baseline experiments like MINOS [35] and NOvA [36]; or atmospheric experiments like Super-Kamiokande [37] and Ice-Cube [38,39], have limited the available parameter space in the "3+1" model.
The Tokai to Kamioka (T2K) experiment is a longbaseline accelerator neutrino experiment in Japan which primarily measures muon neutrino disappearance and electron neutrino appearance.
While T2K is designed for studying standard three flavor oscillation at ∆m 2 ∼ 10 −3 eV 2 /c 4 , it also has the potential to search for oscillation signatures due to sterile neutrinos around this ∆m 2 range.Neutral current (NC) neutrino interactions are also collected in the far detector, Super-Kamiokande (SK), which can be used to enhance the sensitivity to sterile mixing as the sterile neutrinos, unlike other active neutrinos, do not interact through CC or NC scattering.We present a long-baseline search for sterile neutrinos in the "3+1" framework, using both the charged-current (CC) ν ( -) µ and ν ( -) e samples and NC samples at the far detector.
Sec. II briefly describes the sterile neutrino mixing model and its effect on the oscillation probability (or oscillation signatures).The T2K experimental setup is outlined in Sec.III, followed by event selection criteria in Sec.IV.Sec.V explains the analysis strategy and Sec.VI presents our search results.Finally, Sec.VII gives a summary and outlook of our sterile neutrino study.

II. "3+1" STERILE NEUTRINO MIXING
In this study, we focus on a "3+1" like model where a single sterile neutrino is added and mixed with the three active states, which is the simplest model with sterile neutrino frequently used in neutrino oscillation analysis.In this model, there is a new flavor state ν s and a new mass state ν 4 with mass m 4 added to the three flavor framework.The relation between the flavor and mass states is given by where |ν α are the flavor states and |ν k are the mass states.The original 3 × 3 PMNS mixing matrix is expanded to a 4 × 4 matrix as: We choose the parameterization as in [40]: where U ij is a unitary rotation matrix of an angle θ ij in the ij-plane.There are therefore three new mixing angles θ 14 , θ 24 , θ 34 and two new CP-violating phases δ 14 , δ 24 .Note that sin 2 θ 14 has been constrained to small values by reactor experiments [41], and T2K has limited sensitivity to θ 14 and the new CP phases.Since there is no significant correlation between them and the other oscillation parameters in this study, we set θ 14 = δ 14 = δ 24 = 0 to simplify the mixing matrix.At the far detector, the ν µ survival probability can be approximated (omitting δ CP terms) as: and the ν e appearance probability as: (5) Thus the CC channels are sensitive to θ 24 and ∆m 2  41 .Similarly, the active neutrino survival probability, which is manifested in the NC channel, is sensitive to θ 24 , ∆m 2  41 , and θ 34 : where A = cos θ 24 sin θ 34 and B = sin 2θ 24 cos θ 34 .The antineutrino oscillation probabilities follow similarly, but there are small differences due to the δ CP terms which are not explicitly written here.Figure 1 shows schematically how the oscillation probabilities are modified with the mixing of sterile neutrinos.The addition of a sterile neutrino state which experiences no matter potential (both CC and NC) complicates the calculation of oscillation in matter.We estimated the possible impacts due to matter effects at T2K baseline and energy with numerical calculation of oscillation probabilities using nuCraft [42], and found at most a few percent changes on the ν ( -) e and NC samples.This is negligible at current precision, so we simply use the vacuum oscillation probabilities (without approximation) in this study.

III. THE T2K EXPERIMENT
The T2K experiment [43] consists of a neutrino beam, a near detector complex, and the water Cherenkov detector Super-Kamiokande [44] as the far detector at a baseline of 295 km.T2K is sensitive to θ 13 and δ CP through FIG. 1: The muon neutrino survival probability (top), electron neutrino appearance probability (middle) and active neutrino survival probability (bottom) as a function of neutrino energy, with and without sterile neutrino, assuming ∆m 2 31 > 0.
the ν ( -) e appearance searches, and can also make precision measurements on θ 23 and mass difference ∆m 2  32 with the ν ( -) µ disappearance samples.
The main ring accelerator in the Japan Proton Accelerator Research Complex (J-PARC) produces a 30 GeV proton beam with spills every 2.48 s that contain eight beam bunches which are 580 ns apart.The proton beam is extracted to the neutrino beamline and collides with a graphite target to produce charged pions and kaons.In the neutrino (antineutrino) beam mode, magnetic horns are used to focus the positively (negatively)-charged pions and kaons into a 96 m long decay volume filled with helium gas.These mesons typically decay into muon neutrinos (antineutrinos).The neutrino beamline is directed at an angle of 2.5 o away from the far detector, so that the off-axis beam at SK has a narrower peak at 0.6 GeV with much less ν ( -) e content than an on-axis beam.The near detector complex located at 280 m from the neutrino production target contains two detectors.The on-axis detector, INGRID, is an array of 16 iron/scintillator detectors which precisely measure the beam direction and intensity [45].The off-axis detector, ND280, is a magnetized tracking detector which constrains the neutrino flux and cross-section model parameters in the oscillation analysis.The far detector, Super-Kamiokande, is located in Gifu prefecture, at a distance of 295 km from the J-PARC neutrino beam.It is a water Cherenkov detector consisting of 50 kt of ultra-pure water.The tank is optically separated into two regions.The inner detector (ID) is a cylinder containing 32 kt of water and is instrumented with 11,129 inward-facing 20 inch photomultiplier tubes (PMTs).The outer detector (OD) extends 2 m outward from the ID and is instrumented with 1885 outwardfacing 8 inch PMTs.The OD serves as an active veto against cosmic-ray muons, and provides passive shielding from radioactivity in the surrounding rock.The expected event rates at far detector are calculated by multiplying the unoscillated neutrino spectra (predicted by near detectors) with the corresponding oscillation probabilities.

IV. EVENT SELECTION AT FAR DETECTOR
The T2K Runs 1-8 data set used in this analysis was collected from Jan. 2010 to May 2017, corresponding to a beam exposure of 14.7 × 10 20 protons on target (POT) in neutrino mode and 7.6 × 10 20 POT in antineutrino mode.Events at the far detector are required to occur within 1 ms of the beam spill time window, and be fully contained in the fiducial volume of the SK ID.For the CC ν ( -) µ and ν ( -) e samples, a new Cherenkov-ring reconstruction algorithm [46] is used to select neutrino events, which improves signal/background discrimination and expands the fiducial volume to increase statistics.
There are five CC analysis samples that are commonly used in the standard three flavor oscillation analysis [47]: µ CC-0π and ν ( -) e CC-0π samples which are enriched in CC quasi-elastic (CCQE) events, and ν e CC-1π + sample where a π + below Cherenkov threshold is produced.The Bayesian analysis in [46]), and sin 2 θ 13 = 0.0219 (taken from [48]).All sterile mixing angles are set to zero.The event rates are summarized in Table I.
In addition to the CC samples, beam-induced NC events are also collected in SK.These events have previously only been used in publications for systematic uncertainties [47] and cross-section [49] studies.In this analysis, NCπ 0 and NC γ-deexcitation samples are used in the oscillation fit to enhance the sensitivity to sterile mixing parameters.
The NCπ 0 samples select neutrino events with single π 0 production, where π 0 → 2γ decay produces two visible Cherenkov rings in the detector.Events with two electron-like Cherenkov rings are selected as candidates for NCπ 0 samples, and those with decay electron candidates (from muons) are rejected.The invariant mass from the two rings is required to be between 85 MeV/c 2 and 135 MeV/c 2 to be consistent with the π 0 mass.From simulations, 68.5% (53.6%) of events originate from a ∆ resonance and 19.1% (34.9%) from coherent pion production for (anti)neutrino mode.The NC single pion resonant (NC1π) production is described by the Rein-Sehgal model [50], while the coherent production is described by a tuned model of Rein-Sehgal [51].In MC, NC events constitute 97.1% (98.5%) of the sample.
The NC γ-deexcitation sample was first reported in the measurement of neutrino-oxygen NC quasi-elastic (NCQE) cross-section [49].The NCQE cross-section is calculated by a spectral function model [52,53] with the BBBA05 form factor parameterization [54], reweighting as a function of neutrino energy to match the theoretical calculations [55].The NCQE interaction can knockout a nucleon, which produces primary γ rays from residual nucleus deexcitation and secondary γ rays when knocked-out nucleons interact with other nuclei in water.The emitted γ rays are 10 MeV per event on average, which is much less energetic than other samples.Momenta and vertex positions are reconstructed using the low-energy tools developed for SK solar neutrino analyses [56].Cuts on reconstructed energy, fiducial volume, event timing, vertex and reconstruction quality, detector pre-activity, and Cherenkov opening angle are applied sequentially to remove beam-unrelated (e.g.radioactivity) and beamrelated CC backgrounds.The sample is estimated to contain 76.9% NCQE and 17.6% NC non-QE events.The NC γ-deexcitation sample is currently available only for T2K Runs 1-4 (from Jan. 2010 to May 2013), corresponding to 6.56 × 10 20 POT.
Figure 3 shows the reconstructed π 0 momentum and γ energy distributions of the NCπ 0 and NC γ-deexcitation samples respectively, with event rates summarized in Table I.Since the event spectra have little information about true neutrino energy, the sensitivity to ∆m 2  41 is limited in the NC channel.

V. ANALYSIS STRATEGY
The overall analysis method is similar to that used in the standard T2K three flavor analysis [47], including the incorporation of off-axis near detector data and treatment of systematic uncertainties.However, this analysis performs a simultaneous fit to the five CC and three newly added NC samples to constrain the sterile mixing parameters in the "3+1" model.Systematic parameters are added to accommodate the possible uncertainties in the NC samples.
The neutrino flux is calculated by a data-driven simulation [57,58], which incorporates the conditions of the proton beam, magnetic horn current and neutrino beamaxis direction.Hadronic interactions are tuned with the thin target measurements in the NA61/SHINE experiment [59].Stability of the neutrino flux has been monitored by INGRID throughout the whole data taking period.At the peak energy 0.6 GeV, the (anti)neutrino mode beam contains 97.2% (96.2%) ν ( -) µ , with only 0.42% (0.46%) ν ( -) e contamination, and the flux uncertainty is approximately 9%.Neutrino events at the near and far detectors are generated by the NEUT 5.3.2neutrino interaction generator [60], which accounts for general interaction and crosssection effects.Most of the cross-section and neutrino flux parameters are constrained by ND280.The unoscillated CC candidate events at ND280 are classified into different samples according to the event topology, and is fit with a binned Poisson likelihood to extract the best-fit parameters and correlated uncertainties.The central values and their covariances are then propagated to estimate the far detector flux and cross-section parameters and uncertainty covariance matrix.The fit to ND data was done assuming no oscillation at ND280.This approximation is valid for small ∆m 2  41 , below around 0.3 eV 2 /c 4 .However, NC and ν e interaction parameters are not constrained by the ND280 fit.As a result, an additional uncorrelated 30% normalization uncertainty is used in this analysis for the NC1π and NCQE channels.The values of these uncertainties are conservative estimates determined from a previous cross-section analysis [61] and NCQE theoretical model comparisons [62,63].They therefore dominate the overall cross-section uncertainty in the NC oscillation samples.
At the far detector, there are NEUT parameters (present in ND280 fit but not constrained) that control the final state interactions inside nuclei and secondary interactions with water molecules, altering the event topologies of outgoing particles.For each event topology, the SK selection efficiency and mis-identification rate is separately parameterized.The NC γ-deexcitation sample has separate uncertainties related to primary and secondary γ production [49].The uncertainties in these parameters constitute the SK detector uncertainties.
The effects of systematic uncertainties on the predicted event rates are summarized in Table II.The NC crosssection errors are dominant in the NC samples.
Using the flux and cross-section inputs from ND280, the unoscillated event sample spectra at the far detector are calculated.Oscillation parameters are varied to obtain the best agreement between data and predicted event rates.A joint maximum-likelihood fit to eight far detector samples constrains the sterile mixing parameters sin 2 θ 24 , sin 2 θ 34 and ∆m 2  41 .The log-likelihood is defined as where n i is the number of events in the i-th data bin, and ) is the expected event rate with oscillation parameters θ and systematic parameters f .The last term in Eq. 7 accounts for the systematic penalty with ∆ f being the difference between the systematic parameters and their prior values, related by the covariance matrix V.The oscillation parameters sin 2 θ 23 , |∆m 2  32 | and δ CP are allowed to vary without constraint; θ 12 and ∆m 2  21 are fixed to their PDG values [48]; and a penalty term is used to constrain sin 2 2θ 13 = 0.0857 ± 0.0046.FIG.4: The T2K 90% exclusion limits on sin 2 θ 24 as a function of ∆m 2  41 , with results from other experiments [35,37,38].The areas on the right are excluded.
During the fitting process, at a grid point of (sin 2 θ 24 , sin 2 θ 34 , ∆m 2  41 ), the function in Eq. 7 is minimized with respect to the other oscillation parameters and systematic parameters.We use Wilks's theorem to estimate the confidence levels (C.L.) [64].The results are cross-checked with Gaussian CL s contours [65] to ensure no significant bias due to the physical limit of sin 2 θ 24 ≥ 0 and sin 2 θ 34 ≥ 0.

VI. RESULTS
We consider the parameter space of ∆m 2 41 > ∆m 2   21   which is most sensitive in T2K.Two categories of fits are done for normal (NH, ∆m 2 31 > 0) and inverted (IH, ∆m 2  31 < 0) neutrino mass hierarchies respectively.The case of ∆m 2  41 < 0 is very similar and can be obtained by flipping the hierarchy.The "3+1" best-fit differs from the standard three flavor best-fit by ∆χ 2 = 1.0 (4.7) for NH (IH).From 2500 sets of MC studies with statistical fluctuations, this level of disagreement is expected with the standard three flavor hypothesis in 50% (30%) of the studies.
In the (sin 2 θ 24 , ∆m 2 41 ) parameter plane, sin 2 θ 24 is scanned from 10 −3 to 1, and ∆m 2 41 from 10 −4 to 0.3 eV 2 /c 4 .For larger values of ∆m 2  41 , oscillations would be also seen at the near detectors, which is beyond the scope of this analysis.Figure 4 shows the T2K 90% exclusion limits together with results from other experiments.We have set the most stringent limit on sin 2 θ 24 for ∆m 2  41 < 3 × 10 −3 eV 2 /c 4 .The limit is weaker at larger ∆m 2  41 due to the lack of high energy events, resulting from the sharply peaked off-axis neutrino flux.The difference between NH and IH comes from the ∆m 2  43 oscillation term.It becomes particularly important when ∆m 2  41 ∼ ∆m 2  31 , as this results in very different values of ∆m 2  43 in NH and IH.In partly-degenerate cases where ∆m 2  41 and ∆m 2 31 are in integer multiples, the ν ( -) e and NC samples are important in resolving ambiguities.In cases where one of the ∆m 2 values is very small, matter effects can significantly alter which mass states are involved in the oscillation, but the overall ν ( -) e appearance probability is not affected by more than few percent so this does not significantly modify the exclusion limits.
The NC samples allow us to constrain θ 34 in conjunction with θ 24 .Because these samples have low statistics and large cross section uncertainties, we have limited sensitivity, but our results are consistent with other measurements.Figure 5 shows that we constrain sin 2 θ 24 < 0.1 and |U τ 4 | 2 = cos 2 θ 24 sin 2 θ 34 < 0.5 at 90% C.L. if ∆m 2 41 = 0.1 eV 2 /c 4 is assumed.At smaller ∆m 2 41 values, the limits are different between NH and IH.

VII. CONCLUSIONS
Data collected by the T2K experiment between 2010 and 2017 (T2K Runs 1-8) have been used to search for oscillation signatures due to light sterile neutrinos in the "3+1" model.The sterile mixing parameters (sin 2 θ 24 , sin 2 θ 34 , ∆m 2  41 ) are constrained by performing a joint fit of the five CC samples: ν ( -) µ CC-0π, ν ( -) e CC-0π, ν e CC-1π + , and the three new NC samples: ν ( -) NCπ 0 , NC γ-deexcitation, selected at the far detector.Systematic uncertainties on the neutrino flux and CC interaction cross-section are constrained by the ND280 data, while NC cross-section uncertainties are determined from a comparison of theoretical models and external data.The data are consistent with the standard three flavor oscillation hypothesis.Limits have been set on the sterile mixing parameters, with the world's best constraint on sin 2 θ 24 for 10 −4 eV 2 /c 4 < ∆m 2  41 < 3 × 10 −3 eV 2 /c 4 .The data related to the measurement and results presented in this paper can be found in [66].
Our current precision is restricted by statistics and the uncertainty on the NC interaction cross-section.Apart from future updates of the analysis as we take more data, dedicated systematic studies are required for further improvements to the precision.Another possible extension is to perform a joint analysis of near and far detector data that would expand the range of constraint to ∆m 2 41 1 eV 2 /c 4 with additional data at smaller L/E.

FIG. 3 :
FIG.3: Reconstructed π 0 momentum spectra of NCπ 0 samples (top and middle), and reconstructed gamma energy spectrum of NC γ-deexcitation sample (bottom).The numbers of events in data and MC prediction are shown in the legend.

TABLE I :
Number of events expected and observed in the eight oscillation samples used in this analysis.Three flavor oscillation is assumed in expected rate.

TABLE II :
Percentage systematic uncertainty on far detector event yields.