New Constraints and Prospects for sub-GeV Dark Matter Scattering off Electrons in Xenon

We study in detail sub-GeV dark matter scattering off electrons in xenon, including the expected electron recoil spectra and annual modulation spectra. We derive improved constraints using low-energy XENON10 and XENON100 ionization-only data. For XENON10, in addition to including electron-recoil data corresponding to about $1-3$ electrons, we include for the first time events with $\gtrsim 4$ electrons. Assuming the scattering is momentum independent, this strengthens a previous cross-section bound by almost an order of magnitude for dark matter masses above 50 MeV. The available XENON100 data corresponds to events with $\gtrsim 4$ electrons, and leads to a constraint that is comparable to the XENON10 bound above 50 MeV. We demonstrate that a search for an annual modulation signal in upcoming xenon experiments (XENON1T, XENONnT, LZ) could substantially improve the above bounds even in the presence of large backgrounds. We also emphasize that in simple benchmark models of sub-GeV dark matter, the dark matter-electron scattering rate can be as high as one event every ten (two) seconds in the XENON1T (XENONnT or LZ) experiments, without being in conflict with any other known experimental bounds. While there are several sources of backgrounds that can produce single- or few-electron events, a large event rate can be consistent with a dark matter signal and should not be simply written off as purely a detector curiosity. This fact motivates a detailed analysis of the ionization-only ("S2-only") data, taking into account the expected annual modulation spectrum of the signal rate, as well as the DM-induced electron-recoil spectra, which are another powerful discriminant between signal and background.

Currently, the most stringent direct-detection constraint on DM as low as a few MeV comes from XENON10, a two-phase xenon time projection chamber (TPC). When a DM particle scatters off an electron and ionizes a xenon atom in the liquid target, the recoiling electron can ionize other atoms if it has sufficient energy. An electric field accelerates the ionized electrons through the liquid, across a liquid-gas interface, and through a xenon gas region in which interactions between the electrons and xenon atoms create a scintillation ("S2") signal that is proportional to the number of extracted electrons and detected by photomultiplier tubes. XENON10 [22] has taken data consisting of events that have an S2 signal corresponding to one or more electrons, without an observable prompt scintillation signal ("S1"). The data corresponding to events with three electrons or less (n e 3) were analyzed in [2] and shown to constrain DM as low as a few MeV.
The main factor limiting the sensitivity of XENON10 is the large number of observed S2-only events and the absence of a background model (to set a constraint, all events are conservatively assumed to originate from DM). Plausible origins of these events include the photodissociation of negatively charged impurities; ionized electrons that were initially created by highly ionizing background events, but then became trapped in the liquid-gas interface and spontaneously released at a later time; and field emission from the cathode [2,22]. However, more study is needed to understand and characterize these events.
In this letter, we derive new constraints from XENON10, including events with n e 4. The rate of such events is lower than for n e 3, leading to significantly improved constraints for DM masses m χ 50 MeV. We also analyze S2-only data from XENON100, containing n e 4 [23]. We derive the expected recoil spectra for the event rate and the annual and daily modulation amplitude, and show the expected event rates and the implications for a few benchmark DM models.
Other experiments, using semiconductor targets such as germanium (Ge) and silicon (Si), currently have a higher electron-recoil energy threshold and are thus less sensitive by several orders of magnitude than XENON10/100 [5]. Dramatic improvements in sensitiv-exposure [kg-yrs] fiducial mass [kg] XENON10 [22] 0.041 1.2 XENON100 [23] 29.8 48.3 LUX [27] 119 145 XENON1T [28] 2,000 1,000 LZ, XENONnT [28,29] 15,000 5,600 ity in the near future are likely with SuperCDMS [5,24], SENSEI [25], and possibly other experiments. Nevertheless, these experiments will initially have target masses of only O(1 kg), far less than current and future xenon experiments (Table I). Understanding the S2-only events in two-phase TPCs could thus lead to dramatic improvements in cross-section sensitivity and, as we will show, probe simple and predictive benchmark models. The large exposures will also allow for an annual modulation analysis [26], which can significantly improve upon the current limit even if the background rates are high. Theoretical Rates and Recoil Spectra. To calculate the DM-electron scattering rate in liquid xenon, we follow the procedure in [2] (see appendices for more details).
We treat the target electrons as single-particle states of an isolated atom, described by numerical RHF bound wave functions from [30] [31]. The velocity-averaged differential ionization cross section for electrons in the (n, l) shell is is the inverse mean speed for a given velocity distribution as a function of the minimum velocity, v min , required for scattering. We assume a standard Maxwell-Boltzmann velocity distribution with circular velocity v 0 = 220 km/s and a hard cutoff of v esc = 544 km/s [32,33]. σ e is the DM-free electron scattering cross section at fixed momentum transfer q = αm e , while the q-dependence of the matrix element is encoded in the DM form-factor F DM (q). |f nl ion (k , q)| 2 is the ionization form factor of an electron in the (n, l) shell with final momentum k = √ 2m e E er . We calculate this form factor using the given bound wave functions and unbound wave functions that are obtained by solving the Schrödinger equation with a potential that reproduces the bound wave functions. We consider electrons in the following shells (listed with binding energies in eV): 5p 6 (12.4), 5s 2 (25.7), 4d 10 (75.6), 4p 6 (163.5), and 4s 2 (213.8). The differential ionization rate is where N T is the number of target atoms and ρ χ = 0.4 GeV/cm 3 is the local DM density.
Top (bottom): Spectrum of expected number of events for DM-electron scattering in xenon, for mχ = 100 MeV and 1000 kg-years for FDM = 1 (α 2 m 2 e /q 2 ). For the left axes, we set σe to the maximum allowed values by current constraints for two popular benchmark models; for the right axes, the indicated σe produces the correct relic abundance. Colored lines show individual contributions from various xenon electron shells, while the gray band encompasses the spectrum when varying the secondary ionization model. See text for details.
We follow [2] to model the conversion from E er to electron yield, n e . The recoiling electron will ionize and excite other atoms, producing n (1) = Floor(E er /W ) additional "primary quanta", either observable electrons or (unobservable) scintillation photons. For fiducial values, we choose the probability for the initial electron to recombine with an ion to be f R = 0, W = 13.8 eV, and the fraction of primary quanta observed as electrons to be f e = 0.83. To capture the uncertainty in the fiducial values, we vary these parameters in the range 0 < f R < 0.2, 12.4 < W < 16 eV, and 0.62 < f e < 0.91. In addition to primary quanta, if DM ionizes an inner-shell electron, n (2) = Floor((E i − E j )/W ) secondary quanta can be created by photons produced in the subsequent outerto-inner-shell electron transitions with binding energies E i,j . The number of secondary electrons produced follows a binomial distribution with n (1) + n (2) trials and success probability f e .
In Fig. 1, we show the recoil spectra as a function of n e for a hypothetical xenon detector with 1000 kg-years of exposure for F DM = 1 (top) and F DM = α 2 m 2 e /q 2 (bottom). The colored lines show individual contributions from different shells, while the black line shows their sum (for fiducial values). Gray bands show the variation away from the fiducial values discussed above.
To emphasize the importance of studying electron re- Observed number of events versus photoelectrons (PE) in XENON10 (top) [22] and XENON100 (bottom) [23]. DM spectra are shown for mχ = 10 MeV (blue) & 1 GeV (red) with a cross section fixed at our derived 90% C.L. limit (we assume fiducial values for the secondary ionization model).
Insets show spectra in bins of 27PE (20PE), the mean number of PE created by one electron in XENON10 (XENON100).
coil events at current and upcoming xenon experiments, we have fixed σ e to specific values that are allowed by simple and predictive benchmark models [1,5,[34][35][36][37][38][39][40] and further below. We consider the DM (a Dirac fermion or complex scalar χ) to be charged under a broken U (1) D gauge force, mediated by a kinetically-mixed dark photon, A , with mass m A . The A mediates DM-electron scattering, and F DM (q) = 1 (α 2 m 2 e /q 2 ) for a heavy (ultralight) dark photon. The left axis for top (bottom) plot of Fig. 1 shows the event rate for σ e fixed to the maximum value allowed by current constraints for m A = 3m χ (m A keV), while the right axis of the top (bottom) plot fixes σ e so that scalar (fermion) DM obtains the correct relic abundance from thermal freeze-out (freeze-in). Clearly, a large number of DM events could be seen in upcoming detectors. These results are easily rescaled to other DM models that predict DM-electron scattering. New XENON10 and XENON100 bounds. We now recalculate the bounds from XENON10 data [2] (15 kg- 3. 90% C.L. limit on the DM-electron scattering cross section from XENON10 data (blue) and XENON100 data (red) for FDM = 1 (top) & FDM = α 2 m 2 e /q 2 (bottom). Dotted black lines show XENON10 bounds from [2]. days), including for the first time events with n e 4, as well as from XENON100 data [23] (30 kg-years). Since the experimental observable is the number of photoelectrons (PE) produced by an event, we convert n e to PE. An event with n e electrons produces a gaussian distributed number of PE with mean n e µ and width √ n e σ, where µ = 27 (19.7) and σ = 6.7 (6.2) for XENON10 (XENON100). We multiply the signal with the trigger and acceptance efficiencies from [2,23] and then bin both the signal and data in steps of 27PE (20PE), starting from 14PE (80PE) for XENON10 (XENON100). The first bin for the XENON100 analysis is 80-90PE, corresponding to roughly half an electron. We require that the resulting signal is less than the data at 90% C.L. in each bin. For XENON10, the 90% C.L. upper bounds on the rates (after unfolding the efficiencies) are r 1 < 15.18, r 2 < 3.37, r 3 < 0.95, r 4 < 0. 35 Fig. 2 shows the two data sets in PE and two sample DM spectra. Fig. 3 shows the strongest XENON10 and XENON100 limit combined across all bins, and a comparison with the XENON10 bound derived in [2]. In the Appendices, we show cross-section bounds for the individual PE bins, taking into account the systematic uncertainties from the secondary ionization model. For F DM = 1, the inclusion of the high-PE bins in XENON10 significantly improves upon the bound from [2] for m χ 50 MeV (small differences at lower masses are from the limit-setting procedure). The new XENON10 and XENON100 bounds are comparable for m χ 50 MeV. For F DM = α 2 m 2 e /q 2 , the low PE bins determine the bound, and XENON100 is therefore not competitive due to its high analysis threshold. Modulation. A useful discriminant between signal and background is the annual modulation of the signal rate [26] due to the Sun's motion through the DM halo. Fig. 4 shows f mod versus n e , where f mod = Rmax−Rmin 2Ravg is the modulation amplitude, derived by calculating the rates for the average Earth velocity and varying it by ±15.0 km/s. The f mod spectrum is distinctive, which should provide a helpful discriminant between signal and background. The significance of a signal S over a flat background B is then given by sig = f mod S √ S+B . To demonstrate the power of an annual modulation search, we imagine that a future detector with 1000 kgyears of exposure observes the same S2-only event rate and spectrum as observed in XENON10 data, R Xe10 . Requiring the signal rate to be less than the observed event rate yields the same constraints as with XENON10 data, σ e,Xe10 . However, an annual modulation analysis would potentially see a signal of high statistical significance, and in the absence of one a fraction of the observed event rate must be background. Requiring the significance of the annual modulation signal to be less than sig, the expected sensitivity is Sensitivity reach from an annual modulation analysis with a hypothetical 1000 kg detector and 1-year exposure, assuming the observed spectrum and data rate are the same as in XENON10 [22] (solid blue) or XENON100 [23] (solid red). DM-electron scattering event rates assuming a 1-electron (4electron) threshold are shown in dashed (dotted) green. Blue (red) shaded regions show our XENON10 (XENON100) limits. These lines/regions are overlaid on several simple and predictive benchmark models for DM (χ) scattering off electrons via a dark photon A . Top: (FDM = 1) A complex scalar obtains the correct relic density from thermal freezeout (light orange), while a fermion, which obtains its correct relic abundance from an initial asymmetry, must have σe above the dark brown line (assuming no additional annihilation channels) to avoid indirect-detection constraints [41][42][43]. Bottom: (FDM = α 2 m 2 e /q 2 ) Fermion DM coupled to an ultralight mediator A obtains the correct relic density from freeze-in (thick brown line). Gray regions show constraints as in [5], updated on the top plot with data from MiniBooNE [44] and BaBar [45]. Due to earth-scattering effects [46], no XENON10/100 limit exists in the top right region.
We calculate σ mod e for bins of n e = 0.5 − 1.5, 1.5 − 2.5, . . . and show with a blue line the best sensitivity across all bins in Fig. 5 for sig = 1.65 (90% CL) (see Appendices for sensitivities from each bin). Similarly, a red solid line shows σ mod e assuming the future observed rates/spectrum correspond to the current XENON100 rate/spectrum. We overlay these lines on the DM benchmark models discussed above. While hypothetical, this analysis emphasizes the power of an annual modulation analysis. Large Event Rates. To further emphasize the importance of understanding the electron recoil events in xenon TPCs, we show the expected event rates in Fig. 5 for a 1000 kg detector for two thresholds, n e ≥ 1 and n e ≥ 4. We see that the benchmark DM models predict large event rates. For example, Dirac fermion DM coupled to the A that obtains its abundance from an initial asymmetry could produce about one event every two seconds at LZ. This underscores the point that while there are several sources of backgrounds that can produce singleor few-electron events, a large event rate can be consistent with a DM signal and should not be simply written off as a detector curiosity. Conclusions. We derived new constraints on DMelectron scattering, improving upon the previous bound, and showed spectra for the expected number of electrons and the modulation amplitude. While there are several possible detector-specific origins of the observed XENON10/100 events, in principle almost all the observed events could originate from DM-electron scattering without coming into conflict with other existing DM constraints. This is not the case when interpreting these events as arising from few-GeV DM recoiling elastically off nuclei [22,23], which is excluded by existing results from e.g. LUX [27] and CDMSlite [47]. Moreover, simple and predictive DM benchmark models predict large event rates in current and future xenon TPCs. An expanded and dedicated effort by the xenon collaborations to understand the origin of their low-energy electron recoil data is thus imperative and well worth the effort.

ACKNOWLEDGMENTS
We would like to thank especially Aaron Manalaysay and Peter Sorensen for many insightful discussions. We also thank Ran Budnik, Daniel McKinsey, Matt Pyle, and Jingke Xu for useful discussions. We are also very grateful to Jeremy Mardon for contributions at the beginning of this project as well as many useful discussions.

APPENDIX
Here we provide additional details to the calculations described in the main text. We also show spectra plots for additional DM masses, as well as the XENON10/XENON100 limits and the prospects for an annual modulation analysis from each PE bin. For completeness, we also show the expected daily modulation of the signal rate due to the Earth's rotation.
Theoretical Rates. We first quote additional formulas that are required for the rate calculation (see also [2,5]). The velocity-averaged differential ionization cross section for electrons in the (n, l) shell is given in Eq. (1). The full expression for v min is where E nl binding is the binding energy of the shell and q is the momentum transfer from the DM to the electron. The form factor for ionization of an electron in the (n, l) shell with final momentum k = √ 2m e E er is given by where [· · · ] is the Wigner 3-j symbol and j L are the spherical Bessel functions. We solve for the radial wavefunctions R k l (r) of the outgoing unbound electrons taking the radial Schrödinger equation with a central potential Z eff (r)/r. This central potential is determined from the initial electron wavefunction by assuming that it is a bound state of the same potential. We include the shells listed in Table II.
Electron and Photoelectron Yields. We provide additional details to convert the recoiling electron's recoil energy into a specific number of electrons. The relevant  15   TABLE II. Xenon shells and energies. "Photon energy" refers to energy of de-excitation photons for outer-shell electrons deexciting to lower shells. This photon can subsequently photoionize, creating additional quanta. The range of additional quanta takes into account that the higher energy shell may have more than one available lower energy shell to de-excite into. For our limits, we take the minimum of this range. quantities are E er = (n γ + n e )W , E er is the amount of deposited energy from the primary electron, which results in a number of observable electrons, n e , unobservable scintillation photons, n γ , and heat. W is the energy needed to produce a single quanta (photon or electron). We take W = 13.8 ± 0.9 eV as the average energy [48]. E er can create both a number of ions, N i , and a number of excited atoms N ex , where N ex /N i 0.2 at energies above a keV [48,49]. We take into account the possibility that the primary electron and secondary ions can recombine, which is described by a modified Thomas-Imel recombination model [50], and assume that the fraction of ions that can recombine, f R , is effectively zero at low energy. This implies that n e = N i and n γ = N ex . The fraction of initial quanta observed as electrons is given by f e = (1 − f R )/(1 + N ex /N i ) 0.83 [51]. To capture the uncertainty in f R , W , and N ex /N i , we calculate the rates and limits varying these parameters over the ranges 0 < f R < 0.2, 0.1 < N ex /N i < 0.3, and 12.4 < W < 16 eV. For our fiducial values, we set f e = 0.83, f R = 0, W = 13.8 eV.
For each primary electron with energy E er , we assume that there are additional n (1) =Floor(E er /W ) quanta created. Furthermore, we assume that the photons associated with the de-excitation of the next-to-outer shells (5s, 4d, 4p, 4s), which have energies (13.3, 63.2, 87.9, 201.4) eV, can photoionize to create an additional n (2) =(n 5s , n 4d , n 4p , n 4s )=(0, 4, 6-10, 3-15) quanta, respectively (see Table II). The range in values for the 4p and 4s shells takes into consideration that there may be more than one outer-shell electron available that can de-excite down to them. For example, if the 4d shell de-excites to 4p, 6 additional quanta are created, while if the 5s shell de-excites to 4p, it would create 10 additional quanta. For our fiducial values, we take the lower number of quanta to be conservative. However, the choice of the number of additional quanta only affects n e > 6, and even here the difference in event rate is smaller than the uncertainties due to the modeling of the secondary ionization.
The total number of electrons is given by n e = n e +n e , where n e is the primary electron and n e are the secondary electrons produced. n e = 0 or 1 with probability f R or 1 − f R , respectively, while n e follows a binomial distribution with n (1) + n (2) trials and success probability f e . Given this conversion from E er into n e , we can calcu- late the differential rate as a function of number of electrons. In addition to the m χ = 100 MeV spectra shown in the main text, we show the spectra for m χ = 500 MeV and 1 GeV in Fig. 6.
XENON10 and XENON100 constraints for individual photoelectron bins. In the main text, we show the cross-section limits from the XENON10 and XENON100 data using the fiducial values above. In Figs. 7, 8, we show the individual limits for each PE bin as well as the uncertainty bands due to the secondary ionization model.
Modulation. In Fig. 5, we showed how an annual modulation analysis of a hypothetical xenon detector with an exposure of 1000 kg-years could significantly improve on current constraints even if the background rates are significant. In Fig. 5, we only showed the best constraints across all individual n e bins. In Fig. 9, we show the individual n e bins. Furthermore, for completeness, we also show the daily modulation amplitude due to the Earth's rotation with respect to the DM wind. The daily modulation is calculated by modifying the average earth velocity by ±0.23 km/s to obtain the maximum and minimum rates. We show the daily modulation fraction in Fig. 10, where we see that the daily modulation fraction is about an order of magnitude smaller than that of the annual modulation.