# Measuring the Sterile Neutrino Mass in Spallation Source and Direct
Detection Experiments

D. Alonso-González ID<EMAIL_ADDRESS>Instituto de Física Teórica,
IFT-UAM/CSIC, 28049 Madrid, Spain  Departamento de Física Teórica, Universidad
Autónoma de Madrid, 28049 Madrid, Spain  D.W.P. Amaral ID
<EMAIL_ADDRESS>Department of Physics and Astronomy, Rice University,
Houston, TX 77005, USA  A. Bariego-Quintana ID<EMAIL_ADDRESS>Instituto de Física Corpuscular (CSIC - Universitat de València), 46980
Paterna, Valencia, Spain  D. Cerdeño ID<EMAIL_ADDRESS>Instituto de
Física Teórica, IFT-UAM/CSIC, 28049 Madrid, Spain  Departamento de Física
Teórica, Universidad Autónoma de Madrid, 28049 Madrid, Spain  M. de los Rios
ID<EMAIL_ADDRESS>Instituto de Física Teórica, IFT-UAM/CSIC, 28049
Madrid, Spain  Departamento de Física Teórica, Universidad Autónoma de Madrid,
28049 Madrid, Spain

###### Abstract

We explore the complementarity of direct detection (DD) and spallation source
(SS) experiments for the study of sterile neutrino physics. We focus on the
sterile baryonic neutrino model: an extension of the Standard Model that
introduces a massive sterile neutrino with couplings to the quark sector via a
new gauge boson. In this scenario, the inelastic scattering of an active
neutrino with the target material in both DD and SS experiments gives rise to
a characteristic nuclear recoil energy spectrum that can allow for the
reconstruction of the neutrino mass in the event of a positive detection. We
first derive new bounds on this model based on the data from the COHERENT
collaboration on CsI and LAr targets, which we find do not yet probe new areas
of the parameter space. We then assess how well future SS experiments will be
able to measure the sterile neutrino mass and mixings, showing that masses in
the range $\sim 15-50$ MeV can be reconstructed. We show that there is a
degeneracy in the measurement of the sterile neutrino mixing that
substantially affects the reconstruction of parameters for masses of the order
of 40 MeV. Thanks to their lower energy threshold and sensitivity to the solar
tau neutrino flux, DD experiments allow us to partially lift the degeneracy in
the sterile neutrino mixings and considerably improve its mass reconstruction
down to 9 MeV. Our results demonstrate the excellent complementarity between
DD and SS experiments in measuring the sterile neutrino mass and highlight the
power of DD experiments in searching for new physics in the neutrino sector.

††preprint: IFT-UAM/CSIC-23-89

##  I Introduction

The neutrino sector remains one of the most promising places to look for new
physics beyond the Standard Model (SM). Amongst the most obvious open
problems, the SM offers no explanation for the origin of neutrino masses. A
generic prediction of new physics models for neutrino masses is the presence
of new sterile neutrino states, which have very small interactions with the SM
ones. The masses of these new exotic states depend on the actual mechanism by
which neutrinos acquire a mass, but an interesting range of values is the MeV.

The search for sterile neutrinos involves different types of experimental
probes and the constraints depend strongly on the mass range of the new
states. For example, sterile neutrinos have been widely searched for in meson
decays, where masses of up to hundreds of MeV in peak searches of pion and
kaon decays have been probed Britton:1992xv ; Britton:1993cj ; Aguilar-
Arevalo:2019owf ; NA62:2020mcv ; T2K:2019jwa ; NA62:2020mcv , and heavier
steriles have been searched for in neutrino beam dump experiments
Bergsma:1985is ; CooperSarkar:1985nh ; Abreu:1996pa ; Vaitaitis:1999wq . In
our regime of interest (tens of MeV), bounds can be derived through their
possible direct production processes. This could be observed in solar neutrino
data Goldhagen:2021kxe , atmospheric neutrino data Dentler:2018sju , or
neutrino beam experiment data Forero:2021azc  like MINOS/MINOS+ MINOS:2017cae
. In addition, the presence of an extra sterile neutrino may have a non-
negligible impact on different cosmological observations depending on its mass
and couplings Abazajian:2012ys . For example, long-lived sterile neutrinos
with masses of the order of MeV may alter Big Bang nucleosynthesis and the
expansion rate of the universe Sabti:2020yrt ; Abdullahi:2022jlv . Moreover,
sterile neutrinos decaying before recombination may affect the cosmic
microwave background anisotropies Vincent:2014rja ; Boyarsky:2021yoh .

Experiments situated at spallation source (SS) facilities have recently become
excellent probes of new neutrino physics. Most notably, the COHERENT
collaboration COHERENT:2015mry  has been able to observe, for the first time,
a very rare SM phenomenon: the coherent elastic scattering of neutrinos with
nuclei (CE$\nu$NS). The results from both the first run on a CsI target
COHERENT:2017ipa  and a second run that employed LAr in the CENNS-10 detector
COHERENT:2020iec  are compatible with the SM prediction Freedman:1973yd ;
Drukier:1984vhf . This has been used to derive limits on new physics in the
neutrino sector (see, for example, Refs. Kosmas:2017zbh ; Miranda:2020syh ;
Abdullah:2018ykz ; Bauer:2018onh ; Miranda:2020zji ; Bauer:2020itv ), with
particular attention to what future detectors can achieve. Planned experiments
include CENNS610 Akimov:2020pdx  (an extension of CENNS-10 LAr
COHERENT:2019kwz ), CCM ccm , and efforts in the European Spallation Source
facility Baxter:2019mcx . The bounds from COHERENT and the sensitivity of the
planned detectors are generally interpreted in models with low-mass mediators
(or using an effective description in terms of non-standard neutrino
interactions), which alters the SM prediction for CE$\nu$NS Miranda:2020tif ;
Abdullah:2018ykz ; Banerjee:2018eaf ; Papoulias:2019txv ; Khan:2019cvi .
Likewise, they are applicable to inelastic processes that involve the up-
scattering to a heavy neutrino state, for example through the presence of a
nonzero neutrino transition magnetic moment Bolton:2021pey ; Miranda:2021kre ;
DeRomeri:2022twg , or even to a dark fermion Candela:2023rvt .

In parallel, underground experiments searching directly for dark matter
particles have become increasingly sensitive. Planned detectors, especially
those based on liquid noble gases, feature extremely clean, ton-scale targets
with excellent background discrimination that will soon enable them to measure
CE$\nu$NS from solar neutrinos. Although this would constitute a serious
background for dark matter searches, it also offers the unique possibility to
test new neutrino physics Cerdeno:2016sfi ; Dutta:2017nht ; Gelmini:2018gqa ;
Essig:2018tss ; Amaral:2020tga ; Dutta:2020che ; Amaral:2021rzw ;
Munoz:2021sad ; deGouvea:2021ymm  in a way that is complementary to that of
dedicated neutrino detectors. The main advantages of these direct detection
(DD) experiments are that they can probe both electron and nuclear recoils,
which makes them a perfect complement to SS and oscillation experiments
Amaral:2023tbs , and that they are also sensitive to the tau neutrinos in the
solar flux.

The sensitivity of DD experiments to observe heavy neutrino states was studied
in Ref. Shoemaker:2018vii  for the particular case of the neutrino dipole
portal, showing that current xenon-based detectors could significantly improve
existing astrophysical bounds. The neutrino dipole portal was considered to
account for the apparent excess in the low-energy data from electronic recoils
in the XENON1T experiment XENON:2020rca ; Shoemaker:2020kji . However, this
solution was seriously limited by other experimental constraints Brdar:2020quo
, and the excess was not reproduced by XENONnT XENON:2022ltv . Since the
coupling of a sterile neutrino to the leptonic sector is in general severely
limited by experimental searches, in this article we will focus on the
potential interactions with the quark sector. These are more difficult to
probe, but they could lead to changes in the predicted nuclear recoil rates in
DD and SS experiments that could be accessible in near future experiments. For
concreteness, in this work we set up to study the sterile baryonic neutrino
(SBN) Pospelov:2011ha  as an example of models in which the active neutrinos
can up-scatter to heavy states.

More specifically, in this article we study the potential of DD and SS
experiments to not only detect the sterile neutrino but also reconstruct its
parameters—namely, its mass and mixings with the active neutrinos. Our main
goal is to determine the conditions under which the sterile neutrino mass can
be unambiguously measured (distinguished from zero).

In Section II, we introduce an effective construction based on the sterile
baryonic neutrino model and determine the new inelastic contribution to
neutrino-nucleus scattering. In Section III, we address the prospects for
upcoming SS experiments. In Section IV, we extend the analysis to include
future xenon-based DD experiments. Finally, in Section V, we study the
complementary role of DD and SS experiments. We present our conclusions in
Section VI.

##  II The Sterile Baryonic Neutrino

We introduce a dark sector consisting of a new vector mediator, $Z^{\prime}$,
stemming from a broken $\mathrm{U}(1)_{B}$ gauge symmetry and a new baryonic
sterile neutrino, $\nu_{b}$, that is also charged under this new symmetry
Pospelov:2011ha . For the purpose of this work, we regard this model as an
effective theory, and we do not address its possible anomaly-free UV
completion. The relevant part of our Lagrangian is given by

$\mathcal{L}\supset\frac{m_{Z^{\prime}}^{2}}{2}Z^{\prime\mu}Z^{\prime}_{\mu}+g_{b}Z^{\prime\mu}\overline{\nu}_{b}\gamma_{\mu}\nu_{b}+\frac{1}{3}g_{q}Z^{\prime\mu}\sum_{q}\overline{q}\gamma_{\mu}q\,.$
(1)

Here, $m_{Z}^{\prime}$ is the mass of the new boson, $g_{b}$ is its gauge
coupling to the baryonic neutrino and $g_{q}$ to the quarks, and the sum runs
over all quark flavours $q$. In this model, a generic flavour eigenstate,
$\ket{\nu_{\alpha}}$, can then be written as a linear combination of mass
eigenstates, $\ket{\nu_{i}}$, as

$\ket{\nu_{\alpha}}=\sum_{i=1}^{4}U_{\alpha i}^{*}\ket{\nu_{i}}\,,$ (2)

where $\ket{\nu_{4}}$ is the new mass eigenstate with mass $m_{4}$, and
$\alpha\in\\{e,\,\mu,\,\tau,\,b\\}$.

From Eq. 1, and defining the coupling $g_{Z^{\prime}}\equiv\sqrt{g_{b}g_{q}}$,
the neutrino-nucleus up-scattering process $\nu_{\alpha}A\rightarrow\nu_{4}A$
has amplitude

$\mathcal{M}_{\alpha
4}=\frac{g_{Z^{\prime}}^{2}}{q^{2}-m_{Z^{\prime}}^{2}}l^{\mu}h_{\mu}\,,$ (3)

where $q^{2}$ is the square-momentum exchange with the nucleus, $h^{\mu}$ is
the nucleus transition amplitude for the nuclear ground state $A$, and
$l^{\mu}$ is the leptonic transition amplitude. Using Eq. 2 to re-write the
dark baryonic current in terms of the mass eigenstates, we have that

$\begin{split}l^{\mu}\equiv\bra{\nu_{4}}\overline{\nu}_{b}\gamma_{\mu}\nu_{b}\ket{\nu_{\alpha}}&=\sum_{ijk}\bra{\nu_{4}}U_{\alpha
k}^{*}U_{bi}^{*}U_{bj}\overline{\nu}_{j}\gamma_{\mu}\nu_{i}\ket{\nu_{k}}=\sum_{i}U_{\alpha
i}^{*}U_{b4}^{*}U_{bi}\bra{\nu_{4}}\overline{\nu}_{4}\gamma_{\mu}\nu_{i}\ket{\nu_{i}}\\\
&\simeq U_{\alpha
4}^{*}\bra{\nu_{4}}\overline{\nu}_{4}\gamma_{\mu}\nu_{i}\ket{\nu_{i}}\,,\end{split}$
(4)

where, in the last step, we have assumed that
$\left|U_{bi}\right|\ll\left|U_{b4}\right|$ for $i\not=4$ and that
$\left|U_{b4}\right|^{2}\simeq 1$ Arguelles:2022lzs . The differential
neutrino-nucleus up-scattering cross section then follows:

$\diff{\sigma_{\alpha 4}}{E_{R}}=\frac{g_{Z^{\prime}}^{4}A^{2}\left|U_{\alpha
4}\right|^{2}m_{A}}{2\pi
E_{\nu}^{2}{\left(2m_{A}E_{R}+m_{Z^{\prime}}^{2}\right)}^{2}}\left[4E_{\nu}^{2}-2E_{R}\left(m_{A}-E_{R}+2E_{\nu}\right)-\frac{m_{4}^{2}}{m_{A}}\left(m_{A}-E_{R}-E_{\nu}\right)\right]F^{2}(E_{R}),\,$
(5)

where $m_{A}$ is the mass of the target nucleus, $E_{\nu}$ is the energy of
the incoming neutrino, and $E_{R}$ is the nuclear recoil energy. For the
nuclear form factor $F^{2}(E_{R})$, which arises from the hadronic part of the
amplitude, we use the Helm form factor Helm:1956zz  with the parametrisation
introduced in Ref. Lewin:1995rx . This new inelastic scattering process
provides an extra contribution to the usual SM elastic neutrino-nucleus
scattering, which takes place through CE$\nu$NS and has the following
differential cross section,

$\diff{\sigma_{{\rm CE}\nu{\rm
NS}}}{E_{R}}=\frac{G_{F}^{2}}{4\pi}Q_{\nu}^{2}m_{A}\left(1-\frac{m_{A}E_{R}}{2E_{\nu}^{2}}\right)F^{2}(E_{R})\
,$ (6)

where $G_{F}$ is the Fermi constant, and $Q_{\nu}\equiv
N-\left(1-4\sin^{2}\theta_{W}\right)Z$ is the SM coherence factor in terms of
the Weinberg angle, $\theta_{W}$, and the number of neutrons, $N$, and
protons, $Z$.

Note that, for the characteristic recoil energies at SS experiments
($E_{R}\lesssim$100\text{\,}\mathrm{keV}$$) and DD experiments
($E_{R}\lesssim$10\text{\,}\mathrm{keV}$$), the cross section in Eq. 5 can be
interpreted as being proportional to the effective coupling
$g_{Z^{\prime}}^{4}|U_{\alpha 4}|^{2}/m_{Z^{\prime}}^{4}$. As both of these
types of experiments are sensitive to this product of model parameters, they
are only able to make inferences on this effective coupling. Since the focus
of our analysis is the physics underlying the baryonic neutrino, we choose to
fix the parameters related to the new vector mediator to
$m_{Z^{\prime}}=$1\text{\,}\mathrm{GeV}$$ and $g_{Z^{\prime}}=4\times
10^{-3}$, taking into account the constraints found in Ref. Foguel:2022ppx .
Thus, without loss of generality, for as long as $m_{Z^{\prime}}^{2}$ remains
greater than the momentum transfer at these experiments, our results can
simply be rescaled by the factor $g_{Z^{\prime}}^{4}/m_{Z^{\prime}}^{4}$. We
therefore consider a four-dimensional parameter space
$(m_{4},\,|U_{e4}|^{2},|U_{\mu 4}|^{2},|U_{\tau 4}|^{2})$ and Table 1 shows
some representative benchmark points used in this work.

|  $m_{4}$ [GeV] |  $\left|U_{e4}\right|^{2}$ |  $\left|U_{\mu 4}\right|^{2}$ |  $\left|U_{\tau 4}\right|^{2}$
---|---|---|---|---
BP1a |  $2\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  0
BP1d |  $2\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  $9\times 10^{-3}$
BP2a |  $9\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  0
BP2b |  $9\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  $9\times 10^{-4}$
BP2c |  $9\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  $4\times 10^{-3}$
BP2d |  $9\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  $9\times 10^{-3}$
BP3a |  $20\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  0
BP4a |  $40\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  0
BP5a |  $60\times 10^{-3}$ |  0 |  $9\times 10^{-3}$ |  0
Table 1: Benchmark points used in this work.

##  III Spallation source experiments

Neutrino experiments at spallation sources have become an extremely useful
tool to explore new neutrino physics associated with neutrino-nucleus
scattering. The neutrino flux arriving on-target has three components, shown
in Fig. 1. The prompt decay of the initially produced pions,
$\pi^{+}\rightarrow\mu^{+}\nu_{\mu}$, induces a monochromatic beam of muon
neutrinos with energy
$E_{\nu_{\mu}}=(m_{\pi}^{2}-m_{\mu}^{2})/{2m_{\pi}}\simeq$30\text{\,}\mathrm{MeV}$$.
The delayed decay $\mu^{+}\rightarrow e^{+}\nu_{e}\bar{\nu}_{\mu}$ gives rise
to a flux of muon antineutrinos and electron neutrinos with continuous energy
distributions. The corresponding fluxes are given by (see, e.g., Ref.
Coloma:2017egw )

$\displaystyle\diff{\phi_{\nu_{\mu}}}{E_{\nu}}$
$\displaystyle=\xi\delta\left(E_{\nu}-\frac{m_{\pi}^{2}-m_{\mu}^{2}}{2m_{\pi}}\right)\,,$
(7) $\displaystyle\diff{\phi_{\bar{\nu}_{\mu}}}{E_{\nu}}$
$\displaystyle=\xi\frac{64}{m_{\mu}}\left[\left(\frac{E_{\nu}}{m_{\mu}}\right)^{2}\left(\frac{3}{4}-\frac{E_{\nu}}{m_{\mu}}\right)\right]\,,$
$\displaystyle\diff{\phi_{\nu_{e}}}{E_{\nu}}$
$\displaystyle=\xi\frac{192}{m_{\mu}}\left[\left(\frac{E_{\nu}}{m_{\mu}}\right)^{2}\left(\frac{1}{2}-\frac{E_{\nu}}{m_{\mu}}\right)\right]\,,$

where, from kinematics, $E_{\nu}\in\left[0,m_{\mu}/2\right]$ for the
continuous spectra of $\bar{\nu}_{\mu}$ and ${\nu}_{e}$. The constant
${\xi\equiv rR_{\textup{POT}}/(4\pi L^{2})}$ accounts for the luminosity of
the experiment. Here, $r$ is the number of neutrinos of any given flavour
produced per proton collision, $R_{\textup{POT}}$ is the number of protons on
target per unit time, and $L$ is the total length of the experimental
baseline. Given the promising sensitivity of the configurations planned to run
at the European Spallation Source, in this article we will consider it as a
paradigmatic example of a realistic future experiment. Two different setups
can be considered Baxter:2019mcx : a small (10 kg) but extremely sensitive
detector with an energy threshold of $E_{\rm th}=0.1$ keV (which we refer to
as ESS10), and a large detector (1 ton) but with a higher energy threshold of
$E_{\rm th}=20$ keV (which we refer to as ESS). For both configurations, the
baseline is $L=20$ m, $R_{\textup{POT}}=2.8\times 10^{23}$ yr-1, and $r=0.3$.
Despite the great advantage of its extremely low threshold, the small target
size of ESS10 makes it insufficient to explore new regions of the parameter
space of sterile neutrino models, and, for this reason, we will concentrate on
ESS assuming 1 yr of operation. In our analysis, we consider a bin energy
resolution of $5$ keV. For the quenching factor, we have extrapolated that of
COHERENT-LAr COHERENT:2020iec , $Q_{F}=0.246+7.8\times 10^{-4}E_{R}$, whereby
$E[$\mathrm{keV}$\textsubscript{ee}{}]=Q_{F}\,E_{R}$. Following the treatment
in Ref. Amaral:2021rzw , we approximate the efficiency as
$\epsilon(E_{R})=0.5\left(1+\tanh\left(E_{R}-E_{\rm th}\right)\right)/E_{\rm
width}$, where we take $E_{\rm width}=1$ keV for ESS.

Figure 1: Normalised (without the $\xi$ factor) neutrino fluxes at spallation
source facilities.

To compute the differential rate of nuclear recoil events, we integrate each
neutrino flux, $\alpha^{\prime}\in\left\\{e,\mu,\bar{\mu}\right\\}$, taking
into account both SM CE$\nu$NS and new physics up-scattering processes, from
Eq. 6 and Eq. 5, respectively. The differential scattering rate is given by

$\diff{R_{\alpha^{\prime}}}{E_{R}}=\frac{1}{m_{A}}\left(\int_{E_{\nu}^{{\rm
min,CE}\nu{\rm NS}}}^{E_{\nu}^{\rm
max}}\diff{\phi_{\nu_{\alpha^{\prime}}}}{E_{\nu}}\diff{\sigma_{{\rm CE}\nu{\rm
NS}}}{E_{R}}\,\dl E_{\nu}\,+\int_{E_{\nu}^{\rm
min,\alpha^{\prime}4}}^{E_{\nu}^{\rm
max}}\diff{\phi_{\nu_{\alpha^{\prime}}}}{E_{\nu}}\diff{\sigma_{\alpha^{\prime}4}}{E_{R}}\,\dl
E_{\nu}\,\right),$ (8)

where $1/m_{A}$ is the total number of targets per unit mass in a given
experiment,
$\mathrm{d}{\sigma_{\bar{\mu}4}}/\mathrm{d}{E_{R}}=\mathrm{d}\sigma_{\mu
4}/{\mathrm{d}E_{R}}$, and $E_{\nu}^{\rm max}=m_{\mu}/2$ is the maximum
allowed neutrino energy. The minimum neutrino energy required to produce a
recoil of energy $E_{R}$ differs for the elastic and inelastic processes. For
usual SM CE$\nu$NS, it is given by

$E_{\nu}^{{\rm min,\,CE}\nu{\rm
NS}}=\frac{1}{2}\left(E_{R}+\sqrt{E_{R}^{2}+2m_{A}E_{R}}\right)\simeq\sqrt{\frac{m_{A}E_{R}}{2}}\,.$
(9)

However, for the inelastic up-scattering process, the minimum energy must be
high enough to produce the massive sterile neutrino, leading to

$E_{\nu}^{\rm
min,\,\alpha^{\prime}4}=\left(1+\frac{m_{4}^{2}}{2m_{A}E_{R}}\right)E_{\nu}^{\min,\,\mathrm{CE\nu
NS}}\,.$ (10)

Finally, the total number of nuclear recoils in each energy bin is computed by
integrating the differential rate over the experimental range of recoil
energies (given by the specific experimental setup) weighted by the
corresponding energy-dependent efficiency function, $\epsilon(E_{R})$,

$N_{\rm SS}=\varepsilon\sum_{\alpha^{\prime}}\int_{E_{R}^{\rm
min}}^{E_{R}^{\rm max}}\diff{R_{\alpha^{\prime}}}{E_{R}}\epsilon(E_{R})\,\dl
E_{R}\ .$ (11)

where $\varepsilon$ is the experiment exposure: the product of its total mass
and its live time. For the ESS configuration that we are considering,
$\varepsilon=1\,\mathrm{ton}\,\mathrm{yr}$.

Figure 2: Predicted energy spectra at ESS for some representative benchmark
points in Table 1, featuring different values for $m_{4}$ and mixings with the
active neutrinos. The vertical grey dotted (dashed-dotted) line shows the
projected ESS (ESS10) threshold. The quantity $\Delta_{\mu}$ is defined as the
energy difference between the endpoint of the SM spectrum and the contribution
from the monochromatic $\nu_{\mu}$ flux. The dashed green line in the upper-
left panel shows the expected neutrino NSI spectrum with
$\varepsilon_{\mu\mu}^{u}=0.4$.

Fig. 2 shows the differential spectrum for each contribution in Eq. 8 and for
four representative benchmark points (BP1a, BP2a, BP3a, and BP5a with
parameters specified in Table 1), where the sterile neutrino mass is varied
for the same choice of couplings. The inelastic contribution only switches on
above a certain recoil energy, leading to a characteristic bump with energies
in the range

$\displaystyle E_{R}^{\textup{bump}}\in$
$\displaystyle\left[\frac{1}{2m_{A}}\left(2(E_{\nu}^{\rm
max})^{2}-m_{4}^{2}-2E_{\nu}^{\rm max}\sqrt{(E_{\nu}^{\rm
max})^{2}-m_{4}^{2}}\right)\right.\,,$ (12)
$\displaystyle\left.\frac{1}{2m_{A}}\left(2(E_{\nu}^{\rm
max})^{2}-m_{4}^{2}+2E_{\nu}^{\rm max}\sqrt{(E_{\nu}^{\rm
max})^{2}-m_{4}^{2}}\right)\right]\,,$ (13)

where we have made the approximation $E_{\nu}/m_{A}\ll 1$. In the event of a
future observation, this ‘bump’ could be used to determine the mass of the
sterile neutrino, thus helping to discriminate this model from other potential
new physics contributions in the neutrino sector. In practice, this could
confirm the existence of a sterile neutrino (with mass different from zero).
Notice that the lower end of the energy bump takes place at very small values
of the recoil energy, well below the reach of current and future detectors.
For this reason, the sterile neutrino mass reconstruction mostly relies on
determining the upper end of the bump, which is displaced from the end of the
SM CE$\nu$NS spectrum. The contribution from muon neutrinos is particularly
interesting for this purpose. As their flux is monochromatic, the energy bump
in their spectrum is more easily distinguishable from the SM prediction. The
difference of the endpoint in the SM CE$\nu$NS spectrum and the inelastic
contribution from $\nu_{\mu}$ is denoted $\Delta_{\mu}$ in Fig. 2 for each
benchmark point.

To observe this feature, the experimental threshold must be low enough and the
energy resolution of the detector must at least be comparable to
$\Delta_{\mu}$. Since $\Delta_{\mu}$ increases with $m_{4}$ (which we can see
in Fig. 2 or infer from Eq. 13), heavier sterile neutrino masses are easier to
reconstruct. Since the energy thresholds of current and planned experiments at
spallation sources are of the order of $\sim 10$ keV, a measurement of the
sterile neutrino mass is only possible above a certain value of $m_{4}$. In
particular, given the planned characteristics of the ESS experiment, the
signal of both BP1 and BP2 would be indistinguishable from that for $m_{4}=0$.
For reference, the vertical grey dotted (dashed-dotted) lines in Fig. 2
represent the expected energy threshold of both ESS and ESS10 respectively.

It should be emphasized that measuring the sterile neutrino mass—that is,
confirming that $m_{4}=0$ is not within the $2\sigma$ best-fit region—is
crucial to discriminate the signal due to the SBN model from that of a generic
neutrino non-standard interaction (NSI), where no extra neutrinos are
introduced Wolfenstein:1977ue ; Guzzo:1991cp ; Guzzo:1991hi ; Gonzalez-
Garcia:1998ryc ; Bergmann:2000gp ; Guzzo:2000kx ; Guzzo:2001mi ;
GonzalezGarcia:2004wg . Indeed, the spectrum from a particular choice of NSI
can mimic the observed signal in the SBN model when the lower end of the
energy bump is below the experimental threshold. We illustrate this in Fig. 2
for BP1a, where we have generated an NSI spectrum with a pure up-quark
effective NSI parameter of $\varepsilon_{\mu\mu}^{u}=0.4$. For the range of
observable energies, we see that the SBN and NSI spectra almost completely
overlap, making them indistinguishable from one another.

Figure 3: Profile-likelihood results for the BPs with $\left|U_{\mu
4}\right|^{2}=9\times 10^{-3}$, $\left|U_{\tau
4}\right|^{2}=\left|U_{e4}\right|^{2}=0$ and $m_{4}=2$ MeV (upper left panel),
$9$ MeV (upper right panel), $20$ MeV (lower left panel), and $60$ MeV (lower
right panel) using SS data only. Hatched areas correspond to allowed regions
($\Delta\chi^{2}<6.18$). For clarity, we depict as a yellow star the true
values of the analysed benchmark point. The shaded black (grey) regions areas
excluded by the current COHERENT data. The horizontal dashed blue line
corresponds to the upper bound on the sterile neutrino mixing with the muon
sector (Arguelles_2023, ).

To test the reconstruction of the sterile neutrino parameters, we have created
Asimov data sets for each of these benchmark points and attempted to
reconstruct their associated model parameters in the four-dimensional space
$(m_{4},\left|U_{e4}\right|^{2},\left|U_{\mu 4}\right|^{2},\left|U_{\tau
4}\right|^{2})$. In these Asimov sets, our ‘observed’ data are equal to the
theoretically expected number of events for each given benchmark point. The
ensuing limit from such an analysis should asymptotically approach the median
limit arising from many Monte Carlo runs Cowan:2010js . The statistical
details of our analysis can be found in Appendix A. We compute the expected
number of nuclear recoil events from Eq. 11 using an extension of the SNuDD
package SNuDD-code . For each benchmark point, we carry out a profile-
likelihood analysis using the nested sampling algorithm multinest multinest ;
multinest2  via its Python implementation pymultinest .

We show in Fig. 3 the parameter reconstruction corresponding to BP1a, BP2a,
BP3a, and BP5a, assuming the projected configuration of the ESS detector. The
hatched areas correspond to the allowed regions ($\Delta\chi^{2}<6.18$). As we
can see, ESS would be able to observe the first three benchmark points and
measure the coupling $\left|U_{\mu 4}\right|^{2}$. It would also be able to
fully reconstruct the mass of the sterile neutrino in BP3a. Nevertheless, for
BP1a and BP2a, only an upper bound on the sterile neutrino mass can be
extracted (the end-point of the bump cannot be distinguished from the SM
spectrum). Since the sterile neutrino mass for BP5a is above the energy of the
neutrino flux in spallation source experiments, the up-scattering is
kinematically forbidden and hence there will be no observation. For this
benchmark point, we can only obtain an exclusion region.

As a new result, we have derived constraints on the SBN model using current
COHERENT data from the two targets, LAr COHERENT:2020iec  and CsI
COHERENT:2021xmm . To do this, we have used the statistical treatment of
Appendix A. The bounds are represented in Fig. 3 as light and dark grey areas
in the corresponding plots for the LAr and CsI targets, respectively. As we
can see, the excluded areas lie above the upper bound on the sterile neutrino
mixing with the muon sector from Ref. Arguelles_2023  and therefore do not
probe new areas of the parameter space.

It is interesting to note that for sterile neutrino masses above $m_{4}\gtrsim
30$ MeV, the monochromatic $\nu_{\mu}$ flux is not energetic enough to produce
the sterile neutrino and only the $\bar{\nu}_{\mu}$ and $\nu_{e}$ fluxes
contribute in Eq. 8. When this occurs, the characteristic feature
$\Delta_{\mu}$ is no longer present. This makes the mass reconstruction more
difficult and leads to a degeneracy between the mixings with muon neutrinos,
$U_{\mu 4}$, and electron neutrinos, $U_{e4}$. This effect is more pronounced
for $m_{4}\simeq 40$ MeV, where the $\nu_{e}$ and $\bar{\nu}_{\mu}$ fluxes are
comparable. To exemplify this, in Fig. 4 we analyse a benchmark point with
$m_{4}=40$ MeV and $\left|U_{\mu 4}\right|^{2}=9\times 10^{-3}$ (BP4a in Table
1), which we attempt to reconstruct through a profile-likelihood analysis. The
degeneracy on the reconstruction of the mixings (evidenced on the right panel)
induces a similar degeneracy on the sterile neutrino mass (see left and middle
panels of Fig. 4), making measuring $m_{4}$ impossible. This degeneracy is
lifted for sterile neutrino masses $m_{4}\gtrsim$45\text{\,}\mathrm{MeV}$$
(depending on the value of the mixings) when the contributions from the
$\nu_{e}$ and $\nu_{\mu}$ fluxes differ (see Fig. 1).

Figure 4: Profile-likelihood results for BP4a ($m_{4}=40$ MeV, $\left|U_{\mu
4}\right|^{2}=9\times 10^{-3}$, $\left|U_{\tau
4}\right|^{2}=\left|U_{e4}\right|^{2}=0$) using SS data only. For clarity, we
depict as a yellow star the true values of the analysed benchmark point. The
shaded black (grey) regions are excluded by current COHERENT data. The yellow
arrows depict that the benchmark point is outside of the plotting area.

Our analysis so far shows that

  * •

Current limits on the SBN model using COHERENT data do not exclude new areas
of the parameter space, but future experiments like ESS would allow us to
explore regions below current experimental constraints.

  * •

In the event of a positive observation, future SS experiments might be able to
determine the sterile neutrino mass (distinguishing it from the massless case)
for a range $m_{4}\sim 15-50$ MeV. For lighter masses, the observed signal is
indistinguishable from that of a new massless neutrino.

  * •

The sterile neutrino mixing with the electron and muon sectors can, in
general, be disentangled based on the different shapes of the contribution
from the $\nu_{e}$ and $\nu_{\mu}$ fluxes.

  * •

There is, however, a region for sterile neutrino masses around
$m_{4}\sim$40\text{\,}\mathrm{MeV}$$ for which the reconstruction is highly
degenerate and the sterile neutrino mass (and mixing with $\nu_{e}$ and
$\nu_{\mu}$) cannot be measured.

  * •

SS experiments are completely insensitive to the sterile neutrino mixing with
the tau sector, as there is no $\nu_{\tau}$ flux.

In the following sections, we will study how (dark matter) direct detection
experiments can provide complementary information that improves the
reconstruction of the SBN model parameters, partially lifting some of these
degeneracies and considerably improving the mass measurement.

##  IV Direct Detection experiments

While primarily employed in the search for dark matter, direct detection
experiments are becoming so sensitive that they will start observing CE$\nu$NS
from solar neutrinos. Indeed, the sensitivities of xenon-based experiments of
this and future generations—such as LZ LZ:2018qzl , XENONnT XENON:2020kmp ,
and DARWIN DARWIN:2016hyl —are projected to hit the neutrino fog: a region of
the parameter space where a dark matter signal and a neutrino event will be
difficult to disentangle OHare:2021utq . This motivates us to think of these
experiments as neutrino observatories instead of as dark matter detectors,
treating this ‘background’ as a signal to help us learn more about the nature
of both SM and BSM neutrino physics. In this section, we show how these
experiments can use measurements of the solar neutrino scattering rate as a
probe of the SNB model.

In the case of nuclear recoils, the calculation of the differential rate is
similar to that of SS. The key differences are that we instead use the solar
neutrino flux and that we must now account for the oscillation probabilities
as neutrinos propagate to the Earth from the solar core. As we did in Section
III, the SM and new inelastic contributions must be considered separately
since the minimal neutrino energy to produce a nuclear recoil of a given
energy differs. The differential scattering rate, after summing over the
flavours $\alpha\in\\{e,\,\mu,\,\tau\\}$, is ultimately given by111It has
recently been noted that one must be careful when calculating the solar
neutrino scattering rate in the presence of new physics Coloma:2022umy . If
the new physics introduces flavour-changing neutral current processes, then a
more general density matrix formalism must be employed. This was recently done
in the context of DD experiments and general NSI in Ref. Amaral:2023tbs . In
our case, flavour charge is conserved, so we can compute the rate in the usual
manner as we have written.

$\displaystyle\diff{R}{E_{R}}=\frac{1}{m_{A}}$
$\displaystyle\left[\int_{E_{\nu}^{{\rm min,CE}\nu{\rm NS}}}^{E_{\nu}^{\rm
max}}\diff{\phi_{\nu_{e}}}{E_{\nu}}\diff{\sigma_{{\rm CE}\nu{\rm
NS}}}{E_{R}}\,\dl E_{\nu}\,+\sum_{\alpha}\int_{E_{\nu}^{\rm min,\alpha
4}}^{E_{\nu}^{\rm max}}\diff{\phi_{\nu_{e}}}{E_{\nu}}\ P_{e\alpha}\
\diff{\sigma_{\alpha 4}}{E_{R}}\,\dl E_{\nu}\right]\,,$ (14)

where ${\mathrm{d}\phi_{\nu_{e}}}/{\mathrm{d}E_{\nu}}$ is the total
differential solar electron-neutrino flux and $P_{e\alpha}$ is the transition
probability for an electron neutrino to oscillate to the flavour $\alpha$.
Notice that since SM CE$\nu$NS is flavour blind, the transition probabilities
factor out and sum to one. For the new physics contribution, the cross section
is instead flavour dependent, so the probabilities must be retained.

In this work, we consider a multi-ton xenon experiment with an exposure of
$\varepsilon=200\,\mathrm{ton}\,\mathrm{yr}$, a recoil energy threshold of
$E_{\mathrm{th}}=$1\text{\,}\mathrm{keV}$$, and an energy bin resolution of
$1\text{\,}\mathrm{keV}$. This type of experiment has been shown to be a
powerful probe of new physics in the neutrino sector Amaral:2023tbs ;
Amaral:2021rzw ; Amaral:2020tga . When calculating the total number of
expected events, we incorporate experimental effects, folding into Eq. 14 the
energy-dependent efficiency and resolution functions. We do this using

$N_{\mathrm{DD}}=\varepsilon\int_{0}^{E_{\mathrm{max}}}\Bigg{(}\int\diff{R}{E^{\prime}}\epsilon(E^{\prime})\frac{1}{\sigma(E^{\prime})\sqrt{2\pi}}\exp\left[-\frac{\left(E_{R}-E^{\prime}\right)^{2}}{2\sigma^{2}(E^{\prime})}\right]\,\operatorname{d\\!}E^{\prime}\Bigg{)}\operatorname{d\\!}E_{R}\,,$
(15)

where the convolution with the Gaussian resolution function is taken with
respect to the theoretically expected recoil energy, $E^{\prime}$, which is
converted to the observed recoil energy, $E_{R}$. The integral is taken from
$E_{R}=0$, with the threshold of the experiment implicitly incorporated
through the efficiency function, $\epsilon$. Note that it is crucial to
incorporate this convolution with the resolution function, as this smears
lower energy ${}^{8}\mathrm{B}$ events beyond where CE$\nu$NS would be
kinematically forbidden. As experimental thresholds are typically placed near
where this forbidden region occurs, which is useful for dark matter searches,
this smearing allows us to see some events as opposed to almost no events.

To implement Eqs. 14 and 15, we once again make use of the SNuDD package. This
package uses the B16-GS98 standard solar model neutrino flux predictions
Vinyoles:2016djt  and the NuFIT 5.2 oscillation parameter results to compute
the electron neutrino survival and transition probabilities Esteban:2020cvm .
For more information on the package, please see Ref. Amaral:2023tbs  for the
theory and Ref. SNuDD-code  for the code base.

With the existence of the new flavour state $\ket{\nu_{b}}$, it is possible
that the electron neutrinos produced in the Sun can oscillate into baryonic
neutrinos. These neutrinos could then elastically scatter off target nuclei
via the new vector mediator, leading to an observable signal in DD experiments
that could, in principle, dominate over that of our considered inelastic
process Pospelov:2011ha ; Pospelov:2012gm . However, for sterile neutrinos in
the mass range we have considered
($m_{4}\sim$1\text{\,}\mathrm{MeV}$\text{--}$100\text{\,}\mathrm{MeV}$$)
deviations from the unitarity of the PMNS matrix are highly constrained by
flavour and electroweak precision data, as well as direct searches for such
heavy neutrino states Abdullahi:2022jlv . Consequently, we take the liberty of
ignoring transitions to the baryonic neutrino state, neglecting the elastic
scattering process and using the SM prediction for the survival and transition
probabilities.

Figure 5: Predicted DD spectra for four benchmark points with $m_{4}=9$ MeV,
$\left|U_{\mu 4}\right|^{2}=9\times 10^{-3}$, $\left|U_{e4}\right|^{2}=0$, and
$\left|U_{\tau 4}\right|^{2}=0~{}(\textit{upper\ left\ panel}),\ 9\times
10^{-4}~{}(\textit{upper\ right\ panel}),\ 4\times 10^{-3}~{}(\textit{lower\
left\ panel}),\ 9\times 10^{-3}~{}(\textit{lower\ rigth\ panel})$. The SM
spectrum is shown in grey, while the SBN contribution is shown in black. For
completeness, we also show as red dashed-dotted (blue dashed) lines the SBN
contributions arising from the $\nu_{\mu}$ ($\nu_{\tau}$) flux.

Fig. 5 shows the resulting differential spectrum for some representative
benchmark points from Table 1. As in the case of SS experiments, the new
physics contribution from the inelastic process shows a characteristic bump.
There is, however, an important difference. Since the solar neutrino fluxes
are not monochromatic, this feature is not as abrupt as the $\nu_{\mu}$
contribution in SS experiments. Consequently, the reconstruction of the
sterile neutrino mass from a hypothetical future signal in DD experiments is
significantly more challenging. Notice that the lower end of the energy bump
is generally well below the experimental threshold (and is therefore not
observable). Thus, it is difficult to determine a lower bound on the mass of
the sterile neutrino using DD alone. Given the shape of the solar neutrino
flux Bahcall:2004pz , for sterile neutrino masses above
$\sim$2\text{\,}\mathrm{MeV}$$, only the 8B and $hep$ neutrino fluxes
contribute to the inelastic process. Despite this, DD experiments have the
great advantage that they are sensitive to all three flavours of active
neutrinos, thereby conveniently complementing the information from spallation
sources, which lack a tau neutrino flux.

As we did for SS experiments, we can compare the expected number of events for
a given set of model parameters with the simulated data of each benchmark
point detailed in Table 1. Since the expected number of events is
significantly lower than in SS experiments, we model the likelihood as a
product of Poissonian likelihoods for each energy bin. In addition, we
introduce a nuisance parameter to account for the systematic uncertainty on
the ${}^{8}\mathrm{B}$ flux. The full statistical description can be found in
Appendix A. To test how this uncertainty impacts our results, we consider two
cases222These values are motivated by the current uncertainty obtained through
global fits analysis Bergstrom:2016cbh  ($\sigma_{{}^{8}\mathrm{B}}=2\%$) and
the uncertainty to which DUNE will measure $\mathrm{{}^{8}B}$ using a
combination of elastic scattering and charged-current interactions
($\sigma_{{}^{8}\mathrm{B}}=2.5\%$) Capozzi:2018dat .: one with the current
experimental uncertainty of $\sigma_{{}^{8}\mathrm{B}}=4\%$ (SNO:2011hx, ) and
another one with an optimistic uncertainty of $\sigma_{{}^{8}\mathrm{B}}=1\%$.

In Fig. 6, we show as blue hatched regions the parameters that would be
allowed ($\Delta\chi^{2}<6.18$) by a future observation in a multi-ton liquid
xenon experiment with $\sigma_{{}^{8}\mathrm{B}}=1\%$. For comparison, we
include as a blue dashed line the results obtained with
$\sigma_{{}^{8}\mathrm{B}}=4\%$. Given the maximum energy of the 8B solar
neutrino flux, DD experiments will be insensitive to BP3a and BP5a. Hence, DD
experiments can only probe sterile neutrinos with a low mass ($m_{4}\lesssim
20$ MeV) and a large mixing. Regarding the benchmark points of Fig. 6, only
BP1a is observable—while we do observe events for BP2a, the statistics are not
high enough for a reconstruction. For BP2a, BP3a, and BP5a we only obtain an
upper bound on the neutrino mixing. For BP5a, adding DD data leads to a more
constraining upper bound for small sterile neutrino masses. It should be
emphasised that one cannot disentangle the individual contributions from each
of the three neutrino flavours using only DD data, and therefore the
reconstruction of the mixing parameters is completely degenerate (in the
figure, this leads to $|U_{\mu 4}|^{2}$ being unbounded).

##  V The Complementarity of direct detection and spallation source
Experiments

In this section, we forecast the sensitivity that will be achieved by
combining the results of future DD and SS experiments. In particular, we
analyse how their complementarity can be used to break the degeneracies found
in their individual analyses and better determine the parameters of the SBN
model. Since the measurements performed by DD and SS experiments are
independent of one another, we model the total likelihood as the product of
the individual likelihoods described in Appendix A. Using this combined
likelihood, we repeat our previous analysis.

Figure 6: Profile-likelihood results for the BPs with $\left|U_{\mu
4}\right|^{2}=9\times 10^{-3}$, $\left|U_{\tau
4}\right|^{2}=\left|U_{e4}\right|^{2}=0$, and $m_{4}=2$ MeV (upper left
panel), $9$ MeV (upper right panel), $20$ MeV (lower left panel), and $60$ MeV
(lower right panel). The orange (blue) hatched regions show the results using
SS (DD) data only, while the green shaded regions show the results when using
the combination of both types of experiments. For clarity, we depict as a
yellow star the true values of the analysed benchmark points. The shaded black
(grey) regions are excluded by current COHERENT data. Regarding the
uncertainty in the 8B solar neutrino flux, the solid blue line corresponds to
$\sigma_{{}^{8}\mathrm{B}}=1\%$, and the dashed blue line to
$\sigma_{{}^{8}\mathrm{B}}=4\%$.

In Fig. 6, we present the results for the same benchmark points as in Fig. 3,
but now considering the information that DD experiments can contribute. The
blue-shaded areas correspond to the best-fit regions when only DD data are
considered, while green-shaded regions are those that employ the combination
of DD and SS data. Only BP1a is observable by a future multi-ton xenon
experiment. While the corresponding mass of BP1a cannot be determined using DD
alone, the inclusion of DD data leads to a more stringent upper bound on
$m_{4}$. For BP2a, BP3a, and BP5a, DD can only set upper bounds on the mixing
parameters; however, this can still prove to be extremely useful. For example,
when combined with SS results, this can help to exclude regions with small
$m_{4}$. In the case of BP2a, for instance, DD complements the results of SS
and is crucial to better measure the sterile neutrino mass. For BP5a, DD data
improves the exclusion for small values of $m_{4}$.

Figure 7: Profile-likelihood results for BP4a ($m_{4}=40$ MeV, $\left|U_{\mu
4}\right|^{2}=9\times 10^{-3}$, $\left|U_{\tau
4}\right|^{2}=\left|U_{e4}\right|^{2}=0$). The orange (blue) hatched regions
show the results of using SS (DD) data only, while the green-shaded regions
show the results when using the combination of both types of experiment. For
clarity, we depict as a yellow star the true values of the analysed benchmark
point. The shaded black (grey) regions areas are excluded by current COHERENT
data. The yellow arrows depict that the benchmark point is outside of the
plotting area.

A particularly interesting case is that of BP4a. As explained in Section III,
for $m_{4}\simeq$40\text{\,}\mathrm{MeV}$$, the parameter reconstruction using
only data from SS experiments displays a degeneracy in the sterile neutrino
mixings and mass (see Fig. 4). In Fig. 7, we show how this degeneracy is
partially lifted when DD data is included. Although BP4a is not observable in
a future xenon detector because of its large mass, the bounds from DD exclude
the region of the parameter space with small $m_{4}$ and large
$\left|U_{e4}\right|^{2}$, which in turn leads to a good measurement of the
sterile neutrino mass.

Figure 8:  The same as in Fig. 7 but for BP1d.

Another great advantage of combining both types of experiments is that the
solar neutrino flux includes a $\nu_{\tau}$ component due to neutrino
oscillations. This provides an extra handle with which to measure the sterile
neutrino mixing with tau neutrinos. In order to test this, Fig. 8 shows an
analysis of BP1d: a benchmark point with a non-negligible $U_{\tau 4}$ mixing.
Not only is this component measured with DD data, but also the combination
with SS results leads to a better upper bound on the sterile neutrino mass and
an improved reconstruction of $U_{\tau 4}$.

For completeness, Fig. 9 shows a series of examples where both $U_{\mu 4}$ and
$U_{\tau 4}$ are non-vanishing, corresponding to BP2b, BP2c, and BP2d in Table
1. These benchmark points are observable in DD thanks to the $U_{\tau 4}$
component. When the best-fit regions are determined, the upper bound on
$|U_{\mu 4}|^{2}$ from DD data is sensitive to the magnitude of the mixing
with tau neutrinos: for small $|U_{\mu 4}|^{2}$ (e.g., BP2b), the bound on
$|U_{\mu 4}|^{2}$ is less stringent than when $|U_{\mu 4}|^{2}$ increases
(e.g., BP2d). This also makes the combination with SS results less trivial—in
some cases, the excluded regions allow for a better reconstruction of the
sterile neutrino mass (BP2b), whereas in other cases this is not possible
(BP2c and BP2d).

Figure 9:  Profile-likelihood results for benchmark points with $m_{4}=2$ MeV,
$\left|U_{\mu 4}\right|^{2}=9\times 10^{-3}$ and $\left|U_{\tau
4}\right|^{2}=9\times 10^{-4}~{}(\textit{left});\ 4\times
10^{-3}~{}(\textit{middle});\ 9\times 10^{-3}~{}(\textit{right})$. The orange
(blue) hatched regions show the results using SS (DD) data only, while the
green shaded regions show the results when using the combination of both types
of experiments. For clarity, we depict as a yellow star the true values of the
analysed benchmark points. The shaded black (grey) regions are excluded by
current COHERENT data.

###  V.1 How well can we measure the sterile neutrino mass?

As we have demonstrated, the combination of DD data with that from SS
experiments can lead to a better measurement of the sterile neutrino mass.
This can happen even in the cases where DD would not observe a new physics
signal, simply from the effect that the DD exclusions have on the regions of
the parameter space that are consistent with detection in SS experiments.
Reconstructing $m_{4}$ (i.e., confirming that it is non-vanishing) is crucial
to discriminate a sterile neutrino model from other kinds of BSM neutrino
physics (such as NSI on the active neutrinos).

In order to better quantify the relevance of the DD and SS complementary role
in measuring $m_{4}$ and to provide a more general picture, we show in Fig. 10
various projections of the $(m_{4},\,|U_{e4}|^{2},\,|U_{\mu 4}|^{2},\,|U_{\tau
4}|^{2})$ parameter space, indicating the areas where $m_{4}$ can be
reconstructed (i.e., $m_{4}=0$ is not within the $95\%$ CL region). Using the
same colour convention as in previous plots, the orange (blue) areas are those
where $m_{4}$ can be reconstructed solely from SS (DD) data, and green regions
correspond to their combination. From top to bottom, the first row corresponds
to the $(m_{4},\,|U_{\tau 4}|^{2})$ plane with $|U_{e4}|^{2}=0$ and $|U_{\mu
4}|^{2}=4\times 10^{-3}~{}(9\times 10^{-3})$ left (right) column. The second
row shows the $(m_{4},\,|U_{\mu 4}|^{2})$ plane with $|U_{e4}|^{2}=0$ and
$|U_{\tau 4}|^{2}=4\times 10^{-3}~{}(9\times 10^{-3})$ left (right) column. In
the third row, we represent the $(m_{4},\,|U_{e4}|^{2})$ plane for $|U_{\tau
4}|^{2}=0$ and $|U_{\mu 4}|^{2}=4\times 10^{-3}~{}(9\times 10^{-3})$ left
(right) column. The different benchmark points of Table 1 are indicted with
yellow stars.

In all of these figures, we observe a clear synergy between DD and SS
experiments. This is evinced by the green areas extending beyond the union of
the blue and orange ones. In particular, the addition of DD data allows us to
measure smaller values of $m_{4}$. The gap in the orange area of the top right
and lower right panels appears for $m_{4}\simeq 40$ MeV and corresponds to the
regions where the degeneracy between $|U_{e4}|^{2}$ and $|U_{\mu 4}|^{2}$
makes the mass reconstruction impossible for SS experiments alone (see Fig. 4
for BP4a). The addition of DD information is crucial to break this degeneracy
and, hence, allow for a mass reconstruction in this region (as in Fig. 7).

As already mentioned, the performance of DD experiments is extremely sensitive
to the uncertainty in the solar neutrino fluxes. For completeness, in Fig. 10
we show in dashed, dashed-dotted and dotted green lines the results obtained
when combining both types of experiments and considering a ${}^{8}\mathrm{B}$
flux uncertainty of $4\%$, $6\%$ and $12\%$, respectively. As expected, we see
how our results worsen when increasing this uncertainty.

Figure 10: Regions of the parameter space in which we have a $2\sigma$ mass
reconstruction for different solar flux uncertainties. The orange (blue)
hatched regions show the results using SS (DD) data only, while the green
shaded regions show the results when using the combination of both types of
experiments. Upper panel: Results for benchmark points with fix
$\left|U_{e4}\right|^{2}=0$ and $\left|U_{\mu 4}\right|^{2}=4\times
10^{-3}(9\times 10^{-3})$ in left and right column respectively. Middle panel:
Results for benchmark points with fix $\left|U_{e4}\right|^{2}=0$ and
$\left|U_{\tau 4}\right|^{2}=4\times 10^{-3}(9\times 10^{-3})$ in left and
right column respectively. Lower panel: Results for benchmark points with fix
$\left|U_{\tau 4}\right|^{2}=0$ and $\left|U_{\mu 4}\right|^{2}=4\times
10^{-3}(9\times 10^{-3})$ in left and right column respectively. The yellow
arrows depict that the benchmark point is outside of the plotting area.

##  VI Conclusions

In this work, we have analysed the complementarity of direct detection and
spallation source experiments for the study of sterile neutrino physics.
Specifically, we have focused on the sterile baryonic neutrino (SBN) model: an
extension of the SM that incorporates a new gauge boson that couples to
baryons and a sterile neutrino that mixes with the active ones and also
couples to this mediator. Due to this mixing, the sterile neutrino can be
produced through the up-scattering of an active neutrino with the nucleus of a
target material. This inelastic process alters the expected nuclear recoil
spectra for both DD and SS experiments, providing a characteristic signature
that can allow for the measurement of the sterile neutrino mass and mixing
parameters in the event of a future detection.

Using current data from the COHERENT collaboration on CsI and LAr, we have
first derived new constraints on the SBN model, showing that they do not
exclude new areas of the parameter space. Assuming a future SS experiment with
the projected properties of a detector to be installed at the ESS, we have
then assessed how well the sterile neutrino properties would be determined
upon a positive observation. We have shown that the new inelastic contribution
to neutrino-nucleus scattering induces a bump in the nuclear recoil spectrum.
This proves extremely useful to reconstruct the sterile neutrino mass,
conclusively disentangling this model from a generic NSI contribution to the
active neutrinos. We have demonstrated that using only SS data, values in the
range $15-50$ MeV can be measured. However, in a narrow range of masses of the
order of 40 MeV, there is a degeneracy in the measurement of the sterile
neutrino mixing that substantially affects mass reconstruction.

Incorporating future DD data helps in two ways. These detectors have an
excellent energy resolution and generally a lower energy threshold than SS
experiments. Furthermore, DD experiments are sensitive to all three neutrino
flavours, including tau neutrinos, present in the solar neutrino flux. Thus,
they are extremely helpful in removing degenerate solutions in the neutrino
mixing parameter space. Considering the case of a future multi-ton liquid
xenon experiment, we have demonstrated that the combination of future DD and
SS results is crucial to substantially increase the area of the parameter
space where the sterile neutrino mass can be reconstructed (see Fig. 10),
allowing us to measure values as low as $\sim 8$ MeV.

These results strengthen the role of DD experiments as probes of the neutrino
sector and their complementarity with dedicated neutrino detectors.

## acknowledgements

We would like to thank Pilar Coloma, Manuel González-López, Elías López
Asamar, Patrick Foldenauer, Marina Cermeño, Andrés Pérez and Karen Macías for
useful discussions and comments. DAG, DGC and MdlR acknowledge support from
the Comunidad Autonoma de Madrid and Universidad Autonoma de Madrid under
grant SI2/PBG/2020-00005, and by the Spanish Agencia Estatal de Investigación
through the grants PID2021-125331NB-I00 and CEX2020-001007-S, funded by
MCIN/AEI/10.13039/501100011033. DGC also acknowledges support from the Spanish
Ministerio de Ciencia e Innovación under grant CNS2022-135702. DA is supported
by the National Science Foundation under award 2209444.

##  Appendix A Statistical Treatment

In all of our analyses, we consider the profiled log-likelihood-ratio test
statistic, defined as

$q(\bm{\theta;\,\bm{\zeta}_{0}})\equiv-2\ln\left[\frac{\mathcal{L}(\bm{\theta},\,\hat{\bm{\omega}},\,\hat{a};\,\bm{\zeta}_{0})}{\mathcal{L}(\hat{\hat{\bm{\theta}}},\,\hat{\hat{\bm{\omega}}},\,\hat{\hat{a}};\,\bm{\zeta}_{0})}\right]\,,$
(16)

where $\mathcal{L}$ is the likelihood function describing our data given the
model parameters. For later convenience, we have split our model parameters
into three subsets, represented by $\bm{\theta}$, $\bm{\omega}$, and
$\bm{\zeta}_{0}$. The parameters $\bm{\theta}\equiv(m_{4},\,|U_{\alpha
4}|^{2})^{\mathrm{T}}$, for some given flavour index
$\alpha\in\\{e,\,\mu,\,\tau\\}$, are the two parameters we are constraining at
any given time. The parameters $\bm{\omega}\equiv(|U_{\beta
4}|^{2},\,|U_{\gamma 4}|^{2})^{\mathrm{T}}$, with $\alpha\neq\beta\neq\gamma$,
are the two remaining mixings we profile over at a given BP. Finally, as
explained in Section II, we fix the parameters related to the new vector
mediator, denoted by
$\bm{\zeta}_{0}\equiv(g_{Z^{\prime}},\,m_{Z^{\prime}})^{\mathrm{T}}$. We also
introduce a dimensionless pull parameter, $a$, as a nuisance parameter that is
designed to capture systematic uncertainties in the theoretically expected
count. We model this parameter as being Gaussian distributed with a mean of
zero and an experiment-dependent standard deviation. Hatted variables indicate
quantities that maximise the likelihood at a given parameter space point (the
null hypothesis likelihood), while double-hatted variables represent the
quantities that maximise the unconstrained likelihood (that of the alternative
hypothesis).

###  A.1 Spallation Source Experiments

Following Refs. (COHERENT:2017ipa, ; Miranda:2020tif, ; Amaral:2021rzw, ) for
SS experiments, we perform a binned statistical analysis, modelling the
likelihood of each bin $i$ as a Gaussian. In this case, Eq. 16 reduces to the
simpler $\Delta\chi^{2}$ statistic, with

$\chi^{2}(\bm{\theta},\,\bm{\omega},\,a;\,\bm{\zeta}_{0})=\sum_{i=1}^{N_{\mathrm{bins}}}\left(\frac{N^{i}_{\mathrm{obs}}-\left[1+a\right]N_{\mathrm{th}}^{i}(\bm{\theta},\,\bm{\omega};\,\bm{\zeta}_{0})}{\sigma_{\mathrm{stat}}^{i}}\right)^{2}+\left(\frac{a}{\sigma_{\mathrm{sys}}}\right)^{2}\,.$
(17)

Here, $N^{i}_{\mathrm{obs}}$ and
$N_{\mathrm{th}}^{i}(\bm{\theta},\,\bm{\omega};\,\bm{\zeta}_{0})$ are the
numbers of observed and theoretically expected events in the $i^{\mathrm{th}}$
bin, respectively. The quantity $\sigma_{\mathrm{stat}}^{i}$ is the
statistical uncertainty of the observed number of events, which we take to be

$\sigma_{\mathrm{stat}}^{i}\equiv\sqrt{N_{\mathrm{obs}}^{i}+N_{\mathrm{bkg}}^{i}}\,,$
(18)

where $N_{\mathrm{bkg}}^{i}$ is the expected number of background events in
the $i^{\mathrm{th}}$ bin. When performing our analysis of COHERENT data, we
use the backgrounds reported by collaboration (COHERENT:2020iec, ;
COHERENT:2021xmm, ). However, when considering the future ESS experiment, we
instead use the fact that the beam-related neutron (BRN) background represents
an important background in this type of search, with CENNS-10 reporting that
$10\%$ of its measured signal events arose due to this background source
(COHERENT:2019kwz, ). Since we make no assumptions on how well future SS
experiments will handle this background, we take $N_{\mathrm{bkg}}^{i}\equiv
N_{\mathrm{SM}}^{i}/10$, with $N_{\mathrm{SM}}^{i}$ the number of expected
CE$\nu$NS events in the $i^{\mathrm{th}}$ bin under the SM. For the pull
parameter, $a$, we take its uncertainty to be $\sigma_{\mathrm{sys}}=0.05$
(Miranda:2020syh, ; Amaral:2021rzw, ).

To construct the $\Delta\chi^{2}$ for our parameters of interest, we compute
the profiled test statistic

$\Delta\chi^{2}(\bm{\theta};\,\bm{\zeta}_{0})=\chi^{2}(\bm{\theta},\,\hat{\bm{\omega}},\,\hat{a};\,\bm{\zeta}_{0})-\chi^{2}(\hat{\hat{\bm{\theta}}},\,\hat{\hat{\bm{\omega}}},\,\hat{\hat{a}};\,\bm{\zeta}_{0})\,.$
(19)

As explained in Section III, we make use of Asimov data sets throughout our
analyses. This means that our ‘observed’ data are set to the theoretically
expected number of events for each given benchmark point. This leads to two
simplifications. Firstly, as the data are perfectly consistent with a given
BP, we know that the value of the overall minimised $\chi^{2}$ will be zero.
Secondly, the minimisation over $a$ can be done without resorting to numerical
methods for any given $\bm{\theta}$ and $\hat{\bm{\omega}}$. By simply finding
that value of $a$ for which $\partial_{a}(\Delta\chi^{2})=0$, we get the
analytical result

$\hat{a}=\left[\sum_{i}\frac{\left(N^{i}_{\text{obs}}-N^{i}_{\text{th}}\right)N^{i}_{\text{th}}}{\left(\sigma^{i}_{\text{stat}}\right)^{2}}\right]\Bigg{/}\left[\left(\sigma_{\text{sys}}\right)^{-2}+\sum_{i}\left(\frac{N^{i}_{\text{th}}}{\sigma^{i}_{\text{stat}}}\right)^{2}\right]\,.$
(20)

Note that, since $N_{\mathrm{th}}^{i}$ is not a function of $a$, the
minimisation over $a$ and $\bm{\omega}$ can be done separately.

Finally, when drawing our contours for the $95\%$ CL regions, we use the fact
that our $\Delta\chi^{2}$ should be distributed according to a $\chi^{2}$
distribution with $2$ degrees of freedom. This is because, of the $7$
parameters that Eq. 17 depends on, we profile over $3$ of them in Eq. 19,
keeping the remaining $2$, represented by $\bm{\zeta}_{0}$, fixed throughout.
We therefore draw the boundaries of our regions at $\Delta\chi^{2}=6.18$.

###  A.2 Direct Detection Experiments

For DD experiments, we also perform a binned statistical treatment. However,
unlike for SS experiments, we assume that the number of counts in each bin
follows a Poisson distribution due to the lower number of events expected
within the high-energy bins. Inserting a Poisson likelihood for $\mathcal{L}$
in Eq. 16 and once again exploiting our use of Asimov data sets, we get that

$\begin{split}q(\bm{\theta};\,\bm{\zeta}_{0})&=2\left[\sum_{i=1}^{N_{\mathrm{bins}}}(1+\hat{a})N_{\mathrm{th}}^{i}(\bm{\theta,\,\hat{\omega};\,\bm{\zeta}_{0}})-N_{\mathrm{obs}}^{i}+N_{\mathrm{obs}}^{i}\ln\frac{N_{\mathrm{obs}}^{i}}{(1+\hat{a})N_{\mathrm{th}}^{i}(\bm{\theta,\,\hat{\omega};\,\bm{\zeta}_{0}})}\right]+\left(\frac{\hat{a}}{\sigma_{\mathrm{{}^{8}B}}}\right)^{2}\,.\end{split}$
(21)

Note that, as for SS experiments, we have also introduced the pull parameter
$a$ to capture the effect of systematic uncertainties. In the case of DD
experiments searching for CE$\nu$NS, we assume that this is dominated by the
uncertainty in the $\mathrm{{}^{8}B}$ solar neutrino flux,
$\sigma_{\mathrm{{}^{8}B}}$, for which we take different values in the main
text.

As before, we can derive the analytical form for $\hat{a}$; we do this by
solving the equation $\partial_{a}q=0$. We find that

$\hat{a}=\frac{-(1+N_{\mathrm{th}}^{\mathrm{tot}}\sigma_{\mathrm{{}^{8}B}}^{2})+\sqrt{(1+N_{\mathrm{th}}^{\mathrm{tot}}\sigma_{\mathrm{{}^{8}B}}^{2})^{2}-4\sigma_{\mathrm{{}^{8}B}}^{2}(N_{\mathrm{th}}^{\mathrm{tot}}-N_{\mathrm{obs}}^{\mathrm{tot}})}}{2},\,$
(22)

where $N_{\mathrm{obs}}^{\mathrm{tot}}$ and $N_{\mathrm{th}}^{\mathrm{tot}}$
are the total observed and theoretically expected number of events across all
bins, respectively. We note that in Eqs. 21 and 22 we have neglected any
background contribution, as the background ($\mathcal{O}(1)$) in DARWIN is
expected to be much smaller than the expected signal ($\mathcal{O}(10^{2-3})$)
for the majority of bins. Since the pull parameter $a$ only impacts the
signal, the analytical minimisation presented in Eq. 22 is only possible with
zero (or, more generally, constant) background. With a bin-variable background
contribution, the minimisation must instead be done numerically.

To draw our $95\%$ CL limits, we make use of Wilks’ theorem Cowan:2010js .
This tells us that the log-likelihood-ratio test statistic asymptotically
follows a $\chi^{2}$ distribution with number of degrees of freedom equal to
the difference in the number of free parameters between the null and
alternative hypotheses. As previously, this gives us two degrees of freedom.
We therefore draw the boundaries of our regions at $q=6.18$.

## References

  * (1) D. Britton et al., _Improved search for massive neutrinos in $\pi^{+}\to e^{+}\nu$ decay_, _Phys. Rev. D_ 46 (1992) 885–887.
  * (2) D. Britton et al., _Measurement of the $\pi^{+}\to e^{+}$ neutrino branching ratio_, _Phys. Rev. D_ 49 (1994) 28–39.
  * (3) PIENU collaboration, A. Aguilar-Arevalo et al., _Search for heavy neutrinos in $\pi\to\mu\nu$ decay_, _Phys. Lett. B_ 798 (2019) 134980, [1904.03269].
  * (4) NA62 collaboration, E. Cortina Gil et al., _Search for heavy neutral lepton production in $K^{+}$ decays to positrons_, 2005.09575.
  * (5) T2K collaboration, K. Abe et al., _Search for heavy neutrinos with the T2K near detector ND280_ , _Phys. Rev. D_ 100 (2019) 052006, [1902.07598].
  * (6) CHARM collaboration, F. Bergsma et al., _A Search for Decays of Heavy Neutrinos in the Mass Range 0.5-GeV to 2.8-GeV_ , _Phys. Lett. B_ 166 (1986) 473–478.
  * (7) WA66 collaboration, A. M. Cooper-Sarkar et al., _Search for Heavy Neutrino Decays in the BEBC Beam Dump Experiment_ , _Phys. Lett. B_ 160 (1985) 207–211.
  * (8) DELPHI collaboration, P. Abreu et al., _Search for neutral heavy leptons produced in Z decays_ , _Z. Phys. C_ 74 (1997) 57–71. [Erratum: Z.Phys.C 75, 580 (1997)].
  * (9) NuTeV, E815 collaboration, A. Vaitaitis et al., _Search for neutral heavy leptons in a high-energy neutrino beam_ , _Phys. Rev. Lett._ 83 (1999) 4943–4946, [hep-ex/9908011].
  * (10) K. Goldhagen, M. Maltoni, S. E. Reichard and T. Schwetz, _Testing sterile neutrino mixing with present and future solar neutrino data_ , _Eur. Phys. J. C_ 82 (2022) 116, [2109.14898].
  * (11) M. Dentler, A. Hernández-Cabezudo, J. Kopp, P. A. N. Machado, M. Maltoni, I. Martinez-Soler et al., _Updated Global Analysis of Neutrino Oscillations in the Presence of eV-Scale Sterile Neutrinos_ , _JHEP_ 08 (2018) 010, [1803.10661].
  * (12) D. V. Forero, C. Giunti, C. A. Ternes and M. Tortola, _Nonunitary neutrino mixing in short and long-baseline experiments_ , _Phys. Rev. D_ 104 (2021) 075030, [2103.01998].
  * (13) MINOS+ collaboration, P. Adamson et al., _Search for sterile neutrinos in MINOS and MINOS+ using a two-detector fit_ , _Phys. Rev. Lett._ 122 (2019) 091803, [1710.06488].
  * (14) K. N. Abazajian et al., _Light Sterile Neutrinos: A White Paper_ , 1204.5379.
  * (15) N. Sabti, A. Magalich and A. Filimonova, _An Extended Analysis of Heavy Neutral Leptons during Big Bang Nucleosynthesis_ , _JCAP_ 11 (2020) 056, [2006.07387].
  * (16) A. M. Abdullahi et al., _The present and future status of heavy neutral leptons_ , _J. Phys. G_ 50 (2023) 020501, [2203.08039].
  * (17) A. C. Vincent, E. F. Martinez, P. Hernández, M. Lattanzi and O. Mena, _Revisiting cosmological bounds on sterile neutrinos_ , _JCAP_ 04 (2015) 006, [1408.1956].
  * (18) A. Boyarsky, M. Ovchynnikov, N. Sabti and V. Syvolap, _When feebly interacting massive particles decay into neutrinos: The Neff story_ , _Phys. Rev. D_ 104 (2021) 035006, [2103.09831].
  * (19) COHERENT collaboration, D. Akimov et al., _The COHERENT Experiment at the Spallation Neutron Source_ , 1509.08702.
  * (20) COHERENT collaboration, D. Akimov et al., _Observation of Coherent Elastic Neutrino-Nucleus Scattering_ , _Science_ 357 (2017) 1123–1126, [1708.01294].
  * (21) COHERENT collaboration, D. Akimov et al., _First Measurement of Coherent Elastic Neutrino-Nucleus Scattering on Argon_ , _Phys. Rev. Lett._ 126 (2021) 012002, [2003.10630].
  * (22) D. Z. Freedman, _Coherent Neutrino Nucleus Scattering as a Probe of the Weak Neutral Current_ , _Phys. Rev. D_ 9 (1974) 1389–1392.
  * (23) A. Drukier and L. Stodolsky, _Principles and Applications of a Neutral Current Detector for Neutrino Physics and Astronomy_ , _Phys. Rev. D_ 30 (1984) 2295.
  * (24) T. S. Kosmas, D. K. Papoulias, M. Tortola and J. W. F. Valle, _Probing light sterile neutrino signatures at reactor and Spallation Neutron Source neutrino experiments_ , _Phys. Rev. D_ 96 (2017) 063013, [1703.00054].
  * (25) O. G. Miranda, D. K. Papoulias, O. Sanders, M. Tórtola and J. W. F. Valle, _Future CEvNS experiments as probes of lepton unitarity and light-sterile neutrinos_ , _Phys. Rev. D_ 102 (2020) 113014, [2008.02759].
  * (26) M. Abdullah, J. B. Dent, B. Dutta, G. L. Kane, S. Liao and L. E. Strigari, _Coherent elastic neutrino nucleus scattering as a probe of a Z’ through kinetic and mass mixing effects_ , _Phys. Rev. D_ 98 (2018) 015005, [1803.01224].
  * (27) M. Bauer, P. Foldenauer and J. Jaeckel, _Hunting All the Hidden Photons_ , _JHEP_ 07 (2018) 094, [1803.05466].
  * (28) O. G. Miranda, D. K. Papoulias, M. Tórtola and J. W. F. Valle, _Probing new neutral gauge bosons with $CE\nu NS$ and neutrino-electron scattering_, _Phys. Rev. D_ 101 (2020) 073005, [2002.01482].
  * (29) M. Bauer, P. Foldenauer and M. Mosny, _Flavor structure of anomaly-free hidden photon models_ , _Phys. Rev. D_ 103 (2021) 075024, [2011.12973].
  * (30) COHERENT collaboration, D. Akimov et al., _First Measurement of Coherent Elastic Neutrino-Nucleus Scattering on Argon_ , _Phys. Rev. Lett._ 126 (2021) 012002, [2003.10630].
  * (31) COHERENT collaboration, D. Akimov et al., _Sensitivity of the COHERENT Experiment to Accelerator-Produced Dark Matter_ , _Phys. Rev. D_ 102 (2020) 052007, [1911.06422].
  * (32) “Coherent captain mills.” https://p25ext.lanl.gov/%7Elee/CaptainMills/.
  * (33) D. Baxter et al., _Coherent Elastic Neutrino-Nucleus Scattering at the European Spallation Source_ , _JHEP_ 02 (2020) 123, [1911.00762].
  * (34) O. G. Miranda, D. K. Papoulias, G. Sanchez Garcia, O. Sanders, M. Tórtola and J. W. F. Valle, _Implications of the first detection of coherent elastic neutrino-nucleus scattering (CEvNS) with Liquid Argon_ , _JHEP_ 05 (2020) 130, [2003.12050]. [Erratum: JHEP 01, 067 (2021)].
  * (35) H. Banerjee, P. Byakti and S. Roy, _Supersymmetric gauged U(1) ${}_{L_{\mu}-L_{\tau}}$ model for neutrinos and the muon $(g-2)$ anomaly_, _Phys. Rev. D_ 98 (2018) 075022, [1805.04415].
  * (36) D. K. Papoulias, _COHERENT constraints after the COHERENT-2020 quenching factor measurement_ , _Phys. Rev. D_ 102 (2020) 113004, [1907.11644].
  * (37) A. N. Khan and W. Rodejohann, _New physics from COHERENT data with an improved quenching factor_ , _Phys. Rev. D_ 100 (2019) 113003, [1907.12444].
  * (38) P. D. Bolton, F. F. Deppisch, K. Fridell, J. Harz, C. Hati and S. Kulkarni, _Probing active-sterile neutrino transition magnetic moments with photon emission from CE $\nu$NS_, _Phys. Rev. D_ 106 (2022) 035036, [2110.02233].
  * (39) O. G. Miranda, D. K. Papoulias, O. Sanders, M. Tórtola and J. W. F. Valle, _Low-energy probes of sterile neutrino transition magnetic moments_ , _JHEP_ 12 (2021) 191, [2109.09545].
  * (40) V. De Romeri, O. G. Miranda, D. K. Papoulias, G. Sanchez Garcia, M. Tórtola and J. W. F. Valle, _Physics implications of a combined analysis of COHERENT CsI and LAr data_ , 2211.11905.
  * (41) P. M. Candela, V. De Romeri and D. K. Papoulias, _COHERENT production of a Dark Fermion_ , 2305.03341.
  * (42) D. G. Cerdeño, M. Fairbairn, T. Jubb, P. A. N. Machado, A. C. Vincent and C. Bœhm, _Physics from solar neutrinos in dark matter direct detection experiments_ , _JHEP_ 05 (2016) 118, [1604.01025]. [Erratum: JHEP09,048(2016)].
  * (43) B. Dutta, S. Liao, L. E. Strigari and J. W. Walker, _Non-standard interactions of solar neutrinos in dark matter experiments_ , _Phys. Lett. B_ 773 (2017) 242–246, [1705.00661].
  * (44) G. B. Gelmini, V. Takhistov and S. J. Witte, _Geoneutrinos in Large Direct Detection Experiments_ , _Phys. Rev. D_ 99 (2019) 093009, [1812.05550].
  * (45) R. Essig, M. Sholapurkar and T.-T. Yu, _Solar Neutrinos as a Signal and Background in Direct-Detection Experiments Searching for Sub-GeV Dark Matter With Electron Recoils_ , _Phys. Rev. D_ 97 (2018) 095029, [1801.10159].
  * (46) D. W. P. d. Amaral, D. G. Cerdeno, P. Foldenauer and E. Reid, _Solar neutrino probes of the muon anomalous magnetic moment in the gauged $\mathrm{U}{(1)}_{L_{\mu}-{L}_{\tau}}$_, _JHEP_ 12 (2020) 155, [2006.11225].
  * (47) B. Dutta, R. F. Lang, S. Liao, S. Sinha, L. Strigari and A. Thompson, _A global analysis strategy to resolve neutrino NSI degeneracies with scattering and oscillation data_ , _JHEP_ 09 (2020) 106, [2002.03066].
  * (48) D. W. P. Amaral, D. G. Cerdeno, A. Cheek and P. Foldenauer, _Confirming $U(1)_{L_{\mu}-L_{\tau}}$ as a solution for $(g-2)_{\mu}$ with neutrinos_, _Eur. Phys. J. C_ 81 (2021) 861, [2104.03297].
  * (49) V. Munoz, V. Takhistov, S. J. Witte and G. M. Fuller, _Exploring the origin of supermassive black holes with coherent neutrino scattering_ , _JCAP_ 11 (2021) 020, [2102.00885].
  * (50) A. de Gouvêa, E. McGinness, I. Martinez-Soler and Y. F. Perez-Gonzalez, _pp solar neutrinos at DARWIN_ , _Phys. Rev. D_ 106 (2022) 096017, [2111.02421].
  * (51) D. W. P. Amaral, D. Cerdeno, A. Cheek and P. Foldenauer, _A direct detection view of the neutrino NSI landscape_ , 2302.12846.
  * (52) I. M. Shoemaker and J. Wyenberg, _Direct Detection Experiments at the Neutrino Dipole Portal Frontier_ , _Phys. Rev. D_ 99 (2019) 075010, [1811.12435].
  * (53) XENON collaboration, E. Aprile et al., _Excess electronic recoil events in XENON1T_ , _Phys. Rev. D_ 102 (2020) 072004, [2006.09721].
  * (54) I. M. Shoemaker, Y.-D. Tsai and J. Wyenberg, _Active-to-sterile neutrino dipole portal and the XENON1T excess_ , _Phys. Rev. D_ 104 (2021) 115026, [2007.05513].
  * (55) V. Brdar, A. Greljo, J. Kopp and T. Opferkuch, _The Neutrino Magnetic Moment Portal: Cosmology, Astrophysics, and Direct Detection_ , _JCAP_ 01 (2021) 039, [2007.15563].
  * (56) XENON collaboration, E. Aprile et al., _Search for New Physics in Electronic Recoil Data from XENONnT_ , _Phys. Rev. Lett._ 129 (2022) 161805, [2207.11330].
  * (57) M. Pospelov, _Neutrino Physics with Dark Matter Experiments and the Signature of New Baryonic Neutral Currents_ , _Phys. Rev. D_ 84 (2011) 085008, [1103.3261].
  * (58) C. A. Argüelles, N. Foppiani and M. Hostert, _Efficiently exploring multidimensional parameter spaces beyond the Standard Model_ , _Phys. Rev. D_ 107 (2023) 035027, [2205.12273].
  * (59) R. H. Helm, _Inelastic and Elastic Scattering of 187-Mev Electrons from Selected Even-Even Nuclei_ , _Phys. Rev._ 104 (1956) 1466–1475.
  * (60) J. D. Lewin and P. F. Smith, _Review of mathematics, numerical factors, and corrections for dark matter experiments based on elastic nuclear recoil_ , _Astropart. Phys._ 6 (1996) 87–112.
  * (61) A. L. Foguel, P. Reimitz and R. Z. Funchal, _A robust description of hadronic decays in light vector mediator models_ , _JHEP_ 04 (2022) 119, [2201.01788].
  * (62) P. Coloma, P. B. Denton, M. C. Gonzalez-Garcia, M. Maltoni and T. Schwetz, _Curtailing the Dark Side in Non-Standard Neutrino Interactions_ , _JHEP_ 04 (2017) 116, [1701.04828].
  * (63) L. Wolfenstein, _Neutrino Oscillations in Matter_ , _Phys. Rev. D_ 17 (1978) 2369–2374.
  * (64) M. M. Guzzo and S. T. Petcov, _On the matter enhanced transitions of solar neutrinos in the absence of neutrino mixing in vacuum_ , _Phys. Lett. B_ 271 (1991) 172–178.
  * (65) M. M. Guzzo, A. Masiero and S. T. Petcov, _On the MSW effect with massless neutrinos and no mixing in the vacuum_ , _Phys. Lett. B_ 260 (1991) 154–160.
  * (66) M. C. Gonzalez-Garcia, M. M. Guzzo, P. I. Krastev, H. Nunokawa, O. L. G. Peres, V. Pleitez et al., _Atmospheric neutrino observations and flavor changing interactions_ , _Phys. Rev. Lett._ 82 (1999) 3202–3205, [hep-ph/9809531].
  * (67) S. Bergmann, M. M. Guzzo, P. C. de Holanda, P. I. Krastev and H. Nunokawa, _Status of the solution to the solar neutrino problem based on nonstandard neutrino interactions_ , _Phys. Rev. D_ 62 (2000) 073001, [hep-ph/0004049].
  * (68) M. M. Guzzo, H. Nunokawa, P. C. de Holanda and O. L. G. Peres, _On the massless ’just-so’ solution to the solar neutrino problem_ , _Phys. Rev. D_ 64 (2001) 097301, [hep-ph/0012089].
  * (69) M. Guzzo, P. C. de Holanda, M. Maltoni, H. Nunokawa, M. A. Tortola and J. W. F. Valle, _Status of a hybrid three neutrino interpretation of neutrino data_ , _Nucl. Phys. B_ 629 (2002) 479–490, [hep-ph/0112310].
  * (70) M. C. Gonzalez-Garcia and M. Maltoni, _Atmospheric neutrino oscillations and new physics_ , _Phys. Rev. D_ 70 (2004) 033010, [hep-ph/0404085].
  * (71) C. A. Argüelles, G. Barenboim, M. Bustamante, P. Coloma, P. B. Denton, I. Esteban et al., _Snowmass white paper: beyond the standard model effects on neutrino flavor_ , _The European Physical Journal C_ 83 (jan, 2023) .
  * (72) G. Cowan, K. Cranmer, E. Gross and O. Vitells, _Asymptotic formulae for likelihood-based tests of new physics_ , _Eur. Phys. J. C_ 71 (2011) 1554, [1007.1727]. [Erratum: Eur.Phys.J.C 73, 2501 (2013)].
  * (73) D. W. P. Amaral, D. G. Cerdeño, A. Cheek and P. Foldenauer, _SNuDD [Computer Software] , available at https://github.com/snudd/snudd.git_, 2023\.
  * (74) F. Feroz, M. P. Hobson and M. Bridges, _MULTINEST: an efficient and robust Bayesian inference tool for cosmology and particle physics_ , _MNRAS_ 398 (Oct., 2009) 1601–1614, [0809.3437].
  * (75) F. Feroz, M. P. Hobson, E. Cameron and A. N. Pettitt, _Importance Nested Sampling and the MultiNest Algorithm_ , _The Open Journal of Astrophysics_ 2 (Nov., 2019) 10, [1306.2144].
  * (76) J. Buchner, A. Georgakakis, K. Nandra, L. Hsu, C. Rangel, M. Brightman et al., _X-ray spectral modelling of the AGN obscuring region in the CDFS: Bayesian model selection and catalogue_ , _A &A_ 564 (Apr., 2014) A125, [1402.0004].
  * (77) COHERENT collaboration, D. Akimov et al., _Measurement of the Coherent Elastic Neutrino-Nucleus Scattering Cross Section on CsI by COHERENT_ , _Phys. Rev. Lett._ 129 (2022) 081801, [2110.07730].
  * (78) LZ collaboration, D. S. Akerib et al., _Projected WIMP sensitivity of the LUX-ZEPLIN dark matter experiment_ , _Phys. Rev. D_ 101 (2020) 052002, [1802.06039].
  * (79) XENON collaboration, E. Aprile et al., _Projected WIMP sensitivity of the XENONnT dark matter experiment_ , _JCAP_ 11 (2020) 031, [2007.08796].
  * (80) DARWIN collaboration, J. Aalbers et al., _DARWIN: towards the ultimate dark matter detector_ , _JCAP_ 11 (2016) 017, [1606.07001].
  * (81) C. A. J. O’Hare, _New Definition of the Neutrino Floor for Direct Dark Matter Searches_ , _Phys. Rev. Lett._ 127 (2021) 251802, [2109.03116].
  * (82) P. Coloma, M. C. Gonzalez-Garcia, M. Maltoni, J. P. Pinheiro and S. Urrea, _Constraining new physics with Borexino Phase-II spectral data_ , _JHEP_ 07 (2022) 138, [2204.03011]. [Erratum: JHEP 11, 138 (2022)].
  * (83) N. Vinyoles, A. M. Serenelli, F. L. Villante, S. Basu, J. Bergström, M. C. Gonzalez-Garcia et al., _A new Generation of Standard Solar Models_ , _Astrophys. J._ 835 (2017) 202, [1611.09867].
  * (84) I. Esteban, M. C. Gonzalez-Garcia, M. Maltoni, T. Schwetz and A. Zhou, _The fate of hints: updated global analysis of three-flavor neutrino oscillations_ , _JHEP_ 09 (2020) 178, [2007.14792].
  * (85) M. Pospelov and J. Pradler, _Elastic scattering signals of solar neutrinos with enhanced baryonic currents_ , _Phys. Rev. D_ 85 (2012) 113016, [1203.0545]. [Erratum: Phys.Rev.D 88, 039904 (2013)].
  * (86) J. N. Bahcall, A. M. Serenelli and S. Basu, _New solar opacities, abundances, helioseismology, and neutrino fluxes_ , _Astrophys. J. Lett._ 621 (2005) L85–L88, [astro-ph/0412440].
  * (87) J. Bergstrom, M. C. Gonzalez-Garcia, M. Maltoni, C. Pena-Garay, A. M. Serenelli and N. Song, _Updated determination of the solar neutrino fluxes from solar neutrino data_ , _JHEP_ 03 (2016) 132, [1601.00972].
  * (88) F. Capozzi, S. W. Li, G. Zhu and J. F. Beacom, _DUNE as the Next-Generation Solar Neutrino Experiment_ , _Phys. Rev. Lett._ 123 (2019) 131803, [1808.08232].
  * (89) SNO collaboration, B. Aharmim et al., _Combined Analysis of all Three Phases of Solar Neutrino Data from the Sudbury Neutrino Observatory_ , _Phys. Rev. C_ 88 (2013) 025501, [1109.0763].