First Measurement of Missing Energy Due to Nuclear Effects in Monoenergetic Neutrino Charged Current Interactions

Abstract

We present the first measurement of the missing energy due to nuclear effects in monoenergetic, muon neutrino charged-current interactions on carbon, originating from $K^{+}\rightarrow\mu^{+}\nu_{\mu}$ decay-at-rest ($E_{\nu_{\mu}}=235.5$ MeV), performed with the JSNS² liquid scintillator based experiment. Towards characterizing the neutrino interaction, ostensibly $\nu_{\mu}n\rightarrow\mu^{-}p$ or $\nu_{\mu}$${}^{12}\mathrm{C}$ $\rightarrow\mu^{-}$${}^{12}\mathrm{N}$, and in analogy to similar electron scattering based measurements, we define the missing energy as the energy transferred to the nucleus ($\omega$) minus the kinetic energy of the outgoing proton(s), $E_{m}\equiv\omega-\sum T_{p}$, and relate this to visible energy in the detector, $E_{m}=E_{\nu_{\mu}}~{}(235.5~{}\mathrm{MeV})-m_{\mu}~{}(105.7~{}\mathrm{MeV})-E_{vis}$. The missing energy, which is naively expected to be zero in the absence of nuclear effects (e.g. nucleon separation energy, Fermi momenta, and final-state interactions), is uniquely sensitive to many aspects of the interaction, and has previously been inaccessible with neutrinos. The shape-only, differential cross section measurement reported, based on a $(77\pm 3)$% pure double-coincidence KDAR signal (621 total events), provides an important benchmark for models and event generators at 100s-of-MeV neutrino energies, characterized by the difficult-to-model transition region between neutrino-nucleus and neutrino-nucleon scattering, and relevant for applications in nuclear physics, neutrino oscillation measurements, and Type-II supernova studies.

I Introduction

Models of neutrino-nucleus interactions are a crucial part of all accelerator-based neutrino physics programs performing oscillation measurements. Such models, usually implemented via a neutrino event generator, are used to form multiple aspects of each measurement, including expected signal and background rates, observable event topologies, and detector efficiency, bias and resolution. As we enter an era of precision measurements, uncertainties in neutrino interaction physics need to be controlled at the percent level or better for next generation experiments to achieve their physics goals [1]. In recognition of this challenge, there is an ongoing global program of dedicated neutrino interaction measurements, across many energies and with an array of nuclear targets [2]. Monoenergetic (235.5 MeV) muon neutrinos from charged kaon decay-at-rest (KDAR; $K^{+}\rightarrow\mu^{+}\nu_{\mu}$, with a branching ratio of 63.6% [3]) represent a unique and important part of this program, in particular for informing interaction models in the 100s-of-MeV region [4]. At these energies, the utilization of such models is wide-ranging, with relevance across multiple current and future neutrino oscillation experiments [5, 6, 7, 8, 9, 10, 11] and understanding Type-II supernovae and the neutrino-induced nucleosynthesis processes inside [12].

KDAR neutrinos are an important tool for studying neutrino interactions and the nuclear response because their energy is known, unlike all other relevant sources of neutrinos above the $\nu_{\mu}$ charged current (CC) threshold. These neutrinos are also particularly interesting because their characteristic energy transfer range ($0<\omega<120$ MeV) lies in the difficult-to-model transition region between neutrino-on-nucleus and neutrino-on-nucleon scattering, in which the interaction evolves from inducing collective nuclear excitations among multiple nucleons to quasi-elastic scattering off of individual nucleons. In addition, along with neutrino interaction and nuclear physics, KDAR neutrinos have been proposed as a signature of dark matter annihilation in the sun [13, 14] and as a source for neutrino oscillation searches at short- and long-baseline [15, 16, 17, 18].

Despite this widespread applicability, there exists only one measurement of KDAR neutrinos to date, based on a 3.9$\sigma$ observation of the process and coarsely presented in terms of shape-only muon kinetic energy $\frac{1}{\sigma}\frac{d\sigma}{dT_{\mu}}$ [19]. However, several detailed cross section calculations have been performed for KDAR neutrino interactions on carbon and argon [20, 21, 12, 22] and a number of neutrino event generators can be used to simulate these events. Unfortunately, and as exemplified in Ref. [12], the predictions vary significantly; depending on the generator or model, the total KDAR $\nu_{\mu}$-carbon cross sections range from $0.9-1.8\times 10^{-39}~{}\mathrm{cm}^{2}/\mathrm{neutron}$ and agreement among the anticipated kinematic distributions [e.g. muon energy ($E_{\mu}$) and angle ($\theta_{\mu}$)] is similarly weak. In general, the large disagreement amongst predictions underscores the importance of measurements of this process for both fully elucidating the neutrino-nucleus interaction and precisely studying oscillations involving 100s-of-MeV neutrinos.

The J-PARC Sterile Neutrino Search at the J-PARC Spallation Neutron Source (JSNS²) experiment employs a 48 ton liquid scintillator detector (96 inner PMTs and 24 veto PMTs) at 24 m from the Materials and Life Science spallation neutron source originating from 3 GeV protons focussed onto a mercury target, currently operating at 950 kW. Proton interactions with the target, surrounded mainly by concrete and iron shielding, readily produce pions, muons, and kaons to create the primarily decay-at-rest neutrino source – about 97.8% of neutrino parent $K^{+}$ come to rest before decaying. In particular, JSNS² uses the muon decay-at-rest component ($\mu^{+}\rightarrow e^{+}\overline{\nu}_{\mu}\nu_{e}$) to search for possible signatures of short-baseline oscillations near $\Delta m^{2}\sim 1$ eV² ($\overline{\nu}_{\mu}\rightarrow\overline{\nu}_{e}$; $\overline{\nu}_{e}p\rightarrow e^{+}n$). Most relevant for this paper, the protons produce KDAR neutrinos at the rate of 0.0038 KDAR $\nu_{\mu}$ per proton on target (POT), according to a detailed, high statistics, Geant4 [23] simulation of the target and shielding geometry. Notably, however, kaon production at the source is highly uncertain, with the MARS simulation package [24] predicting nearly a factor of two higher rate.

The JSNS² experiment [25, 26] began taking data in July 2020 and has gathered $4.9\times 10^{22}$ POT as of 6/2024; The detector is described in detail in Ref. [27]. JSNS² is sensitive to the scintillation light produced in neutrino interactions inside of the 17 ton cylindrical target volume, composed of linear alkylbenzene, 3 g/L PPO, 15 mg/L bis-MSB, and loaded at 0.1% concentration with gadolinium. The surrounding 31 ton buffer volume is composed of a similar liquid scintillator, but without gadolinium and with 30 mg/L bis-MSB; the outer buffer volume is instrumented to veto events originating outside the target volume. For the run period analyzed here, an additional solvent di-isopropylnaphthalene (DIN) was added to the target volume at the 8% level to enhance pulse shape discrimination ability, for differentiating fast-neutron-like and electron-like signals.

For KDAR $\nu_{\mu}$ CC events, we associate the reconstructed visible energy, after accounting for quenching, light propagation, and other detector effects, to the interaction products via $E_{vis}=T_{\mu}+\sum T_{p}$, where $T_{\mu}$ and $T_{p}$ represent the kinetic energy associated with the produced muon and proton(s) respectively, noting that multiple protons can be produced via final-state interactions (FSI). Small, ones-of-MeV-scale corrections due to nuclear de-excitation gammas, interactions of possible FSI-induced fast neutrons, and the nuclear recoil energy render this definition of “visible energy” not quite literal, but we account for this in the detector unfolding procedure described below. Naively, assuming that the incoming KDAR neutrino interacts quasi-elastically with a neutron to produce a muon and a proton, $\nu_{\mu}n\rightarrow\mu^{-}p$, one might assume that the calorimetric energy is simply $E_{vis,naive}=E_{\nu_{\mu}}(235.5~{}\mathrm{MeV})-m_{\mu}~{}(105.7~{}\mathrm{MeV})+[m_{n}-m_{p}]~{}(1.3~{}\mathrm{MeV})=131.2~{}\mathrm{MeV}$. However, nuclear physics effects, in particular nucleon separation energy and Fermi momenta, as well as Pauli blocking and FSI, among others, significantly modify this expectation. As shown below, predictions for this modification among generators and models vary considerably.

In this article, we report the first measurement of this missing energy due to nuclear effects using neutrinos, quantified as $E_{m}\equiv\omega-\sum T_{p}=E_{\nu_{\mu}}(235.5~{}\mathrm{MeV})-E_{\mu}-\sum T_{p}=E_{\nu_{\mu}}(235.5~{}\mathrm{MeV})-m_{\mu}~{}(105.7~{}\mathrm{MeV})-E_{vis}$. This is analogous to the missing energy variable commonly used to describe quasi-elastic events in electron scattering, [e.g. Ref. [28] with ¹²C$(e,e^{\prime}p)$]. Note that the neutron-proton mass difference is ignored in this definition for simplicity.

II Reconstruction, Simulation, and Analysis

This analysis is based on $1.3\times 10^{22}$ POT taken from January-June 2021. During this data taking period, the beam and detector performed reliably and were stable as determined by a variety of metrics, including light yield based on muon-decay Michel electrons, regular in-situ LED-based calibrations, neutron-capture on Gd signals, beam timing, singles flash rate, neutrino candidates per proton on target, and slow control monitoring. Most notably, we find that the effective light yield varied by as much as 2.6% over the course of the run, which we correct for and assign a systematic uncertainty to.

Our analysis seeks to identify and reconstruct KDAR $\nu_{\mu}$ CC events, characterized by a “prompt” scintillation light flash originating from the muon and proton(s) ($\sim$20-150 MeV) followed by a “delayed” Michel electron signal (0-53 MeV; $\tau_{\mu}=2.0$ $\mu$s [29]), to form a differential cross section measurement in terms of $E_{m}$:

$$(\frac{d\sigma}{dE_{m}})_{i}\sim\frac{\sum_{j}M_{ij}(d_{j}-b_{j})}{\epsilon_{i}}~{},$$ (1)

where the indices $i$ and $j$ respectively represent true and reconstructed $E_{m}$ (based on the prompt flash detector observable, $E_{prompt}$), $d_{j}$ is the measured data distribution passing selection, $b_{j}$ is the predicted background passing selection, $M_{ij}$ is the unfolding matrix to transform from reconstructed to true $E_{m}$, and $\epsilon_{i}$ is the reconstruction efficiency. Note that we report a shape-only differential cross section, $\frac{1}{\sigma}(\frac{d\sigma}{dE_{m}})_{i}$; The choice to ignore normalization is based on the highly uncertain kaon production at the source previously mentioned.

Towards forming this measurement, we utilize a PMT-charge-based maximum-likelihood reconstruction technique to relate the measured PMT pulses to the position and energy of the event. This reconstruction procedure relies on likelihood functions representing the probabilities to measure detector observables as a function of the underlying event quantities and determined using a detailed detector simulation [30, 31]. The simulation is largely based on Geant4 with light propagation model inputs, including PMT geometry and response, Birks’ constant for electrons/muons, attenuation length, and others, from the Daya Bay and RENO experiments [32, 33] and a number of in-situ calibrations, including those mentioned above and with a ²⁵²Cf source deployed across a range of vertical positions in the detector [31].

In particular, cosmic-induced Michel electrons provide an excellent, high statistics calibration source throughout the detector. We compare Monte Carlo (MC) simulated Michel electron events to detected events by fitting each with an analytic function for the Michel spectrum, convolved with an energy resolution function, with fit parameters for resolution and energy scale. The relative energy resolution is modelled as $\frac{\sigma{(E)}}{E}=\sqrt{\sigma_{ep}^{2}\frac{E_{ep}}{E}+C^{2}}$, where $\sigma_{ep}$ is the resolution at the endpoint energy ($E_{ep}$) and is a free parameter in the fit and $C$ represents the constant minimum possible energy resolution, estimated to be $\sim$2.5%. Notably, while the observed endpoint energy resolution is (3.79$\pm$0.08)%, the MC simulated resolution is determined to be (2.78$\pm$0.06)%; we apply an additional energy smearing to MC events to account for this discrepancy. We also compare the energy scale as a function of position throughout the inner volume between MC simulated Michel electron events and observed events, and correct the MC events’ reconstructed energy to compensate for the observed discrepancy between the two. Figure 1 shows the energy scale discrepancy from the nominal value throughout the detector.

For each correction to the energy scale and/or resolution we estimate its uncertainty, and account for these as systematics on the energy reconstruction. We furthermore account for a potential energy scale non-linearity by comparing the scale observed for neutron capture on gadolinium events ($E\approx 8$ MeV) and for Michel electron events ($E_{ep}\approx 53$ MeV). The observed energy scale non-linearity matches the predicted non-linearity from simulation to within 0.2%. We extrapolate this difference to the KDAR energy region and apply it as an additional 0.45% systematic uncertainty on the energy scale. At energies above the Michel electron endpoint we expect electronics saturation and non-linearity in our PMTs’ charge response to become increasingly prevalent. These effects are modeled within the MC simulation and reconstruction software to mitigate their contribution to the energy scale non-linearity. The final energy scale and energy resolution systematic uncertainty (at the Michel endpoint energy) is 0.68% and 0.31% respectively. Figure 1 also shows the observed Michel electron energy spectrum and the corresponding simulation with systematic uncertainties applied.

Figure 1: (Top) The observed Michel electron energy scale as a function of position relative to the middle-most position bin; the value is indicated both by the histogram color and bin text. The red box indicates the fiducial volume used for energy calibration. (Bottom) The observed (black) and MC simulated (blue) Michel electron energy spectrum across the fiducial volume. The blue error band corresponds to the energy reconstruction systematic uncertainties on the energy scale and resolution applied to the simulated events.

For simulating neutrino events, we utilize the NuWro [34, 35] and GiBUU [36, 37] event generators combined with the detector simulation to find the detection efficiency ($\epsilon_{i}$) and unfolding matrix ($M_{ij}$) described above. The Michel electron calibration results and systematic uncertainties are used to constrain the MC simulation energy response.

Figure 2: Measured distributions of signal-like (black) and background-like (red) events, based on whether or not they pass selection criteria other than timing. The two pulses corresponding to the selection’s proton-on-target arrival windows are shown. The background distribution shows the prompt, signal-like ($20<E_{prompt}<150$ MeV) events that are rejected by the non-timing selection criteria–most of these are due to beam fast neutrons.

Signal-induced prompt events are searched for in a timing window corresponding to the beam-on-target arrival plus time of flight. The beam strikes the target with two adjacent pulses, each reasonably consistent with a Landau-Gaussian convolution ($1\sigma\approx 40$ ns), separated by 600 ns, at 25 Hz. We open two corresponding 150 ns windows around the expected KDAR arrival times for selecting events, also in consideration of the $K^{+}$ lifetime ($\tau_{K^{+}}=12.4$ ns). Figure 2 shows the measured beam timing distributions for signal prompt events passing all KDAR selection criteria other than timing, further discussed below, and background events that are rejected by other selection criteria. The delayed signal due to the muon-decay Michel electron is searched for in a 10 $\mu$s timing window after the prompt beam pulse. In addition to timing selection of the two signal flashes, we also require candidate events to have (1) reconstructed energies $20<E_{prompt}<150$ MeV and $20<E_{delayed}<60$ MeV, with the lower limits set by background mitigation and signal efficiency considerations and the upper limits set by kinematic thresholds; (2) a primary reconstructed “vertex”, or more accurately “mean gamma emission point” ($MGEP$), in a fiducial volume, described with cylindrical coordinates relative to the center of the detector, of $R<1.4$ m and $-1.0<z<0.5$ m, where the boundaries of the target volume are $0<R<1.6$ m and $-1.25<z<1.25$ m; (3) a reconstructed distance between the primary flash and delayed flash $MGEP$ ($\Delta_{MGEP}$) of $<30$ cm; and (4) not be associated with veto activity, noting that a 10 $\mu$s deadtime is imposed after a detected cosmic muon [38].

The selection criteria are imposed to mitigate background contributions from accidental coincidences, cosmics, beam-induced fast-neutrons, and beam-induced non-KDAR neutrino and antineutrino events. Background events from accidental coincidences, primarily from uncorrelated cosmic and/or beam-related activity, are characterized by purposefully producing prompt-delayed pairs with events from different beam spills, treating them as though they were from the same beam spill. True correlated signals do not appear across beam signals, so any apparent correlation produced with this scheme will be accidental. The mis-paired events are subjected to the KDAR selection criteria to estimate the rate and spectral shape of the accidental background. By this method the accidental background rate is estimated to be $<$1.14 events in the KDAR dataset. Backgrounds from true correlated events, produced by cosmogenic activity (or any other ambient source), are characterized using data taken when the beam was off. Applying the KDAR event selection to these data produces the estimated rate and spectral shape for the cosmogenic background; The cosmogenic background rate within the analyzed dataset is estimated to be $<$6.57 events.

Beam produced fast-neutrons are a potential background source since they can interact within the target volume to produce a pion and subsequent decay muon and Michel electron, thus potentially faking the prompt and delayed signature of a KDAR event. However, beam neutrons enter the detector relatively late compared to KDAR events, because most are non-relativistic, so we can constrain the fast-neutron rate by looking for an excess of events with late timing; we find that this background is negligible and consistent with zero.

The most significant background is from beam non-KDAR $\nu_{\mu}$ and $\overline{\nu}_{\mu}$ that are sufficiently energetic to undergo the same CC interaction as the KDAR neutrinos, the dominant source of which is from postively charged pion and kaon decay-in-flight ($\pi^{+}\rightarrow\mu^{+}\nu_{\mu}$, $K^{+}\rightarrow\mu^{+}\nu_{\mu}$). Notably, the timing, $E_{prompt}$, $MGEP$, and $\Delta_{MGEP}$ characteristics of the observed background in a kinematically disallowed KDAR sideband region ($\Delta_{MGEP}>30$ cm and $E_{prompt}>150$ MeV) are highly consistent with non-KDAR $\nu_{\mu}$/$\overline{\nu}_{\mu}$ CC events. We estimate this background using a neutrino flux provided by a Geant4 simulation of the beam-target interaction and the surrounding target geometry. Uncertainty on the neutrino flux is estimated by performing the simulation using both the QGSP-BERT [39] and the QGSP-BIC [40] physics lists within Geant4 – the difference between the two simulation results is used as a 1$\sigma$ flux shape uncertainty. The background neutrino interactions are simulated using both the NuWro and GiBUU event generators, once again taking the difference between their predictions as a systematic uncertainty.

Figure 3: After selection, the reconstructed KDAR prompt energy ($E_{prompt}$) spectrum including the remaining background components. The inset shows the reconstructed energy spectrum of delayed Michel signals among KDAR-candidates in black, and the red shows the scaled MC simulated KDAR delayed spectrum.

The final estimates for the signal and background rates come from a simultaneous fit to the KDAR visible energy spectrum, including sideband constraints, and further forcing the signal rate to zero in the kinematically disallowed region, $E_{prompt}>150$ MeV. The fit estimates the KDAR signal selection to be $(77\pm 3)$% pure with $144.4^{+21.3}_{-21.1}$ background counts in the signal region. The total background contribution is estimated to be 99% from non-KDAR neutrinos. Figure 3 shows the KDAR-like selected events and the results of the background estimate, and the inset shows the spectrum for KDAR candidate delayed Michel electron events compared to the MC prediction.

The background estimate and uncertainties are used to estimate the underlying KDAR spectrum, which is then unfolded from reconstructed energy, $E_{prompt}$, to $E_{m}$. The detector response, including detection efficiency, is estimated using MC simulated KDAR events with NuWro and GiBUU and taking their difference as a source of systematic uncertainty. We use Iterative Bayes (D’Agostini) unfolding [41, 42, 43] with a $p>0.9$ between the re-folded and the observed distribution used as the stopping criteria, resulting in 4 iterations for NuWro (central value) and 6 for GiBUU. We also consider a systematic uncertainty based on the scintillator’s proton Birks’ constant, which sets how much of the proton’s energy deposit gets turned into scintillation light. For the detector simulation we use a proton Birks’ constant value of 0.097 mm/MeV, from Ref. [44] and based on a consistent LAB scintillator admixture, but without DIN; we assign a systematic uncertainty to this value based on measurements with pure DIN [45]. We conservatively take the difference between the two Birks’ constant values to be the magnitude of the systematic error, $kB_{proton}=0.097\pm 0.011$ mm/MeV.

Each of the systematic uncertainties are propagated by modifying the MC simulated KDAR events used to produce the unfolding matrix and those used for the non-KDAR neutrino background estimate. The systematically varied response matrix and background spectral shape estimates are used for the background subtraction and unfolding. The difference between the central value unfolded result and the systematically varied result is taken as the corresponding systematic uncertainty on the differential cross section.

Figure 4 shows the observed differential cross section in terms of $E_{m}$ and $E_{vis}$ for the KDAR $\nu_{\mu}$ CC interaction on carbon with statistical and systematic error bars compared to predictions from NuWro (with carbon spectral function), GiBUU, and a relativistic mean field (RMF) prediction [46, 47] using the ACHILLES event generator to simulate FSI [48, 49]. While statistical errors dominate in most of the measurement bins, uncertainties associated with the energy scale and resolution, proton Birks’ constant, and generator uncertainty (affecting the unfolding matrix) contribute significantly in the peak region. The $\chi^{2}=\sum_{i,j}(d_{i}-p_{i})V_{ij}^{-1}(d_{j}-p_{j})$ comparisons between the data ($d$) and predictions ($p$), using the full covariance matrix ($V_{ij}$) are reported in the figure as well. For simplicity, we consider the higher-statistics region of $E_{m}=$585$$ MeV across 16 bins ($i/j$) when forming the $\chi^{2}$. As can be seen, the generator and model predictions disagree with each other substantially and none of them is able to adequately match the data well. While it is difficult to isolate a single cause for this discrepancy across all measurement bins considered, it is much more likely that multiple uncertain nuclear effects are contributing, including FSI, the momentum and separation energy distributions of the nucleons, short-range correlations among nucleons, and even aspects of the weak interaction itself. Notably, however, the data do favor a significant FSI contribution, with a large $\chi^{2}$ improvement after turning on this effect within the NuWro event generator, which best matches the data overall. In general, this measurement serves as a unique and first-of-its-kind benchmark, simultaneously sensitive to multiple aspects of the event, towards elucidating low energy neutrino-nucleus interactions.

Figure 4: The KDAR $\nu_{\mu}$ CC missing energy, $E_{m}$, shape-only differential cross section measurement compared to several neutrino event generator/model predictions. The top x-axis provides the corresponding $E_{vis}$ for each $E_{m}$ value.

III Conclusion

We have presented the first measurement of missing energy for monoenergetic $\nu_{\mu}$ CC interactions based on a $(77\pm 3)$% pure sample of 621 total collected KDAR candidate events. The shape-only differential cross section extracted disagrees significantly with the generator and model predictions considered. This result serves as a known-neutrino-energy, standard candle for improving our understanding of difficult to model low energy neutrino interactions at the transition between on-nucleon and on-nucleus scattering. In the future, JSNS² and JSNS²-II [50] are poised to follow-up on this initial measurement with vastly improved statistics, better control over systematic uncertainties, in particular the proton’s Birks’ constant and energy scale/resolution, and the ability to distinguish KDAR $\nu_{\mu}$ CC events with and without neutrons.

IV Acknowledgements

We deeply thank the J-PARC for their support, especially for the MLF and the accelerator groups to provide an excellent environment for this experiment. We acknowledge the support of the Ministry of Education, Culture, Sports, Science, and Technology (MEXT) and the JSPS grants-in-aid: 16H06344, 16H03967 23K13133, 24K17074 and 20H05624, Japan. This work is also supported by the National Research Foundation of Korea (NRF): 2016R1A5A1004684, 17K1A3A7A09015973, 017K1A3A7A09016426, 2019R1A2C3004955, 2016R1D1A3B02010606, 017R1A2B4011200, 2018R1D1A1B07050425, 2020K1A3A7A09080133, 020K1A3A7A09080114, 2020R1I1A3066835, 2021R1A2C1013661, NRF-2021R1C1C2003615, 2022R1A5A1030700, RS-2023-00212787 and RS-2024-00416839. Our work has also been supported by a fund from the BK21 of the NRF. The University of Michigan gratefully acknowledges the support of the Heising-Simons Foundation. This work conducted at Brookhaven National Laboratory was supported by the U.S. Department of Energy under Contract DE-AC02-98CH10886. The work of the University of Sussex is supported by the Royal Society grant no. IESnR3n170385. We also thank the Daya Bay Collaboration for providing the Gd-LS, the RENO collaboration for providing the LS and PMTs, CIEMAT for providing the splitters, Drexel University for providing the FEE circuits and Tokyo Inst. Tech for providing FADC boards.