Final results on neutrino oscillation parameters from the OPERA experiment in the CNGS beam

The OPERA experiment has conclusively observed the appearance of tau neutrinos in the muon neutrino CNGS beam. Exploiting the OPERA detector capabilities, it was possible to isolate high purity samples of $\nu_{e}$, $\nu_{\mu}$ and $\nu_{\tau}$ charged current weak neutrino interactions, as well as neutral current weak interactions. In this Letter, the full dataset is used for the first time to test the three-flavor neutrino oscillation model and to derive constraints on the existence of a light sterile neutrino within the framework of the $3+1$ neutrino model. For the first time, tau and electron neutrino appearance channels are jointly used to test the sterile neutrino hypothesis. A significant fraction of the sterile neutrino parameter space allowed by LSND and MiniBooNE experiments is excluded at 90% C.L. In particular, the best-fit values obtained by MiniBooNE combining neutrino and antineutrino data are excluded at 3.3 $\sigma$ significance.

IIHE, Université Libre de Bruxelles, B-1050 Brussels, Belgium The OPERA experiment has conclusively observed the appearance of tau neutrinos in the muon neutrino CNGS beam. Exploiting the OPERA detector capabilities, it was possible to isolate high purity samples of νe, νµ and ντ charged current weak neutrino interactions, as well as neutral current weak interactions. In this Letter, the full dataset is used for the first time to test the three-flavor neutrino oscillation model and to derive constraints on the existence of a light sterile neutrino within the framework of the 3 + 1 neutrino model. For the first time, tau and electron neutrino appearance channels are jointly used to test the sterile neutrino hypothesis. A significant fraction of the sterile neutrino parameter space allowed by LSND and MiniBooNE experiments is excluded at 90% C.L. In particular, the best-fit values obtained by MiniBooNE combining neutrino and antineutrino data are excluded at 3.3 σ significance.

INTRODUCTION
The OPERA experiment [1] was designed to conclusively observe the appearance of tau neutrino in the high purity muon neutrino CNGS beam [2,3], in the parameter region indicated by Super-Kamiokande [4] and MACRO [5] to explain the zenith dependence of the atmospheric neutrino deficit. The average neutrino energy in the CNGS beam was designed to be about 17 GeV, theν µ fraction was 2.1% in terms of expected charged current (CC) interactions, the sum of the ν e andν e fractions was below 1%, while the prompt ν τ component was negligible.
The OPERA detector [1] was located in the underground Gran Sasso Laboratory (LNGS), about 730 km away from the neutrino source at CERN. The detector was a hybrid apparatus made of a nuclear emulsion/lead target complemented by electronic detectors (ED) [6]. The target had a total mass of about 1.25 kt and was composed of two identical sections. Each section consisted of 31 walls of Emulsion Cloud Chamber (ECC) bricks, interleaved by planes of horizontal and vertical scintillator strips (Target Tracker) used to select [7] ECC bricks in which a neutrino interaction had occurred. Each ECC brick consisted of 57 emulsion films interleaved with 1 mm thick lead plates, with a (12.7 × 10.2) cm 2 cross section and a total thickness corresponding to about ten radiation lengths. A magnetic muon spectrometer [1,8] was positioned downstream of each target section, and was instrumented by resistive plate chambers and highresolution drift tubes.
The experiment collected data from the CNGS beam from 2008 to 2012, with an integrated exposure of 17.97 × 10 19 protons on target. A total of 19505 neutrino interaction events in the target were recorded by the electronic detectors, of which 5603 were fully reconstructed in the OPERA emulsion films.
In 2015, the OPERA Collaboration reported the discovery of tau neutrino appearance [9]. Five ν τ candidate events [9][10][11][12][13] were observed, with an expected signal-tobackground ratio of ∼10. The analysis of these data lead to the exclusion of the background-only hypothesis with a significance of 5.1 σ.
The final result on ν τ appearance was published in 2018 [14], using an analysis based on a multivariate approach and a looser selection with respect to [9], applied to the complete data set. The achieved significance for tau appearance was 6.1 σ. Additionally, the atmospheric neutrino mass splitting was measured to be ∆m 2 32 = 2.7 +0.7 −0.6 × 10 −3 eV 2 , assuming maximal twoflavor mixing.
A search for electron neutrino CC interactions was performed [15,16]. Unlike ν τ appearance, the sensitivity to the ν µ → ν e oscillation channel is limited mainly by the ν e beam contamination and the systematic uncertainty in the neutrino flux. The collected sample of electron neutrino candidates was in agreement with the expectation from the three-flavor neutrino mixing model and was exploited to obtain a 90% C.L. upper limit on the parameter sin 2 2θ 13 < 0.43 [16].
The data were also used to test the hypothesis of one additional light sterile neutrino state, in the 3 + 1 neutrino mixing model [15][16][17]. OPERA is sensitive to a significant part of the allowed region of the sterile neutrino parameters defined by the LSND [18] and MiniBooNE [19] experiments.
In this Letter, the final samples of tau and electron neutrino CC interactions are used for the first time in a combined analysis to constrain standard oscillation parameters θ 13 , θ 23 , as well as the 3 + 1 mixing model ones. In addition, a search for the ν µ disappearance signal was performed using only electronic detector data, assuming three-flavor neutrino mixing.

DATA SELECTION
Selection of ν τ and ν e CC candidate events using emulsion detector data The emulsion detector data set consists of events in which the neutrino interaction vertex is fully reconstructed in the emulsion films. Tracking capabilities of the OPERA emulsion technique allow isolating samples of tau and electron neutrino CC interactions with a very high signal-to-background ratio [11,20]. In the final analysis, ten tau neutrino candidates have been identified [14], while the expected background amounts to 2.0±0.4 events. The sources of background are: (i) the decay of charmed particles produced in ν µ CC interac-tions, (ii) the reinteraction of hadrons from ν µ events in lead, and (iii) the large-angle scattering (LAS) of muons produced in ν µ CC interactions. Processes (i) and ν µ CC in (ii) represent a background source when the µ − at the primary vertex is not identified.
The search for electron neutrino appearance led to the selection of 35 ν e candidates in the full dataset resulting from a dedicated analysis [16]. The two main sources of background are: (i) ν τ CC interactions followed by a τ → e decay and (ii) muon-less events with a π 0 → γγ decay and an e + e − pair from prompt γ conversion misidentified as an electron. The total expected background amounts to 1.2 ± 0.5 events. Given the prompt ν e and ν e beam components, the number of observed ν e events is consistent both with no-oscillation (31.9 ± 3.3) and standard oscillation (34.3 ± 3.5) hypotheses.
Selection of ν µ CC and ν NC candidate events using only electronic detector data Only electronic detector data were used in a search for ν µ disappearance, since events without a decay topology (most of ν µ CC interactions and ν NC interactions) are not fully reconstructed in the emulsion target [20]. This analysis is complementary to the appearance study, as this dedicated selection of ED data differs from the one used to locate ECC bricks in which ν interactions occurred.
The compositions of CC-like and NC-like samples described below were estimated using standard Monte Carlo simulation validated on data [6], with GENIE v2.8.6 [21,22] as a neutrino interaction generator.
A global selection was applied to ED data: (i) the interaction vertex must be reconstructed inside the active volume defined in [23], (ii) a fiducial volume cut rejecting events in which neutrinos most probably interacted in the 10 brick walls immediately downstream of a spectrometer, and (iii) a cut on the total signal recorded by electronic detectors in order to reject soft non-beam events.
The ν µ CC interaction candidates (CC-like sample) are selected by requiring a muon track with negative charge reconstructed in at least two arms of a single magnetic spectrometer. According to the Monte Carlo simulation, the expected purity of the ν µ CC-like event sample is 99.5%.
The neutrino NC interaction candidates (NC-like sample) are selected by requiring no track identified as a muon. Assuming global best-fit values of the standard oscillation parameters for normal mass hierarchy (NH) [24], the MC expectation is that 75.3% of events are ν and ν NC interactions, 20.2% are ν µ and ν µ CC interactions in which no muon was identified, and 4.5% are CC interactions of other neutrino flavors. The latter are expected to be: (i) prompt ν e and ν e beam components, and (ii) flavors appearing through neutrino oscillations, most significantly in channels ν µ → ν τ and ν µ → ν e .
NC-like and CC-like data samples consist of 1724 and 5782 events, respectively.
It should be noted that the estimated composition of NC-like and CC-like samples quoted above do not explicitly enter the analysis, as the proper oscillation probability is applied to all MC events as described in the next section.

ANALYSIS
Appearance analyses Visible neutrino energy as defined in [14] and reconstructed neutrino energy as defined in [16] are used as observables for tau and electron neutrino samples, respectively. Their distributions are jointly exploited to test neutrino oscillation phenomena both in standard three-flavor and 3+1 mixing models.
The expected distributions are evaluated using GLoBES [25,26]. The smearing matrices, incorporating the detection and selection efficiencies, are calculated by means of MC simulation of the full analysis chain. The simulated CNGS flux [27], default neutrino interaction cross sections evaluated by GENIE v2.8.6 [21,22], and the smearing matrices are used as input. It is worth noting that the oscillations of all CNGS beam components (ν µ ,ν µ , ν e andν e ) are taken into account.
The systematic uncertainty on the expected number of tau neutrino candidates is largely dominated by the limited knowledge of the ν τ CC interaction cross section and detection efficiency, overall estimated to be 20% [14].
The expected number of electron neutrino candidates is affected by systematic uncertainties on: (i) the prompt flavor beam components, (ii) the ν e andν e CC interaction cross sections, and (iii) the neutrino detection efficiency. The systematic error on the expected number of electron neutrino candidates is estimated to be 20% for reconstructed energies below 10 GeV and 10% otherwise [16].
The statistical analysis of the data is based on a maximum-likelihood joint fit across the two samples. The data are binned and the overall likelihood (L) is constructed as the product of the contribution for each data sample (L τ , L e ) and a Gaussian constraint G on ∆m 2 31 as: where ∆m 2 31 and σ 31 are the best fit value and the 1 σ uncertainty from a global fit of the current neutrino oscillation data [24]. The used values for ∆m 2 31 ± σ 31 are (2.56 ± 0.05) × 10 −3 eV 2 and (−2.49 ± 0.05) × 10 −3 eV 2 for NH and inverted hierarchy (IH), respectively. The floating parameter µ 31 represents the true value of the squared mass difference.
The negative log-likelihood functions − ln L x are de-fined as: with x = τ or e, N (x) is the number of bins, n i is the number of events in the i -th bin of the data histogram, and µ i = µ i θ ,φ is the expected number of events in the i-th bin from the set of oscillations parametersθ, namely the mixing matrix elements U αβ and the squared mass splittings ∆m 2 ij , and the set of parametersφ which account for systematic uncertainties. The second term of Eq. 2 accounts for the a priori knowledge of the systematic uncertainties, whereφ j and σ φj are respectively the estimated value and the uncertainty on parameter φ j . The non-zero value of ∆m 2 21 is taken into account as well as matter effects assuming a constant Earth crust density estimated with the PREM [28,29] onion shell model, even though they were checked to be negligible for ∆m 2 41 > 1 eV 2 . Disappearance analysis A search for the ν µ disappearance signal is performed in the dedicated ED data sample. Since the experiment was not equipped with a near detector, the uncertainty on the ν µ flux is estimated to be 15 %. On the other hand, the detector was located far from the first oscillation maximum, so that the ν µ CC rate reduction at the detector due to oscillation does not exceed 2% in the standard 3-flavor model. To overcome these difficulties, a NC-like/CC-like event rate ratio is used as an estimator of the disappearance oscillation probability that is independent from the overall flux normalization uncertainty.
Since it is experimentally unfeasible to reconstruct the energy of an incoming neutrino in NC events, NClike/CC-like event rate is evaluated as a function of E T , a measurement of the energy deposited in the target.
In order to constrain the value of the atmospheric mass splitting ∆m 2 32 , three histograms with identical binning in E T are constructed: (i) the predicted NC-like/CClike ratio evaluated using MC simulation as a function of ∆m 2 32 , (ii) the observed number of CC-like events, (iii) the observed number of NC-like events.
While constructing histogram (i), global best-fit values [24] for all oscillation parameters except the parameter of interest ∆m 2 32 are used, separately for NH and IH. For each value of ∆m 2 32 , a proper oscillation probability is calculated for each simulated event in both NC-like and CC-like samples, including disappearance of ν µ , ν µ , ν e and ν e CC interactions and ν µ → ν τ and ν µ → ν e appearance.
A dedicated statistical model is built to properly construct confidence intervals based on a ratio of two Poisson random variables. Given two outcomes l and k of two Poisson distributions with mean values λ l and λ k , respectively, the Bayesian probability of ratio x ≡ λ l /λ k , assuming a flat prior on x, is provided by the formula The likelihood function is constructed as where j is a bin label running on all three histograms; NC j and CC j are the j-th bin content of the NC-like (ii) and CC-like (iii) histograms, respectively; R j (∆m 2 32 ) is the NC-like/CC-like ratio obtained from histogram (i) as a function of oscillation parameters; P NCj /CCj (x) is a probability distribution defined in Eq. 3.
In order to extract the confidence interval, the upper limit test statistic [30] is used, defined by equation where ∆m 2 32 is the value of ∆m 2 32 which maximizes the likelihood function L(∆m 2 32 ). A complete treatment of the statistical model can be found in [31].
The effect of uncorrelated uncertainties in the knowledge of the CNGS flux is estimated by dividing the simulated neutrino events into subsets (bins) according to true neutrino energy. The bins are 10 GeV wide, and the weight of every event in a single bin is multiplied by the same Gaussian random number with mean value of 1 and standard deviation of 0.15, corresponding to the uncertainty on the total neutrino flux. Event weights in different bins are multiplied by different random numbers from the same Gaussian distribution. A set of 1000 different smeared fluxes is generated in this way for each value of ∆m 2 32 . The likelihood function is modified to take into account these multiple smeared fluxes as follows: where M is the number of smeared fluxes, and R (j) i (∆m 2 32 ) is a MC predicted NC-like over CC-like event ratio in the i-th E T bin obtained using j-th smeared input flux.

RESULTS
Standard three-flavor model Appearance channels are jointly used to obtain constraints on three-flavor model parameters θ 23 and θ 13 . All other parameters are treated as nuisance and the profile likelihood ratio is used to remove the dependence of the likelihood (Eq. 1) on them. The 68% C.L. allowed region for NH, obtained using Wilks' theorem, is reported in Fig. 1. Resulting best fit value and 1 σ confidence interval for θ 23 parameter is θ 23 = 0.78 +0. 32 −0.31 rad, while the 1 σ confidence interval for θ 13 is [0, 0.20] rad. The ν µ → ν µ disappearance channel was explored using the test statistic (Eq. 5) based on the likelihood (Eq. 4) to obtain an upper limit on ∆m 2 32 . The values of ∆m 2 32 are sampled on a grid in the interval [0.0, 6.0] × 10 −3 eV 2 for NH, and [−6.0, 0.0] × 10 −3 eV 2 for IH. In both cases, the grid spacing is 6.0 × 10 −5 eV 2 . The distribution of q ∆m 2 32 (Eq. 5) is calculated for each of these grid points using 10 000 pseudo-experiments. This is then used to calculate a p-value for each point for both NH and IH cases, shown in Fig. 2 as "no smearing" curves. For both NH and IH assumptions, the 90 % confidence interval is given as ∆m 2 32 < 4.1 × 10 −3 eV 2 . The p-values obtained using smeared neutrino fluxes by using likelihood (Eq. 6) are shown in Fig. 2 as "with smearing" curves. There is no significant difference w.r.t p-values obtained without smearing.
Sterile 3+1 model The excess of electron neutrino and antineutrino events reported by the LSND [18] and MiniBooNE [19] experiments may be interpreted as due to the presence of light O(1 eV) sterile neutrinos. OPERA can test this hypothesis searching for deviation from the prediction of the three-flavor neutrino model. The oscillation probability in this model depends on a four by four unitary mixing matrix U and three squared mass differences (∆m 2 21 , ∆m 2 31 and ∆m 2 41 ) [32]. Here, for the first time, tau and electron neutrino appearance channels are jointly used to test the sterile neu-   FIG. 2. p-values as a function of ∆m 2 32 , with all other oscillation parameters fixed to central values [24]. Vertical width of the "NH with smearing" curve represents the 90 % C.L. interval on p-value due to the finite numbers of pseudoexperiments; other curves have similar uncertainty which is not shown for clarity. p-values for ∆m 2 32 < 2.5 × 10 −3 eV 2 are equal to 1 in all cases.
trino hypothesis in the framework of the 3 + 1 model. Since the NC-like over CC-like event ratio is poorly sensitive to the effects induced by the presence of a sterile neutrino state, the muon neutrino disappearance channel is not exploited in this analysis.
Defining sin 2 2θ µτ = 4|U τ 4 | 2 |U µ4 | 2 , sin 2 2θ µe = 4|U e4 | 2 |U µ4 | 2 and ∆m 2 41 as the parameters of interest, the profile likelihood ratio is used to remove the dependence of the likelihood (Eq. 1) on the other parameters treated as nuisance. Exclusion regions of ∆m 2 41 versus sin 2 2θ µτ and sin 2 2θ µe are shown in Fig. 3 and Fig. 4, respectively. The result is restricted to positive ∆m 2 41 values since negative values are disfavored by results on the sum of neutrino masses from cosmological surveys [33]. The MINOS+ Collaboration has recently reported an analysis setting limits on the existence of a sterile neutrino state [34].
For ∆m 2 41 > 0.1 eV 2 , the upper limits on sin 2 2θ µτ and sin 2 2θ µe are set to 0.10 and 0.019 both for the case of NH and IH. The limits obtained on the sterile oscillation parameters improve those reported in [16,17].
The values of the oscillation parameters ∆m 2 41 = 0.041 eV 2 , sin 2 2θ µe = 0.92 corresponding to the MiniBooNE combined neutrino and antineutrino best-fit [19] are excluded with a p-value of 8.9 × 10 −4 , corresponding to a significance of 3.3 σ.

CONCLUSIONS
Searches for ν τ and ν e appearance and for ν µ disappearance were performed using the full OPERA data sample.  [35] and CHORUS [36] are also shown.
The data are compatible with the three-flavor neutrino model and constraints on θ 23 and θ 13 were derived jointly for the first time exploiting tau and electron neutrino appearance channels. A best fit value of θ 23 = 0.78 +0. 32 −0.31 rad at 1 σ C.L. is obtained, while θ 13 is constrained to [0, 0.20] rad at 1 σ C.L.
Additionally, a dedicated sample of OPERA electronic detector data is used to perform a search for the ν µ disappearance signal in the CNGS beam. Assuming all other mixing parameters equal to the global fit central values [24], an upper limit ∆m 2 32 < 4.1 × 10 −3 eV 2 at 90% C.L. is obtained.
FIG. 4. The 90% C.L. exclusion region in the ∆m 2 41 and sin 2 2θµe plane is shown for the normal (NH, solid line) and inverted (IH, dashed line) hierarchy of the three standard neutrino masses. The plot also reports the 90% C.L. allowed region obtained by LSND [18] (cyan) and MiniBooNE combining ν andν mode [19] (yellow). The blue and red lines represent the 90% C.L. exclusion regions obtained in appearance mode by NOMAD [37] and KARMEN2 [38], respectively. The 90% C.L. exclusion region obtained in disappearance mode by the MINOS and DayaBay/Bugey-3 joint analysis [39] is shown as green line. The black star ( ) corresponds to the MiniBooNE best-fit values for the combined analysis of ν and ν data.
istry of Science and ICT (Grant No. NRF-2018R1A2B2007757) for Korea; and TUBITAK, the Scientific and Technological Research Council of Turkey for Turkey (Grant No. 108T324). We warmly acknowledge the important contributions given by our deceased colleagues G. Giacomelli, A. Ljubičić, G. Romano and P. Tolun.