Neutrino non-standard interactions and dark matter searches with multi-ton scale detectors

Future dark matter (DM) direct detection searches will be subject to irreducible neutrino backgrounds that will challenge the identification of an actual WIMP signal in experiments without directionality sensitivity. We study the impact of neutrino-quark non-standard interactions (NSI) on this background, assuming the constraints from neutrino oscillations and the recent COHERENT experiment data, which are relevant for NSI mediated by light mediators, $m_{\rm med} \lesssim{\cal O}$(GeV). We calculate the expected number of neutrino-nucleus elastic scattering events in a Xe-based ton-size dark matter detector, including solar neutrino fluxes from the $pp$ chain and CNO cycle as well as sub-GeV atmospheric fluxes and taking into account NSI effects in both propagation and detection. We find that sizable deviations from the standard model expectation are possible, but are more pronounced for flavor-diagonal couplings, in particular for electron neutrinos. We show that neutrino NSI can enhance or deplete the neutrino-nucleus event rate, which may impact DM searches in multi-ton detectors.


Introduction
A wealth of astrophysical and cosmological data support the idea that the dominant form of matter in the universe is made of dark matter (DM). Several DM candidates have been proposed ranging from astrophysical observable objects [1] up to fundamental particles [2]. Among the latter, the most studied candidates (both theoretically and experimentally) are massive particles whose abundance is determined by thermal freeze-out (weakly interacting massive particles or WIMPs).
WIMPs are being searched for using direct detection experiments since the late 80's, following a proposal by Goodman and Witten [3]. These experiments rely on energy deposition measurements from WIMP-induced nuclear recoil processes, which are expected to result as the Earth (detector) moves in the galactic DM halo. The non-observation of such signatures are translated into upper limits for the WIMP-nucleon (spin-independent or spin-dependent) cross section as a function of the DM mass. At present the most stringent bound is given by the XENON1T collaboration, which in the absence of a nuclear recoil signal has derived the limit σ χ−n 7.7 × 10 −47 cm 2 /c 2 for m χ = 35 GeV/c 2 at 90% CL [4]. Searches beyond these values and down to, say, σ χ−n ∼ 10 −49 cm 2 /c 2 require multi-ton target masses, as planned for the DARWIN detector [5].
The reach of these experiments, however, will be limited by an irreducible neutrino background. Indeed, coherent neutrino nucleus elastic scattering (CνNES) induced by solar and, eventually, atmospheric neutrino fluxes will mimic a possible WIMP-nucleus signal [6]. While the precise extent of this background is affected by various factors, including neutrino flux uncertainties and astrophysical parameters [6,7], in the SM, the CνNES cross section is well understood and is not a source of uncertainty. This can potentially change in the presence of new physics. The most obvious way through which this can occur is via neutrino-quark non-standard interactions (NSI), which historically were introduced by Wolfenstein as a possible way to account for the solar neutrino deficit [8].
Neutrino NSI are constrained by several observations, including neutrino oscillation data [9], neutrino deep inelastic scattering experimental data (CHARM [10] and NuTeV [11]) and, potentially, the measurement of the CνNES process itself, recently observed by the COHERENT experiment [12] (see [13] for a review and [14] for the detailed analysis). CHARM and NuTeV set strong constraints on ν µ -related NSI parameters, but their relevance for constraining CνNES depends on the range of application of the NSI effective operators. In particular, these constraints do not apply if the NSI are mediated by a particle with mass m med < O(1) GeV [15]. In what follows we will consider the possibility that NSI are mediated by a particle in the mass range 10 MeV m med 1 GeV, range in which CHARM and NuTeV bounds do not apply, but one is still allowed to use the constraints on NSI as extracted from a combined analysis of COHERENT and oscillation data in ref. [14]. It is worth pointing out that in this mass range bounds are less severe for all flavor-diagonal NSI couplings, and so compared to the case where deep inelastic scattering data apply one expects larger deviations on the neutrino-nucleus event rates.
Specifically, using current constraints on neutrino NSI derived from neutrino oscillation and CO-HERENT data [14], in this paper we quantify the deviations on the CνNES expected number of events in a Xe-based DM detector using an exposure E = 1 (ton · year). For that aim we consider solar neutrino fluxes for which we include neutrinos from the pp chain (pp, hep, 8 B, 7 Be, pep) and CNO cycle ( 13 N, 15 O, 17 F). We consider as well sub-GeV atmospheric neutrino fluxes. For the solar sector analysis we calculate NSI propagation effects in the mass dominance limit (∆m 2 31 → ∞) and for the calculation of the electron neutrino survival probability we take into account the distribution of neutrino production in the Sun according to the standard solar model (SSM) BS05 [16]. We quantify matter-induced effects in the atmospheric sector and find they are negligible (H mat /H vac 0.1) [17], however we include neutrino flavor oscillations in vacuum in the calculation of the number of neutrino-nucleus events.
The study of beyond-the-standard model physics in neutrino-quark interactions in the context of multi-ton DM detectors, has been recently the subject of various analyses. Ref. [18] studied the capability of such detectors to constraint light mediators scenarios, while ref. [19] instead considered scenarios where the new physics couples both neutrinos and DM. Ref. [20] considered the implications of neutrino NSI on neutrino backgrounds. The analysis was done for solar neutrinos and in particular for the 8 B neutrino flux, for which solar neutrino fluxes can mimic a WIMP-nucleus scattering event. We extend upon this study by reanalyzing NSI effects in the case of lighter mediators, and so in the light of the recently released COHERENT data. We furthermore consider NSI effects all over solar neutrino energies and include sub-GeV atmospheric neutrino fluxes. The former is particularly relevant if one aims at setting e.g. background-free sensitivities for light WIMPS.
The rest of this paper is organized as follows. In sec. 2 we introduce neutrino NSI, first assuming they are generated by heavy mediators (q 2 m 2 X ) 1 and review current experimental allowed ranges. In sec. 3 we discuss CνNES in the presence of NSI, calculate electron neutrino survival probabilities in the Sun for various NSI data sets and study atmospheric vacuum conversion probabilities. In sec. 4 we present our main results, namely the maximum and minimum number of CνNES events in ton-size DM detectors (Xe based). We perform a single-parameter analysis, c'est-à-dire we assume all NSI parameters are zero but one and vary that coupling within its experimental allowed range. Although constrained, this analysis enable us to capture the main features of what should be expected from neutrino NSI in near-future DM direct detection experiments. Finally, in sec. 5 we summarize and present our conclusions. 1 They might as well proceed through the exchange of light mediators (q 2 m 2 X ). The results for this case will be presented elsewhere [21].
2 Neutrino-quark non-standard interactions NSI can be parametrized in terms of the following effective Lagrangian (see e.g. [22]): where ν T = (ν e , ν µ , ν τ ) and qV and qA refer to 3×3 matrices in flavor space whose elements correspond to vector and vector-axial couplings. Hereafter we will denote flavor-related matrices in bold-face. This parametrization assumes that the neutrino NSI are related to new physics at a scale m med q, with q the characteristic energy momentum transferred in the relevant process. It can be as well that the new physics is related with a scale way below the typical momentum exchange (light mediator scenario), in which case a phenomenological analysis requires to take into account the dynamics of the mediator. Such models with light mediators have been discussed in, for instance [23][24][25]. Here, we focus on the intermediate possibility where CHARM and NuTeV constraints do not apply, while an effective approach is still valid so that the constraints on NSI parameters from the COHERENT data analyzed in [14] hold. In this regime, the constraints on the flavor-diagonal NSI couplings are substantially weaker. For instance (see e.g. [13]) Accelerator constraints are irrelevant for mediator masses below O(1 GeV). The precise value is not important here, as long as it is larger than about 100 MeV, the maximal momentum exchange scale we will consider. Since [14] assumes m med 10MeV, our analysis will be therefore limited to 10MeV < m med < 1GeV. While a priori this seems a narrow window in parameter space, it is however for that window where the largest possible deviations-compared to SM expectations-on the neutrino-nucleus scattering rate are expected.
The interactions in (1) lead to forward coherent scattering (order G F interaction) and scattering (order G 2 F interaction) processes such as CνNES. The former are responsible for matter potentials in the Sun and in the Earth and are entirely controlled by qV . For the latter, the vector current determines spin-independent processes while the axial-vector current accounts for spin-dependent ones. Reactions rates driven by vector-axial currents are relatively smaller (in particular for large nuclei) and so we will not consider them. From now on then we drop the index V .
Constraints on NSI are derived from neutrino oscillation [9,13,26,27] and COHERENT data [12]. Let us discuss this in more detail. The presence of the vector current in (1) induce additional diagonal and non-diagonal matter potential terms which affect neutrino flavor evolution in the Sun and in the Earth (depending on the neutrino energy). They produce distortions in the neutrino oscillation probabilities that can be constrained by using solar, KamLAND, atmospheric and long-baseline neutrino data [9]. Limits on all flavor non-diagonal NSI couplings and two combinations of diagonal couplings e.g. q ii − q µµ (i = e, τ ) can be derived this way 2 . For comparison we show the 90% CL limits reported in [9] for the large mixing angle (LMA) case 3 in table 1.  [9], while those in the table to the right arise from the combined analysis of oscillation and COHERENT data [14], and are the ones we use in our analysis. In the oscillation data table¯ q ii refers to q ii − q µµ (i = e, τ ).
The vector current modifies as well the CνNES cross section, and so data from the recent measurement of the CνNES process by the COHERENT experiment [12] can be used to set bounds on the vector NSI parameters. The low-energy neutrinos (E ν 50 MeV) used by COHERENT result from accelerator-driven protons that strike a fixed mercury target thus producing pions. Their decay produce a monochromatic muon neutrino flux and muons that afterward decay, thus producing anti-muon and electron continuous neutrino fluxes. These data therefore enable placing bounds on q ei and q µi (i = e, µ, τ ) [12,28]. Overall, bounds can be placed by two independent experimental sources (order G F and order G 2 F processes), but can be as well derived by combining both, resulting in more competitive limits. Such an analysis has been done in ref. [14], and these are the limits we use in the study of the impact of neutrino NSI in near-future DM detectors in sec. 2. For reference they are displayed in table 1 (table to the right).

Coherent neutrino-nucleus elastic cross section and NSI
CνNES is a neutral-current process in which the momentum exchange is of order q 100 MeV. In the SM it is therefore well-described by the effective Lagrangian where g q V and g q A are the SM vector and axial-vector couplings. The cross section for CνNES with spinless nuclei is given by [29,30] Here m N is the nucleus mass and E r the nuclear recoil energy, whose maximum value, E max r , is given by Note that since E ν m N (E ν 100 MeV for the neutrino-nucleus scattering process to be coherent), The atomic number and number of neutrons are encoded in Q SM = N − (1 − 4 sin 2 θ W )Z, from which one can see the N 2 enhancement of the CνNES cross section. Finally, we use the Helm [31] form factor where q can be traded to recoil energy through q = 6.92 is the order-one spherical Bessel function of the first kind and for the different parameters in (6) we adopt the conventions from [32]. Nuclear radius according to with skin thickness given by s = 0.9 fm, a = 0.52 fm and c = (1.23A 1/3 − 0.6) fm. CνNES can take place with non-zero spin nuclei as well. For example, of the six most abundant Xe isotopes ( 129 Xe, 130 Xe, 131 Xe, 132 Xe, 134 Xe, 136 Xe) 129 Xe is a J = 1/2 fermion state whereas 131 Xe is a J = 3/2. To properly account for scattering with non-zero spin nuclei one has to include contributions from the axial-vector current in (3). These contributions, being related with the spin operator, are suppressed compared with those from the vector current. Thus, to a fairly good approximation, the CνNES is well-described by (4) regardless of the spin of the target nuclei.
Once NSI are introduced, the CνNES cross section is no longer flavor diagonal and so different neutrino flavors have different cross sections. For ν i − N it reads [33,34]: where the only difference with respect to (4) arises from Q SM → Q ν i , with Note that for light mediators (m 2 X q 2 ) this is no longer true, and a different recoil energy dependence is expected. If this is the case, the ν − N and WIMP-nucleon recoil spectra might differ, thus alleviating the effects of the solar and atmospheric neutrino backgrounds [21].

Solar neutrino flavor conversion and neutrino-nucleus scattering rates
Solar neutrinos (ν e ) are produced in regions near the solar core (r 0.2R ), with the exact distribution of neutrino production determined by the SSM. Once produced, electron neutrinos are subject to flavor conversion, governed by the vacuum and matter Hamiltonians according to Here |ν T = |ν e , ν µ , ν τ T refers to the neutrino flavor eigenstate basis, r to the neutrino propagation path, U = U (θ 23 )U (θ 13 )U (θ 12 ) is the leptonic mixing matrix (U (θ ij ) is a 3 × 3 rotation matrix), H vac = diag(0, ∆m 2 21 , ∆m 2 31 ) and in the absence of NSI H mat = √ 2G F n e (r)diag(1, 0, 0), with n e (r) the solar electron number density. Neutrino NSI induce additional matter potentials which change the matter Hamiltonian and thus affect neutrino flavor evolution. In full generality the evolution equation in (10) becomes where the NSI coupling matrices ε f read . The up-and down-quark relative abundances can be written in terms of the neutron relative abundance according to with the neutron number density in turn calculated from the 4 He and 1 H mass fractions (metallicity in the Sun amounts to less than 1.5% so it can be safely neglected), namely Since we are interested in CνNES, from now on we set e ij = 0 and following refs. [9,14] we do not consider simultaneous contributions from up-and down-quark couplings. For all the SSM related quantities we use those from the BS05 SSM [16].
To a fairly good approximation, flavor conversion probabilities can be calculated in the mass dominance limit, ∆m 2 31 → ∞. In this limit, and with vanishing CP-violating phases, neutrino propagation takes place in the basis |ν = U 13 T U 23 T |ν , which neglecting terms proportional to sin θ 13 corresponds to the following neutrino eigenstates: In this basis, which we dub as propagation basis,ν e −ν τ andν µ −ν τ mixing is of order ξ ij = G F n q (r) q ij /(∆m 2 31 /E ν ). Taking n q 10 26 cm −3 (up-quark number density at about r = 0.1R ), E ν 20 MeV and ∆m 2 31 = 2.55 × 10 −3 eV 2 (BFPV for normal ordering) [35], one can see that even for order-one NSI couplings ξ ij 1. This means that in the mass dominance limit,ν τ is decoupled and so neutrino flavor evolution can be studied in the two-flavor approximation with electron neutrino flavor conversion determined by 1 − P ee (E ν ) . Thus, in what follows we discuss the calculation of the averaged survival probability P ee (E ν ) . First of all, due to the quark and electron densities P ee depends not only on E ν but also on the neutrino propagation path r. P ee (E ν , r) can be written as [36] P ee (E ν , r) = cos 4 θ 13 P eff (E ν , r) + sin 4 θ 13 , where the r dependence is introduced by the effective probability given by [37] The mixing angle in matter is calculated from the diagonalization of the 2 × 2 Hamiltonian where Note that in the limit q ij = 0 and cos θ 13 = 0, A reduces to the SM term and B vanishes. The parameters ε D and ε N result from the rotation from the flavor to the propagation basis and read [9]: Here we have adopted the notation cos θ ij = c ij and sin θ ij = s ij . From eq. (18) the mixing angle in matter can be straightforwardly calculated: cos 2θ M (r) = ∆m 2 12 cos 2θ 12 − A (∆m 2 12 cos 2θ 12 − A) 2 + (∆m 2 12 sin 2θ 12 + B) Eqs. (16) and (17) combined with eqs. (19)- (22) allow the determination of P ee (E ν , r) in terms of neutrino oscillation parameters, electron and quark number densities and NSI parameters. The averaged survival probability is then obtained by integrating over r taking into account the distribution of neutrino production in the Sun [9]: where Φ α (E ν ) stands for neutrino fluxes and ρ α (r) for the distribution of neutrino production, α = pp, 7 Be (ground and excited state), pep, hep, 8 B, 13 N, 15 O, 17 F. In practice, to calculate P ee (E ν ) we have followed a single-parameter analysis, varying only one NSI coupling at a time in the ranges allowed by the constraints derived from neutrino oscillation and COHERENT data [14]. These ranges have been split according to q ij | n+1 = q ij | n + δ q ij with δ q ij = q ij | max − q ij | min /10. And the neutrino energy, E ν , has been varied from 0.145 MeV to 14.5 MeV. The neutrino oscillation parameters have been fixed to their BFPVs according to ref. [35]. For each NSI parameter we thus have generated 11 data sets (11 averaged survival probabilities), including the SM case ( q ij = 0). Focusing on up-quark NSI couplings (results for down-quark parameters resemble those of the upquark case), the averaged survival probability deviates the most from the SM expectation for flavordiagonal NSI and u eτ , being more pronounced for the former than for the latter. Fig. 1 shows the result for P ee for the cases u ee = 0 and u µµ = 0 (results for u τ τ = 0 are pretty close to those found for u µµ = 0). It can be seen that the range of variation of P ee has about the same amplitude in both cases. However its behavior with increasing values of the NSI couplings is rather different. For u ee = 0, the averaged survival probability decreases with increasing values of u ee and its largest value (as a function of E ν ) is indeed found when ee = 0. For u µµ = 0, it is the other way around. Increasing values of the parameter leads to larger P ee . This behavior can be readily understood as follows. In the SM case limit, u ij = 0, eq. (22) reduces to The presence of the u ee and u µµ couplings changes this equation by shifting either the first term in the denominator (numerator) or both. As can be seen from eqs. (19)- (22) the shift in the numerator of eq. (24) due to u ee is always negative. Thus, as u ee increases cos 2θ M decreases and accordingly P ee decreases (see eqs. (16) and (17)). The shift due to u µµ is instead positive thus leading to larger values of cos 2θ M as u µµ increases. As we have pointed out, in the mass dominance limit solar neutrino flavor conversion can be reduced to a two-flavor oscillation problem, in which a mainly |ν e state oscillates to a |ν µ state which is a superposition of |ν µ and |ν τ flavor eigenstates. Thus, in an Earth-based detector, neutrino-nucleus (left-hand graph) and u µµ (right-hand graph). Solar neutrino oscillations are analyzed in the mass dominance limit and soν µ c 23 ν µ − s 23 ν τ . The shaded regions indicate current allowed ranges for both parameters as required by neutrino oscillation and COHERENT data [14].
scattering will take place with either |ν e or |ν µ and so the calculation of neutrino-nucleus scattering rates has to be done in the propagation basis as well. Rotating the effective Lagrangian in (1) to this basis, after dropping the axial-vector coupling matrix, results in with the NSI couplings in the propagation basis related with the NSI couplings in the flavor basis through the following transformation The differential cross sections forν e − N ν e − N andν µ − N can then be obtained from (8) and (9) by q ij →˜ q ij . Results for the cases u ee = 0 and u µµ = 0 are shown for both neutrino flavors in fig. 2. As can be seen, deviations on either dσ νe /dE r or dσν µ /dE r are possible. The presence of NSI can enhance but also can deplete the flavored cross sections, thus potentially leading to substantial deviations of the expected neutrino-nucleus scattering event rate.
Assuming a 100% efficiency and with the oscillation probabilities and the ν − N differential cross sections calculated in the propagation basis, the ν − N event rate as a function of the recoil energy E r can be written according to [38]: Here N T = N A /A (N A = 6.022 × 10 29 ton −1 ), E min r = m N E r /2 and we include neutrinos emitted in the pp and CNO chains, α =pp, 7 Be (ground and excited state), pep, hep, 8 B, 13 N, 15 O, 17 F, as in (23). The neutrino fluxes are in units of cm −2 year −1 keV −1 , so the differential rates are in units of ton −1 year −1 keV −1 and the total event rates, R ν , are then measured in ton −1 year −1 . The latter are calculated from (27) by integration over E r from E th r to 100 keV, with E th r ⊂ [10 −3 , 10 2 ] keV and varied in steps δE r = (E n+1 r − E n r )/50.

Atmospheric neutrino-nucleus scattering rates
In contrast to solar neutrinos, matter effects for sub-GeV atmospheric neutrinos are negligible and flavor conversion can be well described by vacuum oscillations [17]. This can be readily seen by considering the PREM model [39], for which n ⊕ e = 3.25 × 10 24 cm −3 and Y ⊕ n = 1.137 both in the Earth core. With these values, using the BFPV for ∆m 2 31 [35] and taking E ν = 100 MeV so to minimize the vacuum term contribution one gets (see eq. (13) Note that this is the largest contribution from H mat to the neutrino flavor evolution equation, thus implying that matter effects for sub-GeV neutrinos can be fairly ignored.
Since only vacuum oscillations are relevant, the full three-flavor treatment is straightforward. In this case the oscillation probabilities can be further simplified by considering the limit ∆m 2 21 → 0. They read [40] where L refers to the distance covered by a neutrino from production at a height h (perpendicular to the surface) to detection and is given by: with θ z the zenith angle of the neutrino trajectory. The oscillation probability in (29) is therefore a function of production height h and zenith angle θ z . The distribution of vertical production height for sub-GeV neutrinos was calculated in [41] and its mean was found to be h 18 km. Thus, rather than using the distribution we fix h to this value and integrate over the zenith angle. The oscillation Shaded regions refer to current limits on u ii (i = e, µ, τ ) as required by neutrino oscillation and CO-HERENT data [14]. The blue shaded region ( u ee ) covers values down to 0.028 (see tab. 1). These are the cross sections which deviate the most from the SM value and that therefore lead to the largest deviations on the number of neutrino-nucleus scattering events in the atmospheric sector.
probabilities relevant for the determination of the ν − N event rate can then be written as As expected, neutrino flavor oscillations in this case are only relevant for ν µ − ν τ conversion. The result is shown in fig. 3. Since we set CP-violating phases to zero, the differential event rates associated with ν e +ν e and ν µ +ν µ fluxes can be written according to dR Atm In (33) we assume time reversal invariance. We use the atmospheric neutrino fluxes from ref. [41], derived with the FLUKA Monte Carlo simulation package [42] and which provides fluxes for ν e , ν µ and the corresponding anti-neutrinos. The full differential ν − N event rate is obtained from (32) and (33) and the total event rate as in the solar case, integrating over E r from E th r to 100 keV, with E th r ⊂ [10 −3 , 10 2 ] keV and varied in steps δE r = (E n+1 r − E n r )/50. The flavored cross sections are then given by (8) for which one finds sizable deviations from the SM only for diagonal NSI couplings, as expected given the tight constraints on the non-diagonal, particularly u eµ and u µτ . The result for diagonal couplings is shown in fig. 4, in which all couplings but one have been assumed to be zero.

NSI at ton-size DM Xe detectors
We now turn to the discussion of the impact of neutrino NSI on the expected number of events in a Xebased DM detector. For concreteness we assume an exposure of E = 1 ton × year. For practical reasons, and using the same approach we used in the calculation of the averaged survival probability in the Sun (see sec. 3.3), we perform a single-parameter analysis, considering only up-type couplings (results for down-type couplings resemble those found for the up-quark case). We first calculate the neutrino-nucleus differential event rate for ten equally-spaced values of the corresponding NSI parameter: . From these results we then identify the parameters that maximize/minimize the recoil spectra. And for these parameters we calculate the neutrino-nucleus scattering events. For solar neutrinos fluxes we employ the results derived in sec. 3.3, eq. (27) combined with the averaged survival probabilities. For atmospheric neutrino fluxes we use the results from sec. 3.3. The results of our calculation are displayed in fig. 5, which we regard as the main result of this paper.
One can see that even those parameters that are forced to be small by phenomenological constraints can have a sizable impact on the number of events. Departures from the SM value (horizontal line at 1 in fig. 5) are energy dependent and are particularly pronounced in the solar-atmospheric transition, at E r 5.8 keV when hep neutrinos reach their energy tail. Note that in the absence of NSI propagation effects, the SM and NSI neutrino-nucleus scattering rates differ only by a global numerical factor, Q SM → Q ν i (see eq. (9)). Thus, this recoil energy dependence is a consequence of NSI propagation effects. Neutrino NSI not only lead to enhancements of the total number of events but can as well produce depletions (as firstly noted in ref. [20]), with both being relevant for DM direct detection experiments. The former places stronger limits on the sensitivity one can reach in a multi-ton DM detector, while the latter enables measurements of WIMP-nucleon cross sections that otherwise would be challenging reaching.
Overall, one observes that no matter the coupling a deviation from the SM expectation is always possible. However, diagonal couplings lead to the largest deviations which can either be enhancements or depletions of the event rate. Particularly important for solar neutrinos, are deviations in the energy range where the background is dominated by 8 B, since it is in that range where solar neutrinos can mimic a WIMP-nucleus signals [6]. In there, our findings show that u µµ leads to enhancements of R NSI ν which can exceed the SM value by about 40%. For u τ τ , enhancements are slightly smaller, but still substantial. For u ee they reach values of order 2 or so, being for this coupling the largest possible. This is somehow expected, as is for this coupling that bounds are less severe. Note that compared to what was found in ref. [20] this enhancement is larger, with the difference just due to the range employed. In ref. [20] u ee | max = 0.45, value derived from from oscillation and deep inelastic scattering data [15], while in our case u ee | max = 0.6. This difference, although small, leads to a sufficiently large enhancement of the neutrino-nucleus cross section (see fig. 2) that accounts for the larger deviation we found.
Since the largest deviations are found for diagonal NSI we focus in these cases, discussing in a bit more detail the results for u ee and u µµ . Maximization or minimization of R ν is determined by the way in which the NSI parameters affect detection and propagation. The latter relevant only in the solar sector. For u ee we found R NSI ν | max at u ee | max = 0.6 for both, the region where recoils are controlled by solar neutrinos fluxes (E r 5.8 keV) and the region governed by atmospheric neutrinos (E r 5.8 keV). This can be understood as follows. In the solar region dσν µ /dE r does not depend on u ee and so it keeps always its SM value (see fig. 2). Theν e − N cross section instead reaches its maximum exactly at u ee | max . The survival probability behaves in the opposite way, as u ee increases P ee decreases globally, reaching minimum values exactly at u ee | max (see fig. 1, left graph). This effect tends to deplete the contribution of theν e − N scattering process to the total number of events. However, since for u ee = 0.6 the survival probability decreases at most 20% from its SM value, the overall effect is not sufficiently strong resulting in R NSI ν | max at u ee | max . R NSI ν | min is found at about the point where dσν e /dE r reaches its minimum, u ee 0.2. In the atmospheric sector, due P eµ , P eτ 1 the event rate is mainly controlled by the ν e − N cross section, which can be largely modified by NSI effects (see figs. 3 and 4). This in turn is reflected in the non-depletion of R NSI ν when passing from the solar to the atmospheric regions. The parameter dependence in the case of u µµ is rather different. For solar neutrinos, R NSI ν | max and R NSI ν | min are found where theν µ − N cross section is maximized/minimized (see fig. 2, right graph). When solar fluxes are no longer relevant and atmospheric kick in, the maximum expected number of events decreases. From (32) one can see that contributions from ν e +ν e fluxes are mainly SM, sizable deviations are due only to ν µ +ν µ fluxes. They are however not as large as in the case of u ee due to the sizable value of P µτ which tends to deplete the value of the second term in eq. (33), thus leading to small R NSI ν | max . The parameter dependence in the case u τ τ is rather similar to the one found in the u µµ scenario in the solar sector, but in the atmospheric region the behavior is different. From eqs. (32) and (33) this is somehow expected. Contributions from electron neutrino fluxes barely deviate from SM values, substantial deviations can only arise from muon neutrino fluxes, with the ν τ − N cross section weighted by the ν u − ν τ oscillation probability which is of order 0.35 (see fig. 4). This results in slightly larger deviations that those found in the u µµ case. Non-diagonal NSI couplings are subject to tighter constraints, thus resulting in smaller effects. The exception being the case of the u eτ parameter, for which we found no depletion and enhancements of order 40%, inline with what was reported in [20].
The exact value (range) of the mediator mass for which our analysis is valid depends-of course-on the recoil energy and hence on the neutrino sector, solar or atmospheric. For the u µµ = 0 case ( fig. 5 left top panel), we show the minimum mass that the mediator should have. In the solar sector masses below ∼ 120 MeV suffice, with the exact value depending on the specific recoil energy. For recoil energies associated with atmospheric neutrinos fluxes, m med 120 MeV. Thus, while for solar neutrinos most of the mediator mass range is available, for atmospheric neutrinos the analysis is limited to a narrow window. However, as we have already pointed out is for this window where NSI effects may be expected to be maximized.

Implications for background-free sensitivities
CνNES leads to recoil energy spectra which are closely degenerate with that induced by WIMP-nucleus scattering. As a result, once a direct detection experiment becomes sensitive to the neutrino background there is a limit on the WIMP cross section (σ DM−n ) that such experiment can explore, due to the inability to distinguish a neutrino-from a WIMP-induced signal. This limit defines the so-called neutrino "floor".
With the CνNES process enhanced (diminished) by the neutrino-nucleus NSI, one expects the neutrino floor to be more (less) severe. Thus, in order to quantify the impact of NSI on the reach of next-generation DM detectors, we have quantified the 90% CL sensitivities achievable. To do so, we have followed the method introduced in ref. [6] which enables the determination of the smallest WIMP cross section that can be measured for any WIMP mass in the presence of one-neutrino background event. It is worth pointing out that the limit on the WIMP cross section derived this way is just an estimate of what one should really expect. The precise determination of the saturation point, where the WIMP cross section one can access does not reduce with increasing (feasible) exposure, requires abandoning the background-free paradigm and using instead a likelihood analysis [6]. Here, however, we aim just at a preliminary understanding of the possible impact that NSI might have on direct detection sensitivities, and therefore stick to a simple background-free analysis.
We calculate the WIMP cross section as a function of the WIMP mass by fixing the WIMP-nucleus scattering event rate, R DM , to 2.3 (90% CL sensitivity limits), namely Here the WIMP event rate, which depends upon astrophysical, nuclear physics and particle physics parameters [32] reads where E is the exposure, ρ DM is the DM density in the local halo, σ DM−N is the DM-nucleus cross section, µ N is the DM-nucleus reduced mass, v 0 is the most likely speed of the DM in the halo, F 2 (E r ) For comparison the SM 90% CL background-free sensitivity curve is shown (dotted purple curve). These sensitivities have been derived by adjusting the exposure such that any measurement contains one-neutrino event [6]. It has been calculated for the case u ee = 0.
is the form factor given in (6) and v is the velocity of the DM particles. The velocity distribution is Maxwellian and truncated at the DM escape velocity: |v + v lab | = v esc , with v lab the laboratory (Earth) velocity relative to the DM local halo. The factor k 0 /k is a normalization factor obtained by integrating the Maxwellian DM velocity distribution from 0 to v esc and is given by Rather than integrating numerically the velocity distribution in (35) we have used its exact analytic form, namely [32] vesc v min with v min = 2E r /m DM r ×c and r = 4m DM m N /(m DM +m N ) 2 . We have written as well the DM-nucleus cross section in terms of the DM-nucleon cross section with µ n the DM-nucleon reduced mass. For the different astrophysical quantities we have used: ρ DM = 0.3 GeV/c 2 /cm 3 , v 0 = 220 km/s, v lab = 232 km/s and v esc = 544 km/s. For the integration we varied E th r from 10 −3 to 10 2 keV (with 10 3 partitions) and took with r given after eq. (37). The exposure has been tuned so to assure that each E th r isocontour in the σ DM−n − m DM plane contains one-neutrino event. This has been done by calculating the exposure as follows: E(E th r ) = 1 neutrino event The result is shown in fig. 6, where we have plotted the SM sensitivities as well and included current (future) exclusion limits from XENON1T and PandaX-II [4,43] (XENONnT and PandaX-4T [44]). The sensitivities shown in the plot refer to those obtained by calculating the exposure in (40) for the values of u ee that maximize/minimize dR ν /dE r (lower and upper curves in fig. 5, left-top graph). The result in fig. 6 show that the effects of neutrino-nucleus NSI, if present, work either way. They can worsen sensitivities or improve them. Note, that despite the effect being small it might have consequences. For example, if the NSI-induced depletion is effective (lower purple curve) the background-free sensitivities at XENONnT will be less affected. If instead the NSI-induced enhancement is present, XENONnT or PandaX-4T could start measuring neutrino scattering events sooner than expected.

Conclusions
The identification of a DM signal in next-to-next generation (multi-ton scale) detectors will be seriously challenged by irreducible solar and atmospheric neutrino backgrounds. Identifying all the features of this background is therefore mandatory to determine the implications it will have and the challenges it will pose. With this in mind, in this paper-following [20]-we have entertained the possibility that neutrinos are subject to NSI and have studied the impact they would have on the expected number of neutrino-nucleus scattering events in a Xe-based detector assuming an exposure of 1 ton · year. We have as well determine the background-free sensitivities achievable if these interactions are present, showing that they can either worsen or improve them, although only slightly in both cases.
Despite the recent observation of the CνNES process by the COHERENT experiment and its consistency with the SM expectation, NSI couplings can still be relatively large, even when the COHERENT scattering data is combined with data from neutrino oscillation experiments. Thus, by using current constraints applicable in scenarios where neutrino NSI arise from mediators with masses in the range [10, 10 3 ] MeV, we have calculated the expected number of neutrino-nucleus scattering events in the presence of NSI including propagation and detection effects and extending the analysis over the solar and atmospheric neutrino sectors. In the former case, we have included all pp and CNO solar neutrino fluxes, while in the latter case we have included ν e , ν µ ,ν e andν µ fluxes for neutrino energies up to ∼ 100 MeV. To simplify the computational burden, we have carried out a single-parameter analysis in which all couplings were put to zero but one and repeated the analysis for the remaining NSI couplings. Since we used constraints from [14], we analyzed up-and down-quark couplings independently and presented only the results for the former, as the results for the latter are rather similar.
For the solar sector, we have studied propagation in the mass dominance limit and have taken into account the distribution of neutrino production within the Sun, for the calculation of the averaged neutrino survival probability. In the atmospheric region, Earth matter effects are negligible due to the largest relevant neutrino energies being below ∼ 100 MeV. In this sector then, NSI affect the expected number of neutrino-nucleus scattering events only through detection. We have anyway included vacuum oscillation probabilities (properly integrated over the zenith angle) in the determination of the number of scattering events.
Our main results are given by fig. 5 and fig. 6 and can be summarized as follows. Largest deviations from the SM expectation are more pronounced for flavor-diagonal NSI and can operate either enhancing the number of events or depleting them. For the most relevant recoil energy interval in the solar sector (that in which the neutrino flux is determined by 8 B neutrinos in which the background can mimic WIMP-nucleon scattering events) the largest enhancement is found for u ee and can exceed the SM expectation up to a factor ∼ 2. For u µµ and u τ τ , enhancements are smaller but still sizable. Depletions in that energy interval can amount to about 30-40% for all diagonal NSI parameters. Moving to higher neutrino energies we found that when atmospheric neutrino fluxes dominate the background, enhancements as large as twice the SM expectation are possible when u ee is present, in all other cases they are below this value but still are sizable. Depletions in the atmospheric sector can amount to about 50% if u τ τ is at work. For the off-diagonal couplings u eµ and u eτ we did not find any depletion, but we did for u µτ although rather small (order 5%) and so probably not experimentally measurable. Thus, if the neutrino-nucleus event rate for energies dominated by atmospheric neutrino fluxes is found to be below SM expectations-arguably-that will support the presence of new physics featuring flavor-diagonal NSI, assuming one can disentangle the new interactions above any other effect.
Deviations, though small, might have an effect on the sensitivities achievable. To estimate the extent of such effects we have calculated the background-free sensitivities. We have shown that neutrinonucleus NSI, if present, may worsen sensitivities or improve them. The effect is small, but might have observable consequences. We have stressed that if neutrino NSI are present either XENONnT or PandaX-4T, among others, will start observing neutrino-induced nuclear recoil events sooner than expected.