Constraint on the solar Δm2 using 4000 days of short baseline reactor neutrino data

There is a well-known 2σ tension in the measurements of the solar Δm between KamLAND and SNO/ Super-KamioKANDE. Precise determination of the solar Δm is especially important in connection with current and future long baseline CP violation measurements. Seo and Parke [Phys. Rev. D 99, 033012 (2019)] points out that currently running short baseline reactor neutrino experiments, Daya Bay and RENO, can also constrain solar Δm value as demonstrated by a GLoBES simulation with a limited systematic uncertainty consideration. In this work, the publicly available data, from Daya Bay (1958 days) and RENO (2200 days) are used to constrain the solar Δm. Verification of our method through Δmee and sinθ13 measurements is discussed in Appendix A. Using this verified method, reasonable constraints on the solar Δm are obtained using above Daya Bay and RENO data, both individually and combined. We find that the combined data of Daya Bay and RENO set an upper limit on the solar Δm of 18 × 10−5 eV at the 95% C.L., including both systematic and statistical uncertainties. This constraint is slightly more than twice the KamLAND value. As this combined result is still statistics limited, even though driven by Daya Bay data, the constraint will improve with the additional running of this experiment.


I. INTRODUCTION
Evidence that neutrinos are massive and mix is well established by a significant number of experiments. In this paper, we are interested in the mass squared difference, Δm 2 21 ; the mass squared difference of the two mass eigenstates that have the greatest fraction of electron neutrino, ν 1 and ν 2 . This mass splitting is responsible for the neutrino flavor transformations that occur inside the Sun (hence the name the solar mass squared difference), and for the antineutrino oscillations observed at an L=E ∼ 15 km=MeV.
In this paper, we use publicly available data to follow up a recent paper [1], that Daya Bay [2] and RENO [3], the short baseline (∼1.5 km) reactor antineutrino experiments currently running, have enough data already collected to constrain Δm 2 21 . The combined constraint by Daya Bay and RENO, gives an important consistency check of the standard three neutrino paradigm as well as adding addition information to the size of Δm 2 21 . The ∼2σ tension between the combined Super-Kamiokande (SK) [4] & Sudbury Neutrino Observatory (SNO) [5] solar neutrino measurements and KamLAND [6] reactor experiment (L=E ∼ 50 km=MeV) is not directly addressed by this constraint. However such a combined Daya Bay plus RENO constraint is at a different L=E range (∼0.5 km=MeV) than the above mentioned measurements as well as JUNO [7]. Moreover, the ratio of Δm 2 21 to Δm 2 31 , at an L=E ∼ 0.5 km=MeV, is required for the precision measurement of leptonic CP violation parameter, by NOvA [8], T2K [9] and future Long Baseline (LBL) experiments.
Currently there are two measurements of the solar mass squared difference, Δm 2 21 . One measurement comes from a combined measurement by SNO and SK using the observation of a day-night asymmetry by SK and the nonobservation of the low energy up turn of the 8 B neutrino survival probability by SNO and SK. This combined result is from SNO and SK. Similar results are obtained by Nu-Fit [10]. The other measurement is from KamLAND, the long baseline reactor antineutrino experiment, see [6], at If CPT invariance is a good symmetry of nature then the Δm 2 21 measured from solar neutrinos and reactor antineutrinos is required to give the same value. Currently this important parameter for neutrino physics suffers from a 2σ level tension. This tension could come from new physics, some error in the analysis of one or more of the experiments or a statistical fluctuation.
In the bievent plane for T2K, see Fig 44 of [14], Nðν μ → ν e Þ ¼ 37 and Nðν μ →ν e Þ ¼ 4 is outside the allowed region (by about 1σ). This can be well accommodated by a Δm 2 21 value, approximately twice the KamLAND value. Again, this is probably a statistical fluctuation but with only the KamLAND precision measurement of Δm 2 21 , other possibilities are still viable. The future medium baseline, L=E ∼ 15 km=MeV, reactor experiment JUNO will measure to better than 1% precision Δm 2 21 and sin 2 θ 12 , see [7]. JUNO experiment is currently under construction and their precision measurements of Δm 2 21 and sin 2 θ 12 will not be available until approximately 5 years from now. Later next decade, the proposed experiments Hyper-K & DUNE will also give us precision measurements of Δm 2 21 using 8 B solar neutrinos, see [15,16] respectively.
In Sec. II, we briefly discuss in detail the effects of increasing Δm 2 21 on theν e survival probability. Then in Sec. III Daya Bay and RENO data sets used in this work are discussed followed by Secs. IV, V, and VI for methods and systematic uncertainties, results, and conclusion, respectively.

II. SURVIVAL PROBABILITY
In vacuum, the electron antineutrino survival probability is Pðν e →ν e Þ ¼ 1 − P 12 − P 13 with P 12 ¼ sin 2 2θ 12 cos 4 θ 13 sin 2 Δ 21 ; where the kinematic phases are given by Δ jk ≡ Δm 2 jk L=ð4EÞ and θ 13 ≈ 8°and θ 12 ≈ 33°are the reactor and solar mixing angles respectively. The P 12 term is associated with the solar oscillation scale of 15 km=MeV and the P 13 term is associated with the atmospheric oscillation scale of 0.5 km=MeV. To excellent fractional precision, 2 the P 13 term can be approximated by where Δm 2 ee ≡ cos 2 θ 12 Δm 2 31 þ sin 2 θ 12 Δm 2 32 [17,18], interpreted as the ν e average of Δm 2 31 and Δm 2 32 . Using the fit values given in [10], and an L=E range around the first oscillation minimum (L=E ∼ 0.5 km=MeV), P 12 and P 13 is well approximated by: The P 12 term is almost negligible for all L=E < 1 km=MeV, if Δm 2 21 ¼ 7.5 × 10 −5 eV 2 . For Daya Bay and RENO this covers the full L=E range.
Suppose that Δm 2 21 is 3 times larger than KamLAND value, i.e., 22.5 × 10 −5 eV 2 , then 1 In the rest of this paper, when referring to neutrinos, we mean both neutrinos and/or antineutrinos. 2 The fractional precision is better than 0.05% for L=E < 1 km=MeV. Also, in this L=E range, the exact P 13 is very insensitive to mass ordering provided the value of jΔm 2 ee j is the same for both mass orderings.
Now P 12 is now no longer tiny compared to P 13 at L=E ¼ 0.5 km=MeV, oscillation minimum, and as L=E gets larger than 0.5 km=MeV, P 12 gets bigger, whereas P 13 is getting smaller. At an L=E ¼ 1 km=MeV, P 12 would be approximately equal to sin 2 2θ 13 (0.08) for this value of Δm 2 21 . It is this quadratic rise in P 12 as Δm 2 21 increases that we exploit to place an upper limit on Δm 2 21 . For further details on the survival probability as Δm 2 21 increases see [1].

III. DAYA BAY AND RENO DATA SETS
In this work, 1958 days of Daya Bay data [19] and 2200 days of RENO data [20] are used, where Daya Bay has about five times more inverse beta decay (IBD) events than RENO in their far detectors. Daya Bay data including background estimation, energy response function, and systematic uncertainties are taken from the supplementary material in [19]. RENO data and background estimation are extracted from Fig. 1 in [20] and systematic uncertainties are also taken from [20]. Table I shows summary of the basic parameters, i.e., L eff , IBD rate, and background rate, for near and far detectors of Daya Bay and RENO used in this analysis. Note that there are two near detectors in different sites for Daya Bay.

IV. METHODS AND SYSTEMATIC UNCERTAINTIES
Best fit values on Δm 2 21 and sin 2 2θ 13 are obtained by finding minimum χ 2 values between data and predictions for all possible combination of the two parameters. Far-tonear ratio method is employed in this χ 2 analysis to avoid the spectral shape anomaly around 5 MeV region [21] as well as to reduce systematic uncertainties.
The χ 2 formalism as written below contains a covariance matrix (V stat;ij ) to include statistical uncertainty and pull parameters (ξ α ) to include systematic uncertainties.
For both Daya Bay and RENO, systematic uncertainties on the relative detection efficiency, relative energy scale and the main background contributions are taken into account as summarized in Table II. Besides the systematic uncertainties, additional systematic paddings (adjustment factors) are added in our work to match Daya Bay and RENO results on θ 13 and Δm 2 ee measurements. For Daya Bay a 1.3 adjustment factor to the relative energy scale and Li-He background uncertainties is added. Whereas in RENO a 1.4 adjustment factor is added to the relative detection efficiency uncertainty. More details on the validation of our method and expected event description can be found in Appendixes A and B. The RENO predictions are computed using the Daya Bay detector response function and the relative far-to-near normalization is computed comparing our total number of expected events with the total number of expected events in the RENO Far detector. In order to match the best fit values of θ 13 and Δm 2 ee a 0.984 adjustment factor is added to this normalization of a total event rate for RENO.

V. RESULTS
A 2-dimensional scan over Δm 2 21 and sin 2 2θ 13 is performed to find the best fit value pair at the minimum value of χ 2 described earlier, where in the oscillation probability, the parameter θ 12 is fixed 3 at sin 2 θ 12 ¼0.310. The Δm 2 ee parameter is constrained with a pull parameter, allowing it to vary within a 2σ range of a prior Δm 2 ee value with a penalizing term Δm 2 ee;prior − Δm 2 ee σ 2 The prior Δm 2 ee value and its uncertainty are taken to be which is inferred from the combined measurement on Δm 2 μμ by current long baseline neutrino experiments in [10] through Δm 2 ee ≃ Δm 2 μμ AE cos 2θ 12 Δm 2 21 , see [17], where the þ=− comes from the unknown mass ordering (NO=IO) and ignoring terms proportional to sin θ 13 Δm 2 21 . The unknown mass ordering is treated as an additional uncertainty (4%) to Δm 2 μμ uncertainty (4%) for the Δm 2 ee uncertainty which, therefore, becomes about 6%. The best fit, 1, 2, and 3σ allowed regions of Δm 2 21 vs sin 2 2θ 13 are shown in Fig. 1 with (solid lines) and without (dashed lines) systematic uncertainties for Daya Bay and RENO, separately and combined. Daya Bay's result is better than RENO's due to about five time more statistics at the far detector, see Table I. In Fig. 3, we give the constraints on the three parameter fit, Δm 2 21 , Δm 2 ee and sin 2 2θ 13 , without imposing any constrain on Δm 2 ee , using the combined Daya Bay and RENO data sets. Both statistical and systematic uncertainties are included in this plot. As before θ 12 is fixed at sin 2 θ 12 ¼ 0.310, see [1] for discussion on allowing sin 2θ 12 to also vary. Results with Δm 2 ee fixed or free are obtained for each experiment and for when the data from both experiments are combined. These are described and given in Appendix C. It was found that the effect of free Δm 2 ee is bigger than that of systematic uncertainty, but our representing results are based on constrained Δm 2 ee since it is a reasonably well measured oscillation parameter using LBL experiments.

VI. CONCLUSION
Using the currently available public data from Daya Bay (1,958 days) and RENO (2,200 days), we have provided additional information on the solar Δm 2 . A reasonable upper bound is obtained from a combined analysis of the Daya Bay and RENO data as 18 × 10 −5 eV 2 at 95% CL, where Δm 2 ee was constrained using a pull parameter with input information from LBL experiments. Our combined analysis result is currently limited by statistics and, as expected, Daya Bay data drives the combined analysis results. Our analysis method was validated by reproducing the Δm 2 ee and sin 2 θ 13 contours for each experiment as discussed in Appendix A.
Given that the previous measurements by KamLAND and SK/SNO of the solar Δm 2 are in a 2σ tension and the importance of solar Δm 2 for the determination of CP violation in LBL experiments, it is crucial that we understand the value of the solar Δm 2 better. It is expected by circa 2025 that the JUNO experiment will provide additional, important information on the value of the of solar The agreement between our results and Daya Bay as well as RENO for the measurements of Δm 2 ee vs sin 2 2θ 13 is an excellent validation of the methods and numbers used in our analysis. Therefore, our constraint on Δm 2 21 , using the publicly available data of Daya Bay and RENO, has a firm basis.

APPENDIX B: NUMBER OF EXPECTED EVENTS AND PULL PARAMETERS IN χ 2
The expected numbers of signal events in a detector d in a prompt energy bin i, X d i , is computed as follows up to a common input (e.g., reactor power, total number of protons) which cancels when taking ratios in the χ 2 computation.
where, the indices i, r, d, and iso refers to the ith energy bin, rth reactor, dth detector, and a fissionable isotope ( 235 U, 239 Pu, 238 U, or 241 Pu), respectively, and a d is the detector efficiency. L rd is the baseline between the reactor r and the detector d. E ν and E rec are the neutrino true energy and the reconstructed energy, both related by the detector response function RðE rec ; E ν Þ. The σðE ν Þ is the IBD cross section computed performing the integral in d cos θ of the differential cross section in [22] and the f iso is the averaged fission fraction 4 and the ϕ iso ðE ν Þ is the Huber-Mueller flux prediction [23,24]. P rd ν e →ν e ðE ν Þ is the oscillation probability from reactor r to detector d in the three neutrino oscillation paradigm.
The pull parameters accounting for detection efficiency (ϵ d ) and relative energy scale (η d ) are included in the number of expected events as follows For RENO, the efficiency pull parameter is included in the ratio. The background pull parameters are included in back- acc;i and B d n;i ) represents the number of total (Li-He, accidental and fast neutron) background events in the ith prompt energy bin in the dth detector, and the small b represents the corresponding pull parameter.

APPENDIX C: FIXED VS FREE Δm 2 ee
For the results in the main body of our paper we constrained Δm 2 ee treating it as a pull parameter using LBL experiments input. In this section we show the impact of Δm 2 ee fixed and set free. A 2-dimensional scan over Δm 2 21 and sin 2 2θ 13 is performed to find the best fit value pair at the minimum value of χ 2 described earlier, where in the oscillation probability θ 12 is fixed as sin 2 θ 12 ¼ 0.310 but Δm 2 ee is set free within the range of ½1.55; 3.55 × 10 −3 eV 2 . Results with a fixed Δm 2 ee ¼ 2.45 × 10 −3 eV 2 are also obtained and compared to those with Δm 2 ee set free. Figure 6, left and middle panels, shows the results of Δm 2 ee fixed and free for Daya Bay and RENO. It is observed that the effect of floating Δm 2 ee is bigger than adding systematic uncertainty for both Daya Bay and RENO. For floating Δm 2 ee case, the corresponding FIG. 5. Our validation on Δm 2 ee vs sin 2 2θ 13 fit using the RENO data (2200 days), including systematics and statistics uncertainties in red solid lines, and including statistics only in blue dashed lines, for 1, 2 and 3σ allowed regions. The fit of the RENO collaboration with 2,200 days from [20] is represented in the solid black lines. The agreement between our analysis (solid red lines) and RENO's analysis (solid black lines) is excellent.

Δm 2
ee values for the minimum χ 2 are found to be 2.50 × 10 −3 eV 2 (2.68 × 10 −3 eV 2 ) for Daya Bay (RENO) and it is within 1σ uncertainty of each of their measurements. Figure 6, right panels shows the results with combined analysis. For floating Δm 2 ee case, the corresponding Δm 2 ee value for the minimum χ 2 is found to be 2.54 × 10 −3 eV 2 and it is within 1σ uncertainty of the Daya Bay best fit value, i.e., ½2.52 AE 0.07 × 10 −3 .