COHERENT analysis of neutrino generalized interactions

Effective neutrino-quark generalized interactions are entirely determined by Lorentz invariance, so they include all possible four-fermion non derivative Lorentz structures. They contain neutrino-quark non-standard interactions as a subset, but span over a larger set that involves effective scalar, pseudoscalar, axial and tensor operators. Using recent COHERENT data, we derive constraints on the corresponding couplings by considering scalar, vector and tensor quark currents and assuming no lepton flavor dependence. We allow for mixed neutrino-quark Lorentz couplings and consider two types of scenarios in which: (i) one interaction at the nuclear level is present at a time, (ii) two interactions are simultaneously present. For scenarios (i) our findings show that scalar interactions are the most severely constrained, in particular for pseudoscalar-scalar neutrino-quark couplings. In contrast, tensor and non-standard vector interactions still enable for sizable effective parameters. We find as well that an extra vector interaction improves the data fit when compared with the result derived assuming only the standard model contribution. In scenarios (ii) the presence of two interactions relaxes the bounds and opens regions in parameter space that are otherwise closed, with the effect being more pronounced in the scalar-vector and scalar-tensor cases. We point out that barring the vector case, our results represent the most stringent bounds on effective neutrino-quark generalized interactions for mediator masses of order $\sim 1\,$GeV. They hold as well for larger mediator masses, case in which they should be compared with limits from neutrino deep-inelastic scattering data.


Introduction
The coherent elastic neutrino-nucleus scattering (CEνNS) process has been recently observed by the COHERENT experiment [1], more than 40 years after its first theoretical description [2]. Compared to other neutrino processes at energies below 100 MeV, CEνNS has a large cross section with a value of order 10 −39 cm 2 , due to the enhancement induced by the square of the number of neutrons in the nucleus. However, despite these large values the CEνNS eluded experimental detection for years due to the complicated measurement of the weak nuclear recoil energies (∼ few keV) produced in the interaction. Its measurement became possible thanks to the development of ultra-sensitive technology in other experimental searches namely, rare decays and weakly interacting massive particle dark matter (DM) [3].
CEνNS occurs when the de Broglie wavelength of the scattering process is larger than the nuclear radius (λ = h/q R N , where q refers to the exchanged momentum), which for typical nuclei translates into q 200 MeV. Accordingly, in ν − N scattering processes in which q is sufficiently small the scattering amplitudes on single nucleons add coherently and lead to an enhanced cross section whose value depends upon the number of nucleons within the nucleus. In the standard model (SM) the CEνNS process is well understood and it is determined by Z boson exchange [2]. It receives contributions from vector and axial nuclear currents, with the latter being-of course-relevant only for nuclei with spin J = 0 [4]. However, even in that case, it is well known that the axial contribution is relevant only for light nuclei [4] and negligible for heavy ones, such as Cs and I used in the COHERENT detector [1].
The COHERENT experiment uses neutrinos produced in the spallation neutron source (SNS) at the Oak Ridge National Laboratory. The spallation process starts with negatively charged Hydrogen ions H − which are accelerated at a LINAC. After being accelerated at ∼ 0.9 c, the two electrons in the H − ions are stripped off and the resulting protons are accumulated in a storage ring. Spallation takes place when 60 Hz proton pulses hit a liquid mercury fixed target. In that process not only neutrons but also pions are produced from spallation. The neutrinos used by COHERENT are thus generated In the SM both vector and axial quark currents are present. Contributions from heavy new physics can be parameterized by a larger set of couplings subject only to the condition of Lorentz invariance. The most studied case of such parameterization corresponds to neutrino NSI [31], where the couplings have a SM-like structure, but are controlled by free parameters that "measure" the relative strength of the new interaction to the Fermi interaction (G F ). Neutrino NSI, however, are a subset of a whole set of interactions which include scalar, pseudoscalar, vector, axial and tensor couplings, which we refer to as neutrino generalized interactions (NGI). They may emerge in BSM scenarios in which e.g. neutrinos couple to heavy scalars [32] or in models where neutrinos have non-vanishing electromagnetic couplings [33]. If we were to consider such scenarios additional constraints from the charged lepton sector should be accounted for 1 , but here we do not consider this possibility and rather stick from the very beginning to non-gauge invariant dimension six operators. In doing so, we then place constraints on the new effective couplings by requiring consistency with the COHERENT measurement. In our analysis we focus on the leading contributions, which means that we do not consider pseudoscalar nor axial quark currents. These are spin-dependent interactions that-as in the SM case-lead to suppressed contributions. It is worth emphasizing that our analysis is complementary to those presented in refs. [9,28] and extends upon these studies by including scalar effective interactions, crossed Lorentz structures and simultaneous presence of different nuclear currents.
The rest of this paper is organized as follows. In sec. 2 we provide a short overview of the CO-HERENT experiment and discuss the definitions and conventions used to perform our analysis, such as theoretical neutrino fluxes and calculation of number of events. We also define the binned χ 2 function and the different measured quantities that are involved. In sec. 3 we present the parametrization for NGI starting with neutrino-quark interactions and ending up with neutrino-nucleus couplings. We provide relations between the quark and nucleus couplings for scalar, vector and tensor currents. In sec. 4 we present our results for the differential cross section and the constraints implied on the effective neutrino-quark parameters by COHERENT data. Finally, in sec. 5 we summarize and present our conclusions. In appendix A we provide details of the cross section calculation for neutrinos and antineutrinos in the zero-momentum limit, including the full set of generalized interactions for spin−1/2 nuclei.

CEνNS signal rate at COHERENT
COHERENT uses neutrinos produced in the Spallation Neutron Source at the Oak Ridge National Laboratory [1]. The interaction of a pulsed proton (∼ 1 GeV) beam with a fixed mercury target produces neutrons from spallation and a substantial amount of low-energy neutrinos, which stem from the decay of stopped pions and muons, π + → µ + + ν µ and µ + → e + + ν e +ν µ . Muon neutrinos-being the by-products of a two-body decay-are monochromatic, and their energy is determined by the pion and muon masses: E νµ = (m 2 π − m 2 µ )/2m π 30 MeV. Accordingly, their energy distribution is given by [35] Electron neutrinos and muon anti-neutrinos instead feature continuous spectra. Their energy distributionnormalized to one-can be read off from the µ + (unpolarized) differential rate, namely where the energy distribution functions are given by [35] with the kinematic end point located at E ν = m µ /2 52.8 MeV. The neutrino flux (per flavor) that reaches the CsI[Na] detector, φ α (E να ) (α = ν µ ,ν µ , ν e ), is then determined by the energy distribution functions in eqs. (1)-(3) times the total number of neutrinos per each flavor, N. The latter is fixed by the number of neutrinos produced per proton collision (r = 0.08 per flavor), the distance from the source to the detector (L = 19.3 m) and the number of protons-on-target (POT, n POT ). For the 308.1 live-days of neutrino production, n POT = 1.76 × 10 23 [1]. Thus, since neutrinos are isotropically produced, N = r × n POT 4πL 2 and The COHERENT detector consists of m det = 14.6 kg of CsI[Na], where the sodium dopant is present with a fractional mass of 10 −5 − 10 −4 and so it does not play any substantial rôle as a target. Notice also that since A Cs A I , both Cs and I yield approximately the same nuclear response. The number of target nuclei is therefore given by n N = 2m det m CsI × N A [28], where m CsI = 2.598 × 10 −1 kg/mol is the CsI molar mass and N A is the Avogadro number.
For a given flavor α and taking into account both the Cs and I nuclei, the expected number of events in the i-th recoil energy bin reads Here E min ν = 2m Na E r (E r refers to the nuclear recoil energy and m Na to the nucleus mass), E max ν = m µ /2 and f a are nuclear fractions: f Cs = 51% and f I = 49%. The observed number of photoelectrons (PE) is related to the recoil energy through [1] In terms of n PE the COHERENT signal covers 25 bins starting from n PE = 1 and extending up to n PE = 49, with bin size equal to 2 photoelectrons. Since the acceptance function vanishes for n PE ≤ 5 (see below), the first three bins contain no information on the scattering process. Furthermore, from n PE ≥ 31 the relation between the number of photoelectrons and the nuclear recoil energy in (6) does not hold anymore. Thus, in our analysis we consider only 14 data bins, from n PE = 7 to n PE = 31 2 , assuming that at n PE = 31 eq. (6) is still valid (excluding this bin has no significant impact in our results). In terms of n PE the recoil energy integration limits are (n PE ∓ 1)/1.17. In our calculation we employ the nuclear Helm form factor (see discussion in sec. 4): where j 1 (x) is the order-one spherical Bessel function, q a = 6.92 × 10 −3 √ A a E r fm −1 and the effective nuclear radius is given by r n = (c 2 + 7π 2 a 2 /3 − 5s 2 ) 1/2 , with s = 0.9 fm, a = 0.52 fm and c = (1.23A 1/3 a − 0.6) fm [37]. The acceptance function A(x) is given by 3 where k 1 = 0.6655, k 2 = 0.4942, x 0 = 10.8507 and H is the Heaviside function. The CEνNS differential cross section dσ a α /dE r depends on the nuclear target and in BSM physics scenarios can be flavor dependent. In the SM it arises from the neutral current vector and axial-vector couplings [2], with the axial contribution terms being subdominant [4]. The leading contribution can be written as follows: Here showing that for heavy nuclei the CEνNS cross section is largely enhanced. From eqs. (5), (7), (8) and (9) we calculated the number of CEνNS events predicted by the SM. The result is shown in fig. 1, for the different F X separately (colored histograms) and for the total neutrino flux (black histograms), together with the COHERENT data with their corresponding uncertainties. As can be seen, these data closely follows the SM prediction [1]. However, due to the still large uncertainties, sizable contributions from BSM physics can be present and can therefore be constrained. As we have already pointed out, since the release of the COHERENT result various BSM scenarios have been analyzed. They include neutrino effective NSI [27,28], NSI via light mediators [9,15,28], neutrino four-fermion contact tensor interactions as well as electromagnetic neutrino couplings [9].
To constrain new physics contributions with COHERENT data, we use the following least-squares function 4 where N meas i is the number of events measured in the i-th bin, N NGI i is the number of events predicted by the NGI scenario (determined by the set of parameters P), σ 2 i is the statistical uncertainty on the experimental data in the i−th bin. The nuisance parameters α and β account for uncertainties on

Neutrino generalized interactions
Most studies of the CEνNS process in the presence of new physics are done assuming neutrino NSI [31], which are determined by the following four-fermion effective operator Here q(V,A) ij are free parameters which are constrained by neutrino oscillation and neutrino scattering data and q = u, d quarks [27,35,38,39] (see also ref. [40] for a review). These couplings parameterize the strength of the new interactions (relative to G F ). The operator in eq. (11) is actually more general and encodes other interactions. For example, it describes as well an effective theory involving operators such as (ν L i q R ) (q R ν L i ), as can be checked by Fierz rearrangement of the fermion fields. Nevertheless, eq. (11) is not the most general effective ν − q operator. A more general treatment is possible by  Figure 2: Beam-on background from prompt neutrons as a function of the number of photoelectrons n PE . It follows from the prompt neutron probability distribution function and it is weighted by the energy delivered during the 308.1 live-days of neutrino production, 7.48 GWhr [36]. Only PE bins considered in our analysis (7 ≤ n PE ≤ 31) are shown.
considering all Lorentz invariant non-derivative interactions of neutrinos with first generation quarks, namely (we use the notation employed in [41]) Here Γ X = {I, iγ 5 , γ µ , γ µ γ 5 , σ µν }, where σ µν = i[γ µ , γ ν ]/2 and without loss of generality the parameters C q X and D q X are real [41]. As in the NSI case, they "measure" the relative strength of the new physics and so their size is of order ( where m X is the mass of the exchanged particle and g X the coupling constant. Due to the quark axial current term, these interactions include diagonal and non-diagonal Lorentz structures. For example, Γ P involves pseudoscalar-pseudoscalar as well as pseudoscalar-scalar neutrino-quark couplings. Among the NGI, those that give the most relevant effect, in the sense that can sizeably diminish/exceed the SM contribution, do not involve nuclear spin. Indeed, effective couplings for nuclear spin-dependent interactions are determined by a sum over spin-up and spin-down nucleons, Z ↑ − Z ↓ and N ↑ − N ↓ (for proton and neutrons respectively). Therefore they are suppressed for all nuclei except for light ones [4]. Since our analysis involves heavy CsI nuclei we then drop the pseudoscalar and axial quark currents and we only keep scalar, vector and tensor quark currents. It is worth emphasizing that with this choice only the parameters P = {C q S , D q P , C q V , D q A and C q T } can be constrained. To compute the CEνNS cross section induced by the NGI we assume a fermion nuclear ground state with spin J = 1/2. This is motivated by the fact that nuclear matrix elements for nucleonic currents can in this case be borrowed from nucleon matrix elements for quark currents. Of course with such procedure one has to bear in mind that the corresponding nuclear form factors are different. In our case all the leading-order decompositions will involve the Helm form factor given in eq. (7) 6 . This is somehow expected given that after dropping the pseudoscalar and axial quark currents the remaining interactions become spin-independent and so they add coherently on the nucleons 7 .
To determine the effective neutrino-nuclear Lagrangian, from which we next calculate the cross section in the zero-momentum transfer limit, we start with the quark currents and we end up with nuclear currents following the procedure where O q,n,N refer to quark, nucleon (n = p, n) and nuclear operators, respectively. For step (I) one calculates quark currents in nucleons according to (see e.g. [43,44]) Here p i and p f refer to initial and final state nucleon momenta. The scalar current receives contributions also from heavy quarks (q = c, b, t), which are not of the form given in eq. (14). These contributions however are suppressed by m n /m q and so we do not consider them. Moreover, we neglect as well the contribution from strange quarks and from gluons and we keep only first generation quarks. For vector currents, the coefficients N n q can be understood essentially as the number of quarks within the nucleon, while for tensor currents δ n q represents a tensor charge. The factors f Tq are related with the fraction of the nucleon mass "carried" by a particular quark flavor. They are derived in chiral perturbation theory from measurements of the π − n sigma term [45]. The factors δ n q that we use here are derived from an analysis based on data from azimuthal asymmetries in semi-inclusive deep-inelastic scattering (DIS) and e + e − → h 1 h 2 X processes [46]. More recent values are given in refs. [47][48][49]. In our calculation we use the numerical values [46,50] For step (II) one evaluates the correlators of nucleonic currents in nuclei, which involve nuclear form factors and which can be written following Lorentz invariance, namely 6 We are assuming that the proton and neutron form factors are equal and well described by the Helm form factor, . A more precise approach in which F Z (q 2 ) is described by the Fourier transform of the symmetrized Fermi distribution and F N (q 2 ) by the Helm form factor could be adopted [29]. However, given the uncertainties of COHERENT data our description is precise enough. 7 This is what one finds in DM direct detection analyses: all nuclear interactions but the pseudoscalar and axial add coherently. Accordingly, apart from these two cases, the corresponding cross sections involve only the Helm form factor (see e.g. [37,42]).
where the momenta of the incoming and outgoing nucleus, p 2 and k 2 , define the exchanged momentum q = k 2 − p 2 . Some words are in order regarding these decompositions. F (q 2 ) refers to the Helm form factor in eq. (7) and it is in practice the only one relevant at leading order. The magnetic moment term in the vector current decomposition, as well as the second and the third terms in the tensor decomposition, are suppressed by O(q/m N ) factors. Thus, keeping just the leading terms, step (II) can be carried out and the neutrino-nucleus (ν − N ) effective Lagrangian can be written where the coefficients C X and D X correspond to ν −N effective couplings determined by the parameters C q X and D q X in eq. (12). Notice that from eq. (18) we can calculate the zero-momentum cross section, while the full cross section will involve the nuclear form factor which in turn will encode the momentum dependence (q 2 dependence). The ν − N coefficients are written as follows: The expression for D P is obtained from that of C S by trading C q S for D q P , while for D A from C V by trading C q V for D q A . The relations in eq. (19) allow to translate the constraints on the ν − N coefficients to the parameters of the "fundamental" Lagrangian.

Neutrino oscillations versus neutrino scattering
Before proceeding with the chi-square analysis, it is worth commenting on which other processes may set constraints on the NGI and on the range of validity of our results. As in the NSI case, interactions in eq. (12) contribute-in principle-to forward coherent scattering (order G F at q 2 = 0) and scattering processes (order G 2 F with q = 0). The former are responsible for matter potentials in matter and are related to neutrino oscillation data, while the latter include not only COHERENT but also DIS data from CHARM and NuTeV [51,52].
Matter potential induced by SM vector interactions in the Sun and in the Earth are responsible for resonant neutrino flavor conversion [31,53,54]. Accordingly, new contributions to the vector current are subject to both constraints, oscillation+scattering. This is indeed the case for neutrino NSI, where it is found that the combined analysis of oscillation+scattering data imply more stringent bounds [35,55]. Scalar interactions couple background fermions (nucleons) with different chiralities but same helicity, and so they lead to helicity suppressed matter potentials (m ν / E ν ) [56]. Constraints from neutrino oscillation data on these couplings are thus loose, if existing at all. Transverse tensor interactions, instead, can induce a sizable matter potential as they couple background fermions with different chiralities and opposite helicities (longitudinal tensor interactions are helicity suppressed as well). However, this tensor matter potential is only relevant in a polarized medium and so it does not sizeably affect neutrino propagation in the Sun or even in supernovae [56].
In the NSI case, DIS data places more severe bounds than COHERENT data does [35]. This should apply as well for the remaining interactions in (12). These limits however do not apply for mediators whose masses are below the typical momentum exchange in DIS processes, O(10 GeV). For m 2 X q 2 DIS , the relative value of the new contribution is σ BSM /σ SM ∼ g 4 X /q 4 /G 2 F and amounts to 1% for g X = 10 −2 . The same parameter choice with m X = 10 2 MeV and evaluated at q 2 COH (10 MeV) 2 gives σ BSM /σ SM 1. This means that for mediator masses below 10 3 MeV DIS constraints can be evaded and COHERENT bounds become dominant.
It follows that the constraints we derive here (see sec. 4) are the most stringent (for all interactions except the vector one), in scenarios where the mediator mass is below 10 3 MeV. For heavier mediators, more severe limits from DIS data may apply, but to the best of our knowledge such bounds do not exist.

Constraints from COHERENT data
To address the implications of COHERENT data on NGI, one has to calculate the number of expected events for a certain parameter choice according to eq. (5). This requires the determination of the corresponding cross sections for ν − N andν − N coherent scattering (the former has been derived in [41]). Starting from the Lagrangian in eq. (18) we calculate the zero-momentum differential cross section at leading order, i.e. neglecting O (E 2 r /E 2 ν ) terms: the index a denoting the target material. Details of the full calculation, including pseudoscalar and axial quark currents, are given in app. A. In the previous expression, E max r 2E 2 ν /m Na and the following definitions apply 8 The ξ X parameters defined in eq. (21) depend upon the nucleus, although for the sake of simplicity we have not written this dependence explicitly. The momentum-dependent cross section is then obtained from eq. (20) introducing the nuclear form factor Should an axial nuclear current be present, eq. (20) would contain two additional terms, corresponding to the axial contribution itself and to an interference term between the vector and axial currents. This axial-vector interference term as well as the last term in eq. (20) (proportional to R) are the only two that come with opposite signs in theν − N and ν − N cross sections. In the former (latter) case we find that the vector-axial interference term leads to constructive (destructive) interference. If we neglect pseudoscalar and axial nuclear currents then the neutrino and anti-neutrino elastic scattering  Table 1: Best-fit-point value (second column), 90% CL (∆χ 2 < 2.71, third column) and 99% CL (∆χ 2 < 6.63, fourth column) ranges for the ξ X (X = S, V, T ) parameters as defined in eq. (21). From these results one can then map to the fundamental neutrino-quark parameters using eq. (19). See text for further details.
cross sections differ only in the term proportional to R, which turns out to be relevant only if scalar and tensor interactions are simultaneously present. This term leads to rather suppressed differences and eventually the neutrino and anti-neutrino cross sections can be considered equal.
In full generality, the parametrization introduced in eq. (21) must include the SM as well. The SM limit is recovered when all the couplings but ξ V = C V are set to zero and ξ V = C V = 1 − (1 − 4 sin 2 θ w ) × N/Z. This contribution is of course always present throughout our analysis, and so from now on we will denote by ξ V the BSM contribution to the vector current. Note that the term proportional to ξ V has en extra term, E r /E ν , compared to the SM cross section for J = 0, eq. (9). This term is a consequence of the nuclear ground state spin, J = 1/2 [41]. From eqs. (19) and (21) one can see that even neglecting pseudoscalar and axial quark currents and without considering lepton flavor dependent couplings, the full problem involves 10 free parameters. In order to technically simplify the analysis, rather than considering the whole set, we stick to two kinds of simplified benchmark scenarios which we will discuss in the next subsections.

Single-parameter scenarios
We start our analysis by considering the single-parameter case parameterized in terms of the different ξ X . These "couplings" are related to the neutrino-quark couplings of the "fundamental" Lagrangian through the relations derived in eq. (19). Thus, in reality, by "single-parameter" scenarios we refer to the cases in which only one interaction at the nuclear level is present at a time. This, however, does not mean that the analysis reduces to a single parameter problem. Take for example the vector case. Vector nuclear currents arise from either Γ V Γ V or Γ A Γ A Lorentz structures, as can be seen by the definition of ξ V . As already introduced in sec. 3, in our analysis we only consider first generation quarks. Thus, in general, when considering the case in which all ξ X vanish except for ξ V , one is eventually dealing with a four-parameter problem (C u V , C d V , D u A and D d A ). Similarly, also the scalar interaction involves four fundamental parameters (C u S , C d S , D u P and D d P ), encoded in ξ S . The tensor current instead depends only on two parameters at the quark level, C u T and C d T . To facilitate the numerical analysis, for the scalar and vector cases we will consider only two neutrinoquark parameters at a time, for which there are six possible choices: (1-i) C q X = 0 or D q X = 0, (1-ii) C u X and D d X different from zero (or C d X , D u X = 0), (1-iii) C u X and D u X different from zero (or C d X , D d X = 0). Of these cases, (1-i) and (1-ii) lead to the same constraints over the different parameters. Constraints derived on C q X apply directly on D q X and those derived on C u X and D d X on C d X and D u X . Cases (1-iii) instead result in different constraints over the different couplings. Nevertheless, they should differ only by small values, given that the differences between the up and down couplings and masses are small (this is actually what is found in NSI analyses [35,55]). We thus consider only the first options in (1-i)-(1-iii) and C u T -C d T for the tensor interaction. For all scenarios we fit the COHERENT data by minimizing the least-squares function (eq. (10)) over the systematic nuisance parameters α and β, and then we calculate ∆χ 2 = χ 2 − χ 2 min . From this procedure we obtain the 90% and 99% CL allowed ranges for each ξ X . Our results are shown in tab. 1. Note that while the best fit point values (BFPVs) for ξ S and ξ T are zero, an additional vector current with ξ V = −0.113 (−1.764) (corresponding to the two minima of the ∆χ 2 (ξ V ) function) improves the COHERENT data fit. This is shown in fig. 3, where the black (red) colored histograms refer to the CEνNS number of events from the three neutrino flavors in the SM (SM plus a vector NGI, with ξ V = −0.113). Values for χ 2 min in both cases are also shown. The results in tab. 1 can be translated into the "fundamental" neutrino-quark parameters by using eq. (19). To do so one has to bear in mind that although the number of events receives contributions from Cs and I, the following simplification applies ξ 2 X Cs . As expected, given that f Cs f I , F 2 (q Cs ) F 2 (q I ) and m N Cs m N I (numerically we find ξ X I /ξ X Cs 0.95 for all X). We then derive the allowed 90%, 99% CL regions for the quark parameters for scenarios (1-i)-(1-iii) and for the tensor case in terms of C u T − C d T . Fig. 4 shows the result for the scalar and vector interactions for scenarios (1-i) and (1-ii) (results for scenario (1-iii) closely resemble those from (1-ii) and so we do not display them), while fig. 5 for the tensor couplings. It is worth emphasizing that the C q X , D q X couplings appearing in the different panels in fig. 4 are not independent, as can be seen from the translation of the ξ X parameters into the C q X couplings, eq. (19). Hence, even if we display the CL contour regions in two quark parameter planes, since the initial χ 2 function depends only on one ξ X , we keep using ∆χ 2 < 2.71 and ∆χ 2 < 6.63 to determine the 90%, 99% CL contours, respectively. Figure 4: 90%CL (∆χ 2 < 2.71, dark reddish) and 99%CL (∆χ 2 < 6.63, light reddish) allowed regions in the neutrino-quark couplings parameter space for scalar and vector interactions. Panels in the left column correspond to constraints for scenario (1-i), while those in the second column to scenario (1-ii). Results for scenario (1-iii) resemble those of scenario (1-ii) and so we do not display them. The dashed lines refer to the values determined by the ξ S,V BFPVs (see tab. 1). For ξ V the χ 2 function exhibits two minima and so for this case the result includes two non-overlapping regions.
As can be seen in fig. 4, among the constraints implied by COHERENT data those for scalar-type interactions are the most stringent. This can be understood as follows. Although the data still involves large uncertainties, one can see that it is rather consistent with SM expectations. The scalar interaction involves a cross section that substantially differs from that of the SM model, and so once added it worsens the fit. Furthermore, translation from ξ S to quark parameters involves nucleon mass fractions times nucleon-to-quark mass ratios which altogether amount to values of order 5 (see eq. (19)). Since C S is bounded from the constraint derived on ξ S , consistency demands a sort of cancellation between the up and down couplings contributions, as fig. 4 (top panels) shows. On the other hand, tensor couplings allow for a relatively large freedom even compared with vector parameters 9 , as depicted in fig. 5. This, however, does not mean that tensor interactions provide a better fit to data than vector
do, as demonstrated by ξ BFPV T = 0. It follows from the translation from nucleus to quark parameters, which in the vector case involve larger coefficients and so leads to narrower allowed regions. Finally, the presence of two minima in the ∆χ 2 (ξ V ) function translates into two separate linear bands in the C u V , C d V plane and in two concentric rings in the C u V , D d A plane.

two-parameter scenarios
In this case we allow for the simultaneous presence of two interactions at the nuclear level. Accordingly, we can distinguish three cases corresponding to ξ S − ξ V , ξ S − ξ T and ξ V − ξ T which involve eight and six quark parameters respectively. As in the one-parameter case, here we focus on smaller-though representative-regions of parameter space. To determine at which extent the presence of a second interaction modifies the constraints obtained in the single-parameter analysis we study three scenarios: Figure 6: Results for the two parameter case analysis. Dark (light) reddish regions correspond to 90%CL (∆χ 2 < 4.61) and 99%CL (∆χ 2 < 9.21) bounds for the neutrino-quark couplings. and C d T = 0. The chi-square function for this analysis becomes now a function of two parameters (ξ X 1 and ξ X 2 ). We present in tab. 2 the BFPVs and the 90% and 99% CL ranges for each ξ X i . The CL ranges for the parameter ξ X 1 are obtained minimizing the least-squares function over the nuisance parameters α and β and over the second interaction parameter ξ X 2 . In principle, the parameter R could also be constrained by COHERENT data. However, its contribution to the CEνNS cross section is subdominant with respect to the SM contribution. Moreover, it depends on the product of two of the fundamental quark couplings C q S , C q T . It turns out that COHERENT bounds are not competitive enough to constrain C q S and C q T via the R parameter, they are instead more stringently constrained by the requirement of perturbativity (understood as C q S,T ≤ 1, i.e. the NGI should not exceed G F ). The constraints given in tab. 2 can then be mapped into the parameters of the neutrino-quark Lagrangian in the same way as in the single-parameter analysis. Using the relations given in eq. (19) we present in fig. 6 the allowed regions for the fundamental parameters in scenarios (2-i)-(2-iii). We only show these three particular cases, but the results in tab. 2 and eq. (19) allow to investigate any case in which two nuclear interactions are simultaneously present.
These results imply that the presence of an additional interaction at the nuclear level relaxes the bounds on the fundamental neutrino-quark couplings. Indeed, the addition of an extra free parameter ξ X allows for more freedom in the values of the NGI parameters. Interestingly, COHERENT constraints on the vector interaction parameter ξ V are sizeably relaxed with the addition of an extra scalar or tensor interaction. This can be seen by studying the dependence of the ∆χ 2 function upon ξ V , depicted in fig. 7. The red solid curve shows the ∆χ 2 function in the single-parameter scenario where only ξ V is switched on, while the blue dashed curve refers to the two-parameter scenario with ξ V and ξ T simultaneously present and the black dotted to the two-parameter scenario with ξ V and ξ S both active. In all three cases the ∆χ 2 function has two minima, but the region between them is heavily modified when an extra interaction is added. In the region around −0.95 the extra vector interaction tends to cancel the SM contribution, thus worsening the fit. As can be seen, the extra contribution (either scalar or tensor) improves the fit by increasing the expected number of events in that region.  Figure 7: Dependence of the ∆χ 2 function on the ξ V parameter in three different scenarios: The singleparameter scenario with only ξ V (red solid), the two-parameter scenario with ξ V and ξ T (blue dashed) and the two-parameter scenario with ξ V and ξ S (black dotted).

Conclusions
We have studied a generic set of effective Lorentz invariant non-derivative neutrino-quark interactions (NGI). These interactions contain as a subset well-studied neutrino-quark NSI, but involve additional scalar, pseudoscalar, axial and tensor couplings. In contrast to vector interactions, they induce matter potentials that are either helicity suppressed or vanish in non-polarized media. Accordingly, they are poorly constrained by neutrino oscillation data. They instead contribute to scattering processes which set bounds on their values. We have considered the contributions of NGI to the CEνNS process and we have employed the recent COHERENT data to place constraints on the different effective parameters.
Our analysis includes scalar, vector and tensor quark currents and excludes pseudoscalar and axial quark couplings, which being spin-dependent are expected to be less constrained. We have considered diagonal as well as non-diagonal Lorentz structures, such as (νγ µ γ 5 ν)(qγ µ q) and (νγ 5 ν)(q q) and under the assumption of no lepton flavor dependence and of a spin-1/2 nuclear ground state, we have calculated the full CEνNS cross section for neutrinos and anti-neutrinos. In order to assess the impact that such interactions have on the CEνNS process, we have then carried out a chi-square analysis in two simplified benchmark scenarios. A first one where only one nuclear interaction is present at a time, dubbed singleparameter case, and a second where two are simultaneously present, called two-parameter case.
In the single-parameter case, our findings show that the scalar interaction is the most constrained, with the tightest bound found for the Lorentz mixed pseudoscalar-scalar coupling. In such a case the effective parameters are bounded to be smaller than 0.05 at 90%CL. For scalar-scalar couplings this bound is relaxed and the parameters can be of order one, but still in a rather narrow region in parameter space. Allowed vector NGI are also sizable reaching values as large as 0.85 at 90%CL, but again in two non-overlapping narrow stripes. We find that tensor interactions are the less constrained, with the reason being the translation between the nuclear to quark parameters, which involves "tensor charges" which are small, thus allowing for more freedom. Nevertheless none of these values lead to an improvement in the COHERENT data fit, as the BFPVs found in our analysis demonstrate.
In the two-parameter case, we have found that the presence of an additional interaction at the nuclear level relaxes the bounds on the fundamental neutrino-quark couplings. The addition of an extra free parameter ξ X allows the NGI to span over relatively larger regions in parameter space. In particular, the allowed ranges for the vector parameter ξ V are sizeably modified with the addition of an extra scalar or tensor interaction. In the region where ξ V tends to cancel the SM contribution, thus worsening the fit, the scalar or tensor contribution enables its improvement to values below 1σ. We have pointed out that further and perhaps more severe constraints on NGI can be derived by considering instead DIS scattering data from CHARM and NuTeV, as it turns out to be the case for neutrino NSI [35]. However, whether this is the case depends on the mass of the mediator responsible for the effective interaction. We have stressed that for mediator masses below ∼ 1 GeV our constraints can be regarded as the current most stringent bounds on NGI. For mediator masses above this value our results are still valid but should be confronted with those from an analysis using DIS data, which to our knowledge does not exist. At any rate, improvements on limits on NGI couplings generated by mediators with masses below 1 GeV will require further improvement of COHERENT data.
CEνNS offers a plethora of physics opportunities, allowing for tests of anomalously large neutrino magnetic moments, sterile neutrinos, new light degrees of freedom, among others [17]. The analysis presented in this paper, while revisiting COHERENT constraints on some BSM interactions already considered in the literature, further complements previous works by considering effective NGI with mixed neutrino-quark Lorentz structures and simultaneous presence of various neutrino-quark interactions, for which we have shown that COHERENT data still allows for sizable values. We begin the computation of the CEνNS cross section by describing the kinematics of the process. Incoming (outgoing) neutrino/anti-neutrino and nucleus four-momenta are labeled p 1 and p 2 (k 1 and k 2 ), as shown in fig. 8. In the lab-frame they are written as whereê r =ê z cos φ+ê x sin φ,ê r =ê z cos φ−ê x sin φ and φ is the scattering angle. The outgoing neutrino energy can be calculated as from which the nuclear recoil energy E r = E ν − E ν follows and from which in turn the maximum nuclear recoil energy is obtained for backward scattering: E max r 2E 2 ν /m N (E ν m N ). The matrix elements for the processν(p 1 )+N (p 2 ) J=1/2 →ν(k 1 )+N (k 2 ) J=1/2 and ν(p 1 )+N (p 2 ) J=1/2 → ν(k 1 ) + N (k 2 ) J=1/2 can be written according to M(ν + N →ν + N ) = G F √ 2 a v s (p 1 )P R Γ a v s (k 1 ) ū r (k 2 )Γ a (C a + iγ 5 D a ) u r (p 2 ) , M(ν + N → ν + N ) = G F √ 2 a ū s (k 1 )Γ a P L u s (p 1 ) ū r (k 2 )Γ a (C a + iγ 5 D a ) u r (p 2 ) .
Here s, s , r, r refer to spin indices and we sum over all Lorentz structures. The differential cross section is in general given by [ where we have averaged over final state spins. Implementing the kinematic relations in (23) and using FeynCalc [58,59] we arrive to the following expressions (the result for the tensor interaction was derived, as far as we know, for the first time in [60]) where we have dropped O(E 2 r /E 2 ν ) terms. For neutrinos the third and last terms are positive and negative respectively, while for anti-neutrinos the signs are opposite. Here the following conventions apply and Our result for anti-neutrinos differs from that found in ref. [41] in the vector, axial and mixed vector-axial terms (couplings). The energy dependence of those terms, however, is the same and so the differences are numerically small.