Excess electronic recoil events in XENON1T

We report results from searches for new physics with low-energy electronic recoil data recorded with the XENON1T detector. With an exposure of 0.65 tonne-years and an unprecedentedly low background rate of 76 (cid:1) 2 stat events = ð tonne × year × keV Þ between 1 and 30 keV, the data enable one of the most sensitive searches for solar axions, an enhanced neutrino magnetic moment using solar neutrinos, and bosonic dark matter. An excess over known backgrounds is observed at low energies and most prominent between 2 and 3 keV. The solar axion model has a 3 . 4 σ significance, and a three-dimensional 90% confidence surface is reported for axion couplings to electrons, photons, and nucleons. This surface is inscribed in the cuboid defined by g ae < 3 . 8 × 10 − 12 , g ae g effan < 4 . 8 × 10 − 18 , and g ae g a γ < 7 . 7 × 10 − 22 GeV − 1 , and excludes either g ae ¼ 0 or g ae g a γ ¼ g ae g effan ¼ 0 . The neutrino magnetic moment signal is similarly favored over background at 3 . 2 σ , and a confidence interval of μ ν ∈ ð 1 . 4 ; 2 . 9 Þ × 10 − 11 μ B (90% C.L.) is reported. Both results are in strong tension with stellar constraints. The excess can also be explained by β decays of tritium at 3 . 2 σ significance with a corresponding tritium concentration in xenon of ð 6 . 2 (cid:1) 2 . 0 Þ × 10 − 25 mol = mol. Such a trace amount can neither be confirmed nor excluded with current knowledge of its production and reduction mechanisms. The significances of the solar axion and neutrino magnetic moment hypotheses are decreased to 2 . 0 σ and 0 . 9 σ , respectively, if an unconstrained tritium component is included in the fitting. With respect to bosonic dark matter, the excess favors a monoenergetic peak at ð 2 . 3 (cid:1) 0 . 2 Þ keV (68% C.L.) with a 3 . 0 σ global ( 4 . 0 σ local) significance over background. This analysis sets the most restrictive direct constraints to date on pseudoscalar and vector bosonic dark matter for most masses between 1 and 210 keV = c 2 . We also consider the possibility that 37 Ar may be present in the detector, yielding a 2.82 keV peak from electron capture. Contrary to tritium, the 37 Ar concentration can be tightly constrained and is found to be negligible. DOI: 10.1103/PhysRevD.102.072004


I. INTRODUCTION
A preponderance of astrophysical and cosmological evidence suggests that most of the matter content in the Universe is made up of a rarely interacting, nonluminous component called dark matter [1].Although several hypothetical dark matter particle candidates have been proposed with an assortment of couplings, masses, and detection signatures, dark matter has thus far eluded direct detection.The XENON1T experiment [2], employing a liquid-xenon time projection chamber (LXe TPC), was primarily designed to detect weakly interacting massive particle (WIMP) dark matter.Due to its unprecedentedly low background rate, large target mass, and low-energy threshold, XENON1T is also sensitive to interactions from alternative dark matter candidates and to other physics beyond the Standard Model (SM).Here we report on searches for (1) axions produced in the Sun, (2) an enhancement of the neutrino magnetic moment using solar neutrinos, and (3) pseudoscalar and vector bosonic dark matter, including axion-like particles (ALPs) and dark photons.
The XENON1T experiment operated underground at the INFN Laboratori Nazionali del Gran Sasso (LNGS) from 2016 to 2018, utilizing a dual-phase LXe TPC with a 2.0tonne active target to search for rare processes.A particle interaction within the detector produces both prompt scintillation (S1) and delayed electroluminescence (S2) signals.These light signals are detected by arrays of photomultiplier tubes (PMTs) on the top and bottom of the active volume, and are used to determine the deposited energy and interaction position of an event.The latter allows for removing background events near the edges of the target volume (e.g., from radioactivity in detector materials) through fiducialization.The S2/S1 ratio is used to distinguish electronic recoils (ERs), produced by, e.g., gamma rays (γs) or beta electrons (βs), from nuclear recoils (NRs), produced by, e.g., neutrons or WIMPs, allowing for a degree of particle identification.The ability to determine scatter multiplicity enables further reduction of backgrounds, as signals are expected to have only single energy deposition.
In this paper, we report on searches for ER signals with data acquired from February 2017 to February 2018, a time period referred to as Science Run 1 (SR1) [3].As the vast majority of background comes from ER events, we search for excesses above a known background level.The analysis is carried out in the space of reconstructed energy, which exploits the anticorrelation of S1 and S2 signals by combining them into a single energy scale [4], thus reducing the statistical fluctuations from electron-ion recombination [5].Both S1 and S2 signals are corrected to disentangle position-dependent effects, such as light collection efficiency (LCE) and electron attachment to electronegative impurities.After correcting to the mean LCE across the TPC, S1 is reconstructed using signals from all PMTs (cS1).For the S2 reconstruction, only the bottom PMT array is used ( cS2 b ) because it features a more homogeneous light collection [3].The full energy region of interest (ROI) is (1,210) keV, which is primarily motivated by the search for bosonic dark matter and discussed further in Sec.III A.
The paper is organized as follows.In Sec.II, we present the theoretical background and signal modeling of the beyond-the-SM channels considered in this search.We describe the data analysis in Sec.III, including the data selection, background model, and statistical framework.In Sec.IV, upon observation of a low-energy excess in the data, we present a hypothesis of a new background component, tritium, which may be observable for the first time in a xenon detector due to our unprecedented low background.We then report the results of searches for solar axions, an anomalous neutrino magnetic moment, and bosonic dark matter.We end with further discussion of our findings and a summary of this work in Secs.V and VI, respectively.The presence of the excess motivated further scrutiny of the modeling of dominant backgrounds, the details of which we present in the Appendix.

II. SIGNAL MODELS
This section describes the physics channels we search for in this work.In Sec.II A, we motivate the search of solar axions, presenting their production mechanisms in the Sun and the detection mechanism in LXe TPCs, and summarize two benchmark axion models.In Sec.II B, we introduce the search for an anomalous neutrino magnetic moment, which would enhance the neutrino-electron elastic scattering cross section at low energies.In Sec.II C, we discuss the signals induced by bosonic dark matter including pseudoscalar and vector bosons, examples of which are axion-like particles and dark photons, respectively.Expected energy spectra of these signals in the XENON1T detector are summarized at the end of this section.
For all signal models presented below, the theoretical energy spectra in a LXe TPC were converted to the space of reconstructed energy by accounting for detector efficiency and resolution, summarized in Fig. 1.The efficiency is shown in Fig. 2 and discussed in Sec.III A. For the energy resolution, the theoretical spectra were smeared using a Gaussian distribution with energy-dependent width, which was determined using an empirical fit of monoenergetic peaks as described in [2,4].The energy resolution σ is given by with a ¼ ð0.310 AE 0.004Þ ffiffiffiffiffiffiffiffi keV p and b ¼ 0.0037 AE 0.0003.When building a signal model, the resolution is first applied to the deposited, "true" energy spectrum, and then the smeared distribution is corrected according to the predicted loss due to efficiency.This implies that, near threshold, the mean reconstructed energy is higher than the true energy, as the reduced efficiency at lower energies shifts the mean of the observed distribution upwards.This type of reconstruction bias is fully accounted for in this analysis.

A. Solar axions
As a solution to the strong CP problem in quantum chromodynamics (QCD), Peccei and Quinn postulated a mechanism that naturally gives rise to a Nambu-Goldstone boson, the so-called axion [6][7][8].In addition to solving the strong CP problem, QCD axions are also well-motivated dark matter candidates, with cosmological and astrophysical bounds requiring their mass to be small (typically ≪keV) [9][10][11][12][13].On account of this mass constraint, dark matter axions produced in the early Universe cannot be observed in XENON1T.However, solar axions would emerge with-and in turn deposit-energies in the keV range [14][15][16], the precise energies to which XENON1T was designed to be most sensitive.An observation of solar axions would be evidence of beyond-the-SM physics, but would not by itself be sufficient to draw conclusions about axionic dark matter.
We consider three production mechanisms that contribute to the total solar axion flux: (i) atomic recombination and deexcitation, Bremsstrahlung, and Compton (ABC) interactions [14,17], (ii) a monoenergetic 14.4 keV M1 nuclear transition of 57 Fe [15], and (iii) the Primakoff conversion of photons to axions in the Sun [18,19].The ABC flux scales with the axion-electron coupling g ae as and was taken from [14].The 57 Fe flux scales with an effective axion-nucleon coupling g eff an ¼ −1.19g 0 an þ g 3 an and is given by [20,21] Φ 57 Fe a where g 0=3 an are the isoscalar/isovector coupling constants and k a and k γ are the momenta of the produced axion and photon, respectively.The Primakoff flux scales with the axion-photon coupling g aγ and is given by [22] dΦ Prim a dE a ¼ g aγ GeV −1

2
E a keV 2.481 e −E a =ð1.205 keVÞ where E a is the energy of the axion.All three flux components could be detected in XENON1T via the axioelectric effect-the axion analog to the photoelectric effect-which has a cross section that scales with axionelectron coupling g ae and is given by [20,[23][24][25] where β and E a are the velocity and energy of the axion, respectively, α is the fine structure constant, and m e is the mass of the electron.The energy-dependent photoelectric cross section, σ pe , was obtained from [26] and interpolated between points using the logarithms of both photon energies and cross sections.Combining the production and detection mechanisms, we are able to constrain the values of jg ae j (ABC), jg ae g eff an j ( 57 Fe), and jg ae g aγ j (Primakoff). 1 We consider these three observables independently in the analysis, lest we implicitly assume any particular axion model.Still, it FIG.2. Efficiency as a function of energy.The dashed (dotted) line refers to detection (selection) efficiency, while the blue curve and band illustrate the total efficiency and the associated 1-σ uncertainty, respectively.The detection threshold is indicated by the right bound of the gray shaded region.FIG. 1. Left: expected signal in energy space for ABC solar axions with a coupling g ae ¼ 5 × 10 −12 (blue), for solar axions produced from the deexcitation of 57 Fe with coupling g eff an ¼ 1 × 10 −6 (red), and for solar axions produced from the Primakoff effect with coupling g aγ ¼ 2 × 10 −10 (orange).Right: signature of an enhanced neutrino magnetic moment with magnitude 7 × 10 −11 μ B (green) and a 20 keV=c 2 ALP with coupling constant g ae ¼ 2 × 10 −13 (purple).Both the true deposited energy spectra in a xenon detector without efficiency loss (unshaded) and the expected observed spectra in XENON1T including the specific detector resolution and efficiency (shaded) are shown.is important to note that these values are indeed related to each other and to the axion mass under different models.
For QCD axions, the mass m a is related to the decay constant f a via m a ≃ 6 × 10 6 GeV f a eV=c 2 ; ð6Þ and the axion couplings to matter are mostly model dependent.We describe here two benchmark classes of QCD axion models: Dine-Fischler-Srednicki-Zhitnitsky (DFSZ) [27,28], in which axions couple to electrons at tree level, and Kim-Shifman-Vainshtein-Zhakharov (KSVZ) [29,30], where couplings to leptons occur only at loop level.For this reason, the ABC flux is dominant in DFSZ models, while the Primakoff flux is dominant in KSVZ models.Since the axioelectric cross section scales with the axion-electron coupling, XENON1T is in general more sensitive to DFSZ-type axions.
In DFSZ models, the axion-electron coupling is given by where and X u and X d are the Peccei-Quinn (PQ) charges of the up and down quarks, respectively [20,31,32].The couplings to quarks take on a similar expression with respect to β DFSZ .The axion-nucleon couplings g 0 an and g 3 an are functions of X u , X d , and f a , and can be found in [31,33].For a DFSZ axion, it follows that g ae and g eff an are both nonzero in general, as they are connected via β DFSZ and f a .The axionphoton coupling does not depend on the PQ charges but is directly related to the axion decay constant (and thus the mass), where z ¼ m u =m d , m u=d are the respective masses of the up/down quarks, and E=N represent the model-dependent electromagnetic/color anomalies of the axial current associated with the axion field [34].It is typically assumed that E=N ¼ 8=3 in DFSZ models.
In KSVZ models, the PQ charges of the SM quarks vanish, and there is no β DFSZ -like parameter.The axionelectron coupling strength, induced by radiative corrections, depends on the axial current [31,35], where w ¼ m u =m s , m s is the mass of the strange quark and Λ is the cutoff of the QCD confinement scale.The isoscalar/isovector axion-nucleon couplings g 0=3 an do not depend on the PQ charges and are also found in [31,33].The axion-photon coupling is given by Eq. ( 9).For KSVZ models, a benchmark value of E=N ¼ 0 is often used, but many values are possible [36].
As mentioned above, no particular axion model is assumed in the analysis itself; the three flux components are considered completely independent of each other.Since, in principle, it is possible for all three components to be present at the same time, our solar axion model includes three unconstrained parameters for the different components.Were a signal observed, the results of the threecomponent analysis could then be used to constrain different axion models and possibly infer the axion mass.This approach also implies that the results hold generally for solar axion-like particles, which do not have strict relationships between the couplings, as described in Sec.II C.
The expected spectra from solar axions with g ae ¼ 5 × 10 −12 , g aγ ¼ 2 × 10 −10 GeV −1 , and g eff an ¼ 10 −6 are shown in Fig. 1 (left) with before/after detector effects indicated by unshaded/shaded curves, respectively.The rate of the ABC component is proportional to g ae 4 , the 57 Fe component is proportional to ðg ae g eff an Þ 2 , and the Primakoff component is proportional to ðg ae g aγ Þ 2 .

B. Neutrino magnetic moment
In the SM, neutrinos are massless and therefore without a magnetic dipole moment.However, the observation of neutrino oscillation tells us that neutrinos have mass and the SM must be extended, thus implying a magnetic moment of μ ν ∼ 10 −20 μ B [37][38][39][40], where μ B is the Bohr magneton.Larger values of μ ν have been considered theoretically and experimentally [40][41][42].Interestingly, in addition to providing evidence of beyond-SM physics, the observation of a μ ν ≳ 10 −15 μ B would suggest that neutrinos are Majorana fermions [40].Currently, the most stringent direct detection limit is μ ν < 2.8 × 10 −11 μ B from Borexino [42], and indirect constraints based on the cooling of globular cluster and white dwarfs are an order of magnitude stronger at ∼10 −12 μ B [32,43,44].
An enhanced magnetic moment would increase the neutrino scattering cross sections at low energies (on both electrons and nuclei) and thus could be observable by lowthreshold detectors such as XENON1T.Here we only consider the enhancement to elastic scattering on electrons, given by [45] where E r is the electronic recoil energy and E ν is the energy of the neutrino.Note that Eq. (11)  small corrections need to be made for the electron binding energies at Oð keVÞ energies.We search for an anomalous magnetic moment using solar neutrinos, predominantly those from the protonproton (pp) reaction [46].The expected energy spectrum for μ ν ¼ 7 × 10 −11 μ B is shown in Fig. 1 (right), which was calculated by folding the expected solar neutrino flux [46] with Eq. ( 11) and applying a step-function approximation to account for the electron binding energies.In the energy range considered here, this approximation agrees well with more detailed calculations [47].Note that this signal would be added to the SM neutrino elastic scattering spectrum, which we treat as a background as described in Sec.III B.

C. Bosonic dark matter
ALPs, like QCD axions, are pseudoscalar bosons, but with decay constant and particle mass [Eq.( 6)] decoupled from each other and instead taken as two independent parameters.This decoupling allows for ALPs to take on higher masses than QCD axions; however, it also implies that ALPs do not solve the strong CP problem.
ALPs are viable dark matter candidates [48] and could be absorbed in XENON1T via the axioelectric effect [Eq. ( 5)] like their QCD counterparts.Assuming ALPs are nonrelativistic and make up all of the local dark matter (density ρ ∼ 0.3 GeV=cm 3 [49]), the expected signal is a monoenergetic peak at the rest mass of the particle, m a , with an event rate given by (see [25,50]) where A is the average atomic mass of the detector medium (A ≈ 131 u for xenon).The rate coefficient from our calculation is consistent with [51] for the dark matter density used in this work.
In addition to the pseudoscalar ALPs, XENON1T is also sensitive to vector bosonic dark matter, of which dark photons are a common example.Dark photons can couple weakly with SM photons through kinetic mixing [52] and be absorbed with cross section σ V given by [53] where σ pe , α, and β are the same as in Eq. ( 5), and κ parametrizes the strength of kinetic mixing between the photon and dark photon.Similarly to Eq. ( 12), by following the calculation in [25], the rate for nonrelativistic dark photons in a detector reduces to where m V is the rest mass of the vector boson.Like the pseudoscalar above, absorption of a vector boson would also result in a monoenergetic peak broadened by the energy resolution of the detector, but with a rate that is inversely proportional to the particle mass.The expected spectrum for a 20 keV=c 2 ALP with g ae ¼ 2 × 10 −13 is shown in Fig. 1 (right).Vector bosons have the same signature as ALPs, but the rate scales differently with mass [see Eqs. ( 12), ( 14)].

III. DATA ANALYSIS
This section describes the data-analysis methods employed to search for the aforementioned signals.The event-selection criteria and their overall efficiency, the detection efficiency, as well as the determination of fiducialization and ROI are given in Sec.III A. Section III B details each component of our background model, the predictions of which are consistent with the results of a background-only fit to the data.In Sec.III C, we define the likelihood used for the fitting and discuss the statistical framework.

A. Data selection
The data-selection criteria for this search are similar to [3], with the selections and efficiencies optimized and reevaluated for the different parameter space and extended energy range.For an event to be considered valid, an S1-S2 pair is required.A valid S1 demands coincident signals in at least three PMTs, and a 500 photoelectron (PE) threshold is imposed on the S2 size.This S2 threshold is more stringent than that in [3] in order to reject background events originating from radon daughters on the TPC surface [54].Since signal events are expected to deposit energy only once in the detector, events with multiple interaction sites are removed.A variety of selection criteria are applied to ensure data quality and a correct S1 and S2 pairing, which is detailed in [55].The efficiencies and uncertainties of the selection criteria are estimated in a procedure similar to [55], and the cumulative selection efficiency is determined using an empirical fit of the data.The average cumulative selection efficiency over the (1,210) keV region is ð91.2AE 0.3Þ%.
The combined efficiency of detection and event selection with uncertainties is shown in Fig. 2. The detection efficiency, dominated by the threefold coincidence requirement of S1s, was estimated using both a data-driven method of sampling PMT hits from S1s in the 20-100 PE range and an independent study based on simulation of low-energy S1 waveforms [55].The difference between the two methods (∼3% average relative difference in the dropoff region) was considered as a systematic uncertainty.This efficiency was then converted from S1 to reconstructed energy using the detector-response model described in [54], accounting for additional uncertainties such as the photon yield.The S2 efficiency can be assumed to be unity for the energies considered here [55].
Events with energies between (1, 210) keV are selected for this search, with the lower bound determined by requiring the total efficiency be larger than 10%, shown in Fig. 2, and the upper bound limited by an increasing γ-ray background from detector materials, which is difficult to model due to large uncertainties on its spectral shape.While extending the ROI to 210 keV is primarily motivated by the bosonic dark matter search, all profile likelihood fits use this full energy range, as it also allows for better constraints on the background components.The same 1042 kg cylindrical fiducial volume as in [56] was used to reduce the surface and material backgrounds.After event selection and strict fiducialization, the surface backgrounds, accidental coincidences, and neutrons make up less than 0.003% of the total events (<0.3% below 7 keV) and thus are negligible for this search.Additionally, events within 24 hours from the end of calibration campaigns using injected radioactive sources were removed due to residual source activity.The final effective SR1 live time is 226.9 days and thus the total exposure is 0.65 tonne-years.

B. Background model
Within the (1,210) keV ROI and the 1042 kg fiducial volume, ten different components were used to model the background and fit the data, as listed in Table I and illustrated in Fig. 3.
Six components, numbers i-vi in Table I, exhibit continuous energy spectra and were modeled based on either theoretical predictions or GEANT4 Monte Carlo simulations, and the rest are monoenergetic peaks that were modeled as Gaussian functions of known energies and resolution.The spectrum of each background component considers the detector energy resolution and efficiency loss in the same way as the signal model construction in Sec.II.The rates of the background components are  85 Kr, which is subdominant due to its removal via cryogenic distillation [57,58].The shape of these spectra, particularly at low energies, can be affected by atomic screening and exchange effects, as well as by nuclear structure [59,60].The β decays of 214 Pb and 85 Kr are first forbidden nonunique and first forbidden unique transitions, respectively; however, spectra from the IAEA LiveChart (Nuclear Data Services database) [61] are based on calculations of allowed and forbidden unique transitions, neither of which includes exchange effects [62].Likewise, models from GEANT4 [63] include only the screening effect; however, its implementation displays a nonphysical discontinuity at low energies [62,64].For this work, we performed dedicated theoretical calculations to account for possible low-energy discrepancies from these effects in 214 Pb and 85 Kr spectra.These calculations are described in detail in the Appendix.
The activity of 214 Pb can be constrained using in situ measurements of other nuclei in the same decay chain.These constraints, described in [54], place a lower bound of 5.1 AE 0.5 μBq=kg (from coincident 214 BiPo) and upper bound of 12.6 AE 0.8 μBq=kg ( 218 Po α decays).For this analysis, we leave the normalization of the 214 Pb rate unconstrained and use the fit to extract the activity.The background-only fit results give an event rate of 63.0 AE 1.3 events=ðtonne × year × keVÞ [abbreviated as events=ðt • y • keVÞ for the rest of paper] over the ROI after efficiency correction.With the 11% branching ratio (from [65]) and the spectrum of 214 Pb decay to the ground state (calculated in the Appendix), the 214 Pb activity is evaluated to be ð11.1 AE 0.2 stats AE 1.1 sys Þ μBq=kg throughout SR1 and is well within the upper/lower bounds.The 10% systematic uncertainty is mainly from the aforementioned branching ratio [65].
The 85 Kr decay rate is inferred from dedicated measurements of the isotopic abundance of 85 Kr= nat Kr (2 × 10 −11 mol=mol) and the nat Kr concentration evolution in LXe [66].The same measurements also allow for the time dependence of the 85 Kr decay rate to be taken into account.The average rate of 85 Kr is 7.4 AE 1.3 events=ðt • y • keVÞ over the ROI in SR1.
An additional background arises from γ emissions from radioimpurities in detector materials that induce Comptonscattered electrons; however, this background is subdominant in the ROI due to the strict fiducial volume selection.The rate from materials is constrained by radioassay measurements [67] and predicted by simulations [68] to be 2.7 AE 0.3 events=ðt • y • keVÞ.This background is modeled by a fixed, flat component in the fit.
One of the continuous backgrounds considered was 136 Xe, a 2νββ emitter intrinsic to xenon.This component has an increasing rate as a function of energy over the ROI.It was constrained in the fit according to the predicted rate and associated uncertainties on (i) a 136 Xe isotopic abundance of ð8.49AE 0.04 stat AE 0.13 sys Þ% as measured by a residual gas analyzer [69], (ii) the reported half-life [70], and (iii) the calculated theoretical spectrum [71,72].
The first observation of two-neutrino double electron capture (2νECEC) of 124 Xe was recently reported using mostly the same SR1 data set (but different selection cuts) as used in this analysis [73] and is treated as a background here.In [73], we considered the dominant branching ratio of 2νECEC, the capture of two K-shell electrons inducing a peak at 64.3 keV.It is also possible to capture a K-shell and L-shell electron (36.7 keV) or two L-shell electrons (9.8 keV) with decreasing probability, as calculated in [74].For this analysis, the event selection and consideration of time dependence allow us to include all three peaks in the background model.The predicted rates of the peaks are taken from an updated half-life [75] with fixed branching ratios from [74]; the overall rate was not constrained in the fit since the half-life was derived from the same data set.
Three additional backgrounds were included for neutronactivated isotopes: 133 Xe (β), 131m Xe [internal conversion (IC)], and 125 I [electron capture (EC)].These isotopes were produced during neutron calibrations and decayed away with half-lives of Oð10Þ days.The IC decay of 131m Xe produces a monoenergetic peak at 164 keV [76], which, along with the other monoenergetic backgrounds, has the same signature as a bosonic dark matter signal.It was well constrained using its half-life and known dates of neutron calibration. 133Xe decays to an excited state with a dominant branching ratio and emits an 81 keV prompt γ upon deexcitation [77], resulting in a continuous spectrum starting at ∼75 keV, given the energy resolution.The rate was also constrained in the fit with prediction obtained using time dependence.The third activated isotope 125 I, a daughter of 125 Xe, decays via EC of K-shell, L-shell, and M-shell with decreasing probability and produces peaks at 67.3, 40.4,and 36.5 keV, respectively [78].Similar to 124 Xe 2νECEC, all three peaks of 125 I EC are included in the background model with the fixed branching ratios from [78].The 125 I contribution was constrained using a model based on the time evolution of 125 Xe throughout SR1, as detailed in [73].
During SR1, a background from 83m Kr (IC) was present due to a trace amount of 83 Rb (EC, T 1=2 ∼ 86 days) in the xenon recirculation system, which presumably was caused by a momentary malfunction of the source valve and confirmed using half-life measurements.83m Kr decays via a two-step scheme (second step T 1=2 ∼ 154 ns) [79] resulting in many of these events being removed by the multisite selections mentioned in Sec.III A; however, due to the short half-life of the second step, these decays are often unresolved in time and hence contribute as a monoenergetic peak at 41.5 keV.This component was also constrained using a time-evolution model.
Elastic scattering of solar neutrinos off electrons is expected to contribute subdominantly over the entire ROI.The expected energy spectrum was obtained using the standard neutrino flux in the Large Mixing Angle Mikheyev-Smirnov-Wolfenstein model and cross section given by the SM [46,80].Based on rate calculations of neutrino-electron scattering in xenon as given in [81], a 3% uncertainty was assigned and used to constrain the solar neutrino rate in the fit.
We denote the background model described above as B 0 .This model was used to fit the SR1 data in (1,210) keV by maximizing the likelihood constructed in Sec.III C. The fit results are consistent with predictions, as summarized in Table I.The best fit of B 0 is shown in Fig. 3, where the top panel is the full SR1 data set and the bottom two panels are partitions of SR1, which were fit simultaneously to include the temporal information of several backgrounds (see Sec. III C).This fit gives a background rate of 76 AE 2 events=ðt • y • keVÞ within the (1, 30) keV region after efficiency correction with the associated uncertainty from the fitting.Figure 4 shows a zoom in (0, 30) keV region of Fig. 3 with a finer binning.
In Sec.IV, we raise the possibility of an additional background component, the β decay of tritium, that we did not include while constructing the background model.A validated β-decay spectrum from the IAEA LiveChart [61,82] was used for the 3 H model, as described in the Appendix.We treat the possible tritium contribution separately from B 0 for reasons discussed in Sec.IVA.

C. Statistical method
An unbinned profile likelihood method is employed in this analysis.The likelihood is constructed as where μ s and μ b are the expected total signal and background events.Both μ b and θ are nuisance parameters, where θ includes shape parameters for the efficiency spectral uncertainty (see Fig. 2), as well as peak location uncertainties, specifically for 124 Xe (3 peaks), 83m Kr, and 131m Xe. Having largely subdominant event rates, the three peak locations from 125 I EC are fixed at their expected positions to save computation time.Index i runs over all observed events with the total number of N (42251 events), and E i corresponds to the energy of the ith event.f b and f s are the background and signal probability distribution functions, and index j runs over all the background components.C μ and C θ are constraints on the expected numbers of background events and the shape parameters.Index m runs over backgrounds including 85 Kr, solar neutrino, 136 Xe, 83m Kr, 125 I, 133 Xe, and 131m Xe, while index n is for all six shape parameters.
Due to time-dependent backgrounds, the SR1 data set is divided into two partitions: SR1 a consisting of events within 50 days following the end of neutron calibrations and SR1 b containing the rest, with effective live times of 55.8 and 171.2 days, respectively.Including this time information allows for better constraints on the timeindependent backgrounds and improves sensitivity to bosonic dark matter, especially as the time-dependent background from 133 Xe impacts a large fraction of its search region.The full likelihood is then given by where L a and L b are evaluated using Eq. ( 15) in each partition.Nuisance parameters that do not change with time, along with all of the signal parameters, are shared between the two partitions.The constant nuisance parameters are as follows: (i)  In the following sections, this excess is interpreted under solar axion, neutrino magnetic moment, and tritium hypotheses.
(iii) The solar neutrino rate, which would vary by ∼3% between the two partitions on account of Earth's orbit around the Sun.This is ignored due to the subdominant contribution from this source.(iv) The decay rates of the intrinsic xenon isotopes 136 Xe and 124 Xe, as well as the Compton continuum from materials.The remaining parameters all display time dependencies that are modeled in the two partitions.
The test statistic used for the inference is defined as where ðμ s ; μb ; θÞ is the overall set of signal and nuisance parameters that maximize L, while Lðμ s ; μb ; θÞ is maximized L by profiling nuisance parameters with a specified signal parameter μ s .The statistical significance of a potential signal is determined by qð0Þ.For the neutrino magnetic moment and bosonic dark matter searches, a modified Feldman-Cousins method in [83] was adopted in order to derive 90% C.L. bounds with the right coverage.
We report an interval instead of an upper limit if the global significance exceeds 3σ.For bosonic dark matter, this corresponds to 4σ local significance on account of the lookelsewhere effect, which is not present for the neutrino magnetic moment search.The 3σ significance threshold only serves as the transition point between reporting oneand two-sided intervals and was decided prior to the analysis to ensure correct coverage.A two-sided interval does not necessarily indicate a discovery, which in particle physics generally demands a 5σ significance and the absence of compelling alternate explanations.Since the solar axion search is done in the space of g ae , g ae g aγ , and g ae g eff an , we extend its statistical analysis to three dimensions.For this search, we use a standard profile likelihood construction where the true 90th percentile of the test statistic [Eq.(17)] was evaluated at several points on a three-dimensional (3D) grid and interpolated between points to define a 3D "critical" volume of true 90-percent threshold values.By construction, the intersection of this volume with the test statistic qðg ae ; g ae g aγ ; g ae g an Þ defines a three-dimensional 90% C.L. volume in the space of the three axion parameters.In Sec.IV, we report the twodimensional projections of this volume, found by profiling over the third respective signal component.

IV. RESULTS
When compared to the background model B 0 , the data display an excess at low energies, as shown in Fig. 4. The excess departs slightly from the background model near 7 keV, rises with decreasing energy with a peak near 2-3 keV, and then subsides to within AE1σ of the background model near 1-2 keV.Within this reference region of 1-7 keV, there are 285 events observed in the data compared to an expected 232 AE 15 events from the background-only fit, a 3.3σ Poissonian fluctuation.Events in this energy region are uniformly distributed in the fiducial volume.The temporal distribution of these events is discussed in Sec.IV E.
Several instrumental backgrounds and systematic effects were excluded as possible sources of the excess.Accidental coincidences (AC), an artificial background from detector effects, are expected to be spatially uniform, but are tightly constrained to have a rate of <1 event=ðt • y • keVÞ based on the rates of lone signals in the detector, i.e., S1s (S2s) that do not have a corresponding S2 (S1) [54].Surface backgrounds have a strong spatial dependence [54] and are removed by the fiducialization (1.0 tonne here vs 1.3 tonnes in [3], corresponding to a radial distance from the TPC surface of ≳11 cm) along with the stricter S2 threshold cut.Both of these backgrounds also have well-understood signatures in the (cS1, cS2 b ) parameter space that are not observed here, as shown in Fig. 5.
The detection and selection efficiencies were verified using 220 Rn calibration data.The β decay of 212 Pb, a daughter of 220 Rn, was used to calibrate the ER response of the detector, and thus allows us to validate the efficiency modeling with a high-statistics data set.Similarly to 214 Pb, the model for 212 Pb was calculated to account for atomic screening and exchange effects, as detailed in the Appendix.A fit to the 220 Rn data with this model and the efficiency parameter described in Sec.III C is shown in Fig. 6 for a 1-tonne fiducial volume, where good agreement is observed (p-value ¼ 0.50).Additionally, the S1 and S2 signals of the low-energy events in background data were found to be consistent with this 220 Rn data set, as shown in FIG. 5. Distribution of low-energy events (black dots) in the (cS1, cS2 b ) parameter space, along with the expected surface (purple) and AC (orange) backgrounds (1σ band). 220Rn calibration events are also shown (density map).All the distributions are within the 1-tonne fiducial volume.Gray lines show isoenergy contours for electronic recoils, where 1 and 7 keV contours, the boundaries of the reference region, are highlighted in blue.Fig. 5.This discounts threshold effects and other mismodeling (e.g., energy reconstruction) as possible causes for the excess observed in Fig. 4.
Uncertainties in the theoretical background models were considered, particularly for the dominant 214 Pb component.More details can be found in the Appendix, but we briefly summarize them here.A steep rise in the spectrum at low energies could potentially be caused by exchange effects in β-decay emission; however, this component is accurate to within 1% and therefore negligible with respect to the observed excess.The remaining two components, namely, the endpoint energy and nuclear structure, tend to shift the entire β distribution, rather than cause steep changes over a range of ∼10 keV.Conservatively, the combined uncertainty from these two components is þ6% in the 1-10 keV region, as described in the Appendix.In comparison, a þ50% uncertainty at 2-3 keV on the calculated 214 Pb spectrum, as constrained by the higher energy component, would be needed to make up the excess.
We also considered backgrounds that might in principle be present in trace amounts.First, low-energy x-rays from 127 Xe EC, as seen in [84,85], are ruled out for a number of reasons. 127Xe is produced from cosmogenic activation at sea level; given the short half-life of 36.4 days and the fact that the xenon gas was underground for OðyearsÞ before the operation of XENON1T, 127 Xe would have decayed to a negligible level.Indeed, high-energy γs that accompany these x-rays were not observed, and with their OðcmÞ mean free path in LXe they could not have left the OðmÞ-sized TPC undetected.For these reasons, we conclude that 127 Xe was no longer present during SR1.
Another potential background is 37 Ar, which decays via EC to the ground state of 37 Cl, yielding a 2.82 keV peak with a 0.90 branching ratio [86].It was considered by the LUX Collaboration as a background to explain a possible excess rate at ∼3 keV in their data [87].Its ingress was hypothesized to come from either from an initial amount in the xenon gas or from an air leak during operations; however, no definitive conclusion was drawn based on measurements of both the leakage rate and the 37 Ar concentration in air at the experimental site [88].We consider the two aforementioned possibilities for the introduction of 37 Ar into the xenon target and place quantitative constraints on each source.

37
Ar has a half-life of T 1=2 ¼ 35.0 days [86] and a typical abundance in nat Ar of ∼10 −20 mol=mol [89].Given an initial measured nat Ar concentration of <5 ppm in the xenon inventory [90], 37 Ar decayed to a negligible level, <1 events=ðt • yÞ, by the start of the XENON1T commissioning phase (>400 days).As with krypton, argon is not removed by the getter in the purification system, although it is removed by online 85 Kr distillation (see Sec. III B).This further suppresses its presence prior to SR1.These factors conclusively rule out the presence of 37 Ar from its initial concentration in the xenon inventory.
With respect to an 37 Ar component from a constant air leak, the similarities between krypton and argon noble gases allow us to use nat Kr to constrain the concentration of 37 Ar in the detector.From frequent measurements using rare gas mass spectrometry [66] and its natural abundance [91], the observed increasing concentration of nat Kr of <1 ppt=year gives an upper limit on the leak rate of ≤0.9 liter=year during SR1, following online distillation [92].We make a conservative assumption that the nat Kr increase is due entirely to a leak (neglecting emanation).
The air inside the experimental hall at LNGS, supplied from outside of the laboratory and fully exchanged within 2.5 hours, has an 37 Ar concentration of <3.2 mBq=m 3 , as determined from measurements taken in July 2020 following the methods in [93,94].We set a constraint using a robust upper limit of 5 mBq=m 3 for the 37 Ar equilibrium concentration to account for possible seasonal variations [95,96].The estimate is further refined after considering the differential leak rates of the two noble gases based on their respective viscosities in air, as well as accounting for the relative volatility of argon in liquid/gaseous xenon.Applying these corrections and conservatively assuming that 37 Ar reached an equilibrium activity by the start of SR1, we find that its expected rate is <5.2 events=ðt • yÞ.To explain the excess in XENON1T, the 37 Ar rate is required to be ∼65 events=ðt • yÞ, implying that the deduced upper limit is a factor of 13 too low to account for the excess.This conservative constraint on its presence in SR1 therefore excludes 37 Ar from a constant air leak as an explanation for the excess.The time dependence of a potential 37 Ar background is discussed further in Sec.IV E; however, no clear trend is observed due to low statistics, and any temporal fluctuations are still constrained by the measured krypton concentrations throughout SR1.Given its short half-life, low measured concentration, and strong constraints from the leak hypothesis, we conclude that 37 Ar cannot make up the excess, although it may be present in the detector at a negligible level.
We also considered an additional background that has never been observed before in LXe TPCs: the β emission of tritium, 2 which has a Q value of 18.6 keV and a half-life of 12.3 years [100].Tritium may be introduced from predominantly two sources: cosmogenic activation of xenon during above-ground exposure [101] and emanation of tritiated water (HTO) and hydrogen (HT) from detector materials due to its cosmogenic and anthropogenic abundance.In contrast to 127 Xe and 37 Ar, the tritium hypothesis cannot be ruled out.In Sec.IVA, we consider several possible mechanisms for the introduction of tritium into the detector and the uncertainties involved in its production and reduction processes in an attempt to estimate its concentration.

A. Tritium hypothesis
In order to determine the hypothetical concentration of tritium required to account for the excess, we search for a 3 H "signal" on top of the background model B 0 .When compared to B 0 , the tritium hypothesis is favored at 3.2σ and the fitted rate is 159 AE 51 events=ðt • yÞ (68% C.L.), which would correspond to a 3 H=Xe concentration of ð6.2 AE 2.0Þ × 10 −25 mol=mol.As tritium is expected to be removed by the xenon purification system, this concentration would correspond to an equilibrium value between emanation and removal.The spectral fits under this hypothesis are illustrated in Fig. 7(a).
Due to its minute possible concentration, long half-life with respect to our exposure, and the fact that it decays through a single channel, we are unable to confirm the presence of tritium from SR1 data directly.We therefore try to infer its concentration from both initial conditions and detector performance parameters.
A tritium background component from cosmogenic activation of target materials has been observed in several dark matter experiments at rates compatible with predictions [102], although it has never before been detected in xenon.From exposure to cosmic rays during above-ground storage of xenon, we estimate a conservative upper limit on the initial 3 H=Xe concentration of <4 × 10 −20 mol=mol, based on GEANT4 activation rates [101] and assuming saturation activity.At this stage, tritium will predominantly take the form of HTO, given the measured ppm water impurities in the xenon gas and equilibrium conditions [90,103].Through xenon gas handling prior to filling the detector (i.e., condensation of H 2 O=HTO on the walls of the cooled xenon-storage vessel) and purification via a high-efficiency getter with a hydrogen removal unit [2,104], we expect the concentration to be reduced to <10 −27 mol=mol, thus reaching negligible levels with respect to the observed excess.
Tritium may also be introduced as HTO and HT via their respective atmospheric abundances.Water and hydrogen, and therefore tritium, may be stored inside materials, such as the TPC reflectors and the stainless steel of the cryostat.This type of source is expected to emanate from detector and subsystem materials at a rate in equilibrium with its removal via getter purification.Tritium can be found in water at a concentration of ð5 − 10Þ × 10 −18 atoms of 3 H for each atom of hydrogen in H 2 O [105][106][107].Here we assume the same abundance of 3 H in atmospheric H 2 as for water. 3Using the best-fit rate of tritium and the HTO atmospheric abundance, a combined (H 2 O þ H 2 ) impurity concentration of ≳30 ppb in the LXe target would be required to make up the excess.Since water impurities affect optical transparency, the high light yield in SR1 indicates an Oð1Þ-ppb H 2 O concentration [68,110], thus implying a maximum contribution from HTO to the 3 H=Xe concentration of ∼1 × 10 −26 mol=mol.With respect to H 2 , we currently have no direct or indirect measurements of its concentration in the detector.Instead, we consider that O 2 -equivalent, electronegative impurities must reach subppb levels in SR1, given the achieved electron lifetime of ∼650 μs (at 81 V=cm) [3,111].Thus, for tritium to make up the excess requires a factor ∼100 higher H 2 concentration than that of electronegative impurities.Under the above assumptions, tritium from atmospheric abundance appears to be an unlikely explanation for the excess.However, we do not currently have measurements of the equilibrium H 2 emanation rate in XENON1T, and thus the HT concentration cannot be sufficiently quantified.
In conclusion, possible tritium contributions from cosmogenic activation or from HTO in SR1 appear too small to account for the excess, while it is not possible to infer the concentration of HT.In addition, various factors contribute further to the uncertainty in estimating a tritium concentration within a LXe environment, such as its unknown solubility and diffusion properties, as well as the possibility that it may form molecules other than HT and HTO.Since 2 Tritium in the form of tritiated methane has been used for calibration of LXe TPCs [97][98][99], including XENON100, but was not used as a calibration source in XENON1T.Following the XENON100 tritium calibration, neither the xenon gas nor the materials that came into contact with the tritiated methane were used in XENON1T. the information and measurements necessary to quantify the tritium concentration are not available, we can neither confirm nor exclude it as a background component.Therefore, we report results using the background model B 0 and then summarize how our results would change if tritium were included as an unconstrained background component.All reported constraints are placed with the validated background model B 0 (i.e., without tritium).

B. Solar axion results
We search for ABC, 57 Fe, and Primakoff axions simultaneously.Under this signal model, B 0 is rejected at 3.4σ, a value determined using toy Monte Carlo methods to account for the three parameters of interest in the alternative hypothesis.A comparison of the best fits under the alternative hypothesis (B 0 þ axion) and null hypothesis (B 0 ) can be found in Fig. 7

(b).
A three-dimensional confidence volume (90% C.L.) was calculated in the space of g ae vs g ae g aγ vs g ae g eff an .This volume is inscribed in the cuboid given by g ae < 3.8 × 10 −12 g ae g eff an < 4.8 × 10 −18 g ae g aγ < 7.7 × 10 −22 GeV −1 : While easy to visualize, this cuboid is more conservative (it displays overcoverage) than the three-dimensional confidence volume it encloses and does not describe the correlations between the parameters.The correlation information can be found in Fig. 8, which shows the two-dimensional projections of the surface.For the ABC-Primakoff and ABC- 57 Fe projections (Fig. 8 top and middle, respectively), g ae can be easily factored out of the y-axis to plot g aγ vs g ae (top) and g eff an vs g ae (middle).This is not as straightforward for the 57 Fe-Primakoff projection (Fig. 8 bottom).Also shown in Fig. 8 are constraints from other axion searches [84,85,[112][113][114][115][116] as well as predicted values from the benchmark QCD models DFSZ and KSVZ.region is anticorrelated in this space and-due to the presence of the low-energy excess-suggests either a nonzero ABC component or nonzero Primakoff component.Since our result gives no absolute lower bound on g ae , the limit on the product g ae g aγ cannot be converted into a limit on g aγ on its own; i.e., with g ae g aγ ¼ 7.6 × 10 −22 GeV −1 , g aγ → ∞ as g ae → 0, as shown in Fig. 8 (top).
Figure 8 (middle) is taken from the projection onto the ABC-57 Fe plane.Unlike the ABC-Primakoff case, these two signals are not degenerate; however, they still display anticorrelated behavior.The reason for this is that the test statistic q [Eq.(17)] is relatively large with small g ae , meaning small changes in the 57 Fe rate about the best fit make q cross the 90% threshold value and thus be excluded by our 90% confidence volume.There is no statistical significance (<1σ) for the presence of a 14.4 keV peak from 57 Fe axions.
Last, Fig. 8 (bottom) shows the projection onto the Primakoff-57 Fe plane, where no correlation is observed.The Primakoff and 57 Fe components are both allowed to be absent as long as there is a nonzero ABC component.This means that, of the three axion signals considered, the ABC component is the most consistent with the observed excess.
The three projections of Fig. 8 can be used to reconstruct the three-dimensional 90% confidence volume for g ae , g ae g aγ , and g ae g eff an .Due to the presence of an excess at low energy, this volume would suggest either a nonzero ABC component or a nonzero Primakoff component.However, the coupling values needed to explain this excess are in strong tension with stellar cooling constraints [114][115][116][117][118], with the exception of a minute region in the 3D coupling space which corresponds to small g ae and large g eff an , g aγ .The CAST constraints [112] as shown are valid for axion masses below 10 meV=c 2 while those from XENON1T and similar experiments hold for all axion masses up to ∼100 eV=c 2 .For an axion mass below 10 meV=c 2 , the CAST result prefers the region with large g ae and small g aγ ; however, there is no tension between the CAST result and this result for higher axion masses (m a > 250 meV=c 2 ) due to the limited sensitivity of CAST for high-mass axions.
As described above, we cannot exclude tritium as an explanation for this excess.Thus, we report on an additional statistical test, where an unconstrained tritium component was added to the background model B 0 and profiled over alongside the other nuisance parameters.In this case, the null hypothesis is the background model plus tritium (B 0 þ 3 H) and the alternative includes the three axion signal components (B 0 þ 3 H þ axion), where tritium is unconstrained in both cases.The solar axion signal is still preferred in this test, but its significance is reduced to 2.0σ.The fits for this analysis are shown in Fig. 7(d).The tritium component is negligible in the alternate best fit, but its presence allows for a better fit under-and thus a reduced significance of rejecting-the null hypothesis.
FIG. 8. Constraints on the axion-electron g ae , axion-photon g aγ , and effective axion-nucleon g eff an couplings from a search for solar axions.The shaded blue regions show the twodimensional projections of the three-dimensional confidence surface (90% C.L.) of this work and hold for m a < 100 eV=c 2 .See text for more details on the three individual projections.All three plots include constraints (90% C.L.) from other axion searches, with arrows denoting allowed regions, and the predicted values from the benchmark QCD axion models DFSZ and KSVZ.

C. Neutrino magnetic moment results
When compared to the neutrino magnetic moment signal model, the background model B 0 is rejected at 3.2σ.The best fits of the null (B 0 ) and alternative (B 0 þ μ ν ) hypotheses for this search are shown in Fig. 7(c).
The 90% confidence interval for μ ν from this analysis is given by μ ν ∈ ð1.4; 2.9Þ × 10 −11 μ B and is shown in Fig. 9 along with the constraints from other searches.The upper boundary of this interval is very close to the limit reported by Borexino [42], which is currently the most stringent direct detection constraint on the neutrino magnetic moment.Similar to the solar axion analysis, if we infer the excess as a neutrino magnetic moment signal, our result is in strong tension with indirect constraints from analyses of white dwarfs [119] and globular clusters [44].The result is also compatible with the constraint from XENON1T using the S2-only method, which is able to probe a lower energy region and is further discussed in Sec.IV E. It is important to note that the neutrino flavor does impact the interaction involving the magnetic moment, which in reality is a 3 × 3 matrix due to neutrino mixing.Our result, based on a flavor-insensitive detection of solar neutrinos, is thus directly comparable to Borexino's, but not necessarily to Gemma's (reactor electron antineutrinos) or the astrophysical limits (electron neutrinos).
As in Sec.IV B, we report on an additional statistical test where an unconstrained tritium component was included in both null and alternative hypotheses.In this test, the significance of the neutrino magnetic moment signal is reduced to 0.9σ with the presence of a tritium background.This is the most sensitive search to date for an enhanced neutrino magnetic moment with a dark matter detector and suggests that this beyond-the-SM signal be included in the physics reach of other dark matter experiments.

D. Bosonic dark matter results
For bosonic dark matter, we iterate over (fixed) masses between 1 and 210 keV=c 2 to search for peaklike excesses.The trial factors to convert between local and global significance were extracted using toy Monte Carlo methods.While the excess does lead to looser constraints than expected at low energies, we find no global significance over 3σ for this search under the background model B 0 .We thus set an upper limit on the couplings g ae and κ as a function of particle mass.
These upper limits (90% C.L.) are shown in Fig. 10, along with the sensitivity band in green (1σ) and yellow (2σ).The losses of sensitivity at 41.5 and 164 keV are due to the 83m Kr and 131m Xe backgrounds, respectively, and the gains in sensitivity at around 5 and 35 keV are due to increases in the photoelectric cross section in xenon.The fluctuations in our limit are due to the photoelectric cross section, the logarithmic scaling, and the fact that the energy spectra differ significantly across the range of masses.For most masses considered, XENON1T sets the most stringent direct-detection limits to date on pseudoscalar and vector bosonic dark matter couplings.
Due to the presence of the excess, we performed an additional fit using the bosonic dark matter signal model, with the particle mass allowed to vary freely between 1.7 and 3.3 keV=c 2 .The result gives a favored mass value of ð2.3 AE 0.2Þ keV=c 2 (68% C.L.) with a 3.0σ global (4.0σ local) significance over background.A log-likelihood ratio curve as a function of mass is shown in Fig. 11 (left), along with the asymptotic 1-σ uncertainty and the local significance for each mass.The spectral fit of the 2.3 keV peak is illustrated in Fig. 11 (right).Since the energy reconstruction in this region is validated using 37 Ar calibration data, whose distribution has a mean value within <1% of the expectation at 2.82 keV [86], this analysis can also be used to compare the data to potential monoenergetic backgrounds in this region.

E. Additional checks
Here we describe a number of additional checks to investigate the low-energy excess in the context of the tritium, solar axion, and neutrino magnetic moment hypotheses.
The time dependence of events with energies in (1, 7) keV in SR1 was investigated and found to be FIG.9. Constraints (90% C.L.) on the neutrino magnetic moment from this work compared to experiments Borexino [42] and Gemma [120], along with astrophysical limits from the cooling of globular clusters [44] and white dwarfs [119].The constraint from XENON1T using ionization signal only (S2 only) is also shown (see Sec. IV E).Arrows denote allowed regions.The upper boundary of the interval from this work is about the same as that from Borexino and Gemma.If we interpret the low-energy excess as a neutrino magnetic moment signal, its 90% confidence interval is in strong tension with the astrophysical constraints.
inconclusive.The event rate is slightly higher in the beginning of SR1, but the rate evolution is statistically consistent with (i) a constant rate, (ii) a constant background rate (B 0 ) plus a subtle ∼7% (peak-to-peak) rate modulation from the change in Earth-Sun distance, and (iii) a constant background rate (B 0 ) plus an exponentially decreasing component with a fixed half-life of 35 days ( 37 Ar half-life) or 12.3 years ( 3 H half-life).As another test of time dependence, we split SR1 into three periods with equal exposure and fit the data in each period with the ABC solar axion signal model.Similarly to the (1,7) keV rate evolution, the best-fit signal rate is the highest in the first period of SR1, but is not statistically significant as the signal rate is consistent within uncertainty between the three periods.We therefore conclude that, due to limited statistics, at this time we cannot use time dependence to exclude any of the hypotheses discussed in this work.More detailed time dependence studies will be presented in a forthcoming publication.
Since the excess events have energies near our 1 keV threshold, where the efficiency is ∼10%, we considered higher analysis thresholds to check the impact of this choice on the results.With the excess most prominent between 2 and 3 keV, where the respective detection efficiencies are ∼80% and 94%, changing the analysis threshold has little impact unless set high enough so as to remove the events in question.This is not well motivated, given the high efficiency in the region of the excess.For all thresholds considered (namely, 1.0, 1.6, 2.0, 3.0 keV), the solar axion model gives the best fit to the data.We hence conclude that our choice of analysis threshold impacts neither the presence nor interpretation of the low-energy excess.
We also checked data from Science Run 2 (SR2), an R&D science run that followed SR1, in an attempt to understand FIG. 10.Constraints on couplings for bosonic pseudoscalar ALP (top) and vector (bottom) dark matter, as a function of particle mass.The XENON1T limits (90% C.L.) are shown in black with the expected 1 (2) σ sensitivities in green (yellow).Limits from other detectors or astrophysical constraints are also shown for both the pseudoscalar and vector cases [53,84,85,[121][122][123][124][125][126][126][127][128].FIG.11.Left: the log-likelihood ratio q for different bosonic dark matter masses with respect to the best-fit mass at 2.3 keV=c 2 .At each mass, we show the result for the corresponding best-fit coupling.The green band shows an asymptotic 68% C.L. confidence interval on the bosonic dark matter mass.The local significance for each mass is indicated by the right y-axis.Right: best fit of a 2.3 keV peak and B 0 to the data.A 0.4 keV binning is used for better visualization.the observed excess.Many purification upgrades were implemented during SR2, including the replacement of the xenon circulation pumps with units that (i) are more powerful, leading to improved purification speed, and (ii) have lower 222 Rn emanation, leading to a reduced 214 Pb background rate in the TPC [129,130], which is further decreased by online radon distillation.The resulting increased purification speed and reduced background make SR2 useful to study the tritium hypothesis.If the excess were from tritium (or another non-noble contaminant), we would expect its rate to decrease due to the improved purification; on the other hand, the rate of the signal hypotheses would not change with purification speed.
While the SR2 purification upgrades allowed for an improved xenon purity and a reduced background level, the unavoidable interruption of recirculation for the upgrades also led to less stable detector conditions.Thus, in addition to a similar event selection process as SR1 in Sec.III A, we removed several periods of SR2 for this analysis to ensure data quality.Periods where the electron lifetime changed rapidly due to tests of the purification system were removed to reduce uncertainty in the energy reconstruction.We also removed data sets during which a 83m Kr source was left open for calibration.Data within 50 days of the end of neutron calibrations were also removed to reduce neutron-activated backgrounds and better constrain the background at low energies.After the other selections, this data would have only added ∼10 days of live time; thus, for simplicity, it was removed rather than fit separately like the SR1 data set.With these selections, the effective SR2 live time for this analysis is 24.4 days, with an average ER background reduction of 20% in (1,30) keV as compared to SR1.
A profile likelihood analysis was then performed on SR2 with a similar background model as SR1, denoted as B SR2 .Since we are primarily interested in using this data set to test the tritium hypothesis, we focus on the tritium results.
Similarly to SR1, we search for a tritium signal on top of the background model B SR2 and find that the background-only hypothesis is slightly disfavored at 2.3σ.The SR2 spectrum, along with the fits for the null (B SR2 ) and alternative (B SR2 þ 3 H) hypotheses, can be found in Fig. 12.A log-likelihood ratio curve for the tritium component is given in Fig. 13, which shows that the fitted tritium rate is 320 AE 160 events=ðt • yÞ, higher than that from SR1 but consistent within uncertainties.The rate uncertainty in SR2 is much larger than that in SR1 due to limited statistics.The solar axion and magnetic moment hypotheses give similar results, with significances ∼2σ and best-fit values larger than, but consistent with, the respective SR1 fit results.Thus, these SR2 studies are largely inconclusive.
Last, we also checked these hypotheses in a different energy region using the so-called S2-only approach, where the requirement for an S1 signal is dropped, allowing for a ∼200 eV energy threshold.XENON1T's S2-only analysis [123] was used to place limits on the tritium rate ½<2256 events=ðt • yÞ and g ae (<4.8 × 10 −12 ) that are far greater than, and therefore consistent with, the constraints derived here.The S2-only analysis is not as sensitive to the tritium and axion signals because both spectra peak above 1 keV.On the other hand, many of the predicted signal events from neutrino magnetic moment fall below 1 keV as the rate increases with falling energy, so the S2-only search is more relevant for this hypothesis.It yields a 90% C.L. one-sided limit of μ ν < 3.1 × 10 −11 μ B , consistent with the upper boundary of the 90% confidence interval obtained in Sec.IV C. Therefore, none of the discussed hypotheses are in conflict with the S2-only result.FIG.12.A fit to SR2 data if tritium is treated as a signal.The red (gray) line is the fit with (without) tritium in the background model.FIG. 13.The log-likelihood ratio curve for the tritium rate in SR2.The orange line and band indicate the best fit and 1σ uncertainty for the tritium rate in SR1.The SR2 fit result is consistent with SR1, but with a large uncertainty due to limited statistics.

V. DISCUSSION
We observe an excess of electronic recoil events at low energies in XENON1T data.In the reference region of 1-7 keV, 285 events are observed whereas 232 AE 15 events are expected from the background-only fit to the data.The β decay of tritium is considered as a possible explanation, as it has a similar spectrum to that observed and is expected to be present in the detector at some level.We are unable to independently confirm the presence of tritium at the Oð10 −25 Þ mol=mol concentration required to account for the excess, and so treat it separately from our validated background model.If electronic recoils from tritium decay were the source of the excess, this would be its first indication as an atmospheric source of background in LXe TPCs.The tritium hypothesis clearly represents a possible SM explanation for the excess, but-based on spectral shape alone-the solar axion model is the most favored signal by the data at 3.4σ, albeit at only ∼2σ if one considers tritium as an additional background.
If this excess were a hint of a solar axion, our result would suggest either (i) a nonzero rate of ABC axions or (ii) a nonzero rate of both Primakoff and 57 Fe axions.If we interpret the excess as an ABC axion signal (i.e., take g aγ and g eff an to be zero), the required value of g ae is smaller than that ruled out by other direct searches but has a clear discrepancy with constraints from indirect searches [44,131].These constraints are a factor of ∼5-10 lower than reported here, although subject to systematic uncertainties.It is noteworthy that some of these astrophysical analyses, while their constraints are still stronger than direct searches, do in fact suggest an additional source of cooling compatible with axions [117,131].If the indirect hints and the XENON1T excess were indeed explained by axions, the tension in g ae could be relieved by underestimated systematic uncertainties in, e.g., stellar evolution theory [44] or white dwarf luminosity functions [132], or by a larger solar axion flux than that given in [14].
Although not considered in this work, XENON1T is also directly sensitive to the axion-photon coupling via the inverse Primakoff effect, whereby a solar axion coherently scatters off the effective electric field of the xenon atom, thus producing an outgoing photon and inducing an electronic recoil.This detection channel was first considered only recently for xenon-based detectors in [133,134], which demonstrated that the tension of axion-photon coupling between the XENON1T excess and stellar constraints can be significantly reduced.
Continuing to interpret the excess as a hypothetical QCD axion signal, we can extend the analysis to make statements on the axion mass m a under assumptions of different models, as outlined in Sec.II A. As examples, we consider a DFSZ model with variable β DFSZ and KSVZ model with variable electromagnetic anomaly E (for simplicity, we fix the color anomaly N ¼ 3).Comparing these two classes of models with our 90% confidence surface, we find that both are consistent with our result for a subset of parameters.For the DFSZ model, we find m a ∼ 0.1-4.1 eV=c 2 and cos 2 β DFSZ ∼ 0.01-1 would be consistent with this work.Alternatively, under the KSVZ model, m a ∼ 46-56 eV=c 2 and E ¼ 6 would be similarly consistent.These modelspecific mass ranges are not confidence intervals, as their specific assumptions were not included when constructing Fig. 8.We instead report a single, model-independent confidence region on the couplings to allow comparison with a variety of models, not just the examples mentioned here.
Additionally, we describe a direct search for an enhanced neutrino magnetic moment.This signal also has a similar spectrum to the excess observed, but at 3.2σ displays a lower significance than that from solar axions.We report a confidence interval of μ ν ∈ ð1.4;2.9Þ×10−11 μ B (90% C.L.), the upper boundary of which is very close to the worldleading direct limit reported by Borexino [42].This shows that dark matter experiments are also sensitive to beyond-SM physics in the neutrino sector.Here we only search for an enhanced neutrino-electron cross section due to an anomalous magnetic moment, but a similar enhancement would also occur in neutrino-nucleus scattering [135].With the discrimination capabilities of LXe TPCs to ER and NR events, it would be interesting to consider this channel in future searches as well.
If from an astrophysical source, the excess presented here is different from the result reported by the DAMA experiment, which claims that an observed annual modulation of events between 1 and 6 keV might be due to a dark matter signal [136,137].We present here a leptophilic dark matter model, where WIMPs couple with electrons through an axial-vector interaction [138].This model was used to explain the DAMA signal but was rejected already by the XENON100 experiment [139].Interpreting the modulating source observed by DAMA under this model, the expected signal rate in the XENON1T detector would be more than 2 orders of magnitude higher than the total event rate we observed, as shown in Fig. 14.Consequently, the excess observed in this work is unrelated to the one observed by DAMA.

VI. SUMMARY
We report on searches for new physics using low-energy electronic recoils in XENON1T.In a search for bosonic dark matter, world-leading constraints are placed on the interaction strengths of pseudoscalar and vector particles.An excess is observed at low energies that is consistent with a solar axion signal, a bosonic dark matter signal with a mass of 2.3 keV=c 2 , a solar neutrino signal with enhanced magnetic moment, or a possible tritium background.We are unable to confirm nor exclude the presence of tritium at this time.
In an attempt to understand the low-energy excess, we performed a number of additional studies.The analysis of an additional data set called SR2-which displays a ∼20% lower background rate but only ∼10% statistics compared to SR1-is consistent with the SR1 analysis but largely inconclusive about the nature of the excess.An S2-only search, which is able to probe sub-keV energies, similarly yielded consistent constraints for all the discussed hypotheses.Compared to the excess observed by DAMA, it is much lower in rate and thereby unrelated.
The signals discussed here can be further explored in the next-generation detectors, such as the upcoming PandaX-4T [140], LUX-ZEPLIN [141], and XENONnT [142] experiments.The next phase of the XENON program, XENONnT, featuring a target mass of 5.9 tonnes and a factor of ∼6 reduction in ER background, will enable us to study the excess in much more detail if it persists.Preliminary studies based on the best-fit results of this work suggest that a solar axion signal could be differentiated from a tritium background at the 5σ level after only a few months of data from XENONnT.

APPENDIX: β SPECTRA MODELING
This Appendix briefly describes the different theoretical models used in the present work to compute the β spectra for 214 Pb, 212 Pb, and 85 Kr.

GEANT4 radioactive decay module
The radioactive decay module (RDM) in GEANT4 4simulates the decay of a given radionuclide using the nuclear data taken from an evaluated nuclear structure data file (ENSDF) [143].The required β spectra are generated in a dedicated class using an analytical model.The β spectral shape, i.e., the unnormalized emission probability per electron energy, is derived from Fermi's golden rule as dN dW ∝ pWq 2 FðZ; WÞCðWÞSðZ; WÞ; ðA1Þ with Z the atomic number of the daughter nucleus.Here, W is the total energy of the β particle and is related to its kinetic energy E by W ¼ 1 þ E=m e , with m e the electron rest mass.The maximum energy W 0 is defined identically from the energy of the transition E 0 .The β particle momentum is p ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi W 2 − 1 p and the (anti)neutrino momentum is q ¼ W 0 − W, assuming a massless particle (m ν ¼ 0).
The Fermi function FðZ; WÞ corrects for the static Coulomb effect of the nucleus on the β particle.Considering the Coulomb field generated by a pointlike nucleus, the Dirac equation can be solved analytically and the well-known expression of the Fermi function can be derived.GEANT4 follows the approximate expression of the Fermi function from [144].
The shape factor CðWÞ takes into account the nuclear and lepton matrix elements.Assuming constant values of the lepton wave functions within the nuclear volume, one can demonstrate that allowed and forbidden unique transitions can be calculated without involving the structure of the nucleus.For an allowed transition, the shape factor is constant: CðWÞ ¼ 1.In GEANT4, first, second, and third forbidden unique transitions are calculated following the approximate expressions given in [145] that were established by considering the analytical solutions of the Dirac equation, the same as for the Fermi function.In any other case, the decay is treated as allowed.
The atomic screening effect corresponds to the influence of the electron cloud surrounding the daughter nucleus on the β particle wave function.GEANT4 takes this into account following the most widespread approach set out by Rose in [146] almost a century ago.For a β electron, this effect is evaluated by subtracting from the particle energy W a constant Thomas-Fermi potential V 0 which only depends on Z.This corrected energy W 0 ¼ W − V 0 replaces W in all the quantities required for the calculation of the spectral shape, except in the (anti)neutrino energy q because this neutral particle is not affected by the Coulomb field.The parametrization of the potential used in GEANT4 is close to the prescription given in [147].The screening correction is then given by SðW; ZÞ ¼ It is noteworthy that this correction can only be applied for W ≥ V 0 , which creates a nonphysical discontinuity in the spectrum at W ¼ V 0 , as seen in Fig. 15.

IAEA LiveChart
The β spectra available on the IAEA LiveChart website [61] are produced with the first version of the BetaShape program [148].The required information for each transition is taken from the most recent ENSDF file with results from the latest nuclear data evaluation [59].
The physics model in BetaShape has already been detailed in [62], except for the atomic screening effect.The β spectral shape is described in the Behrens and Bühring formalism [149] with all quantities as defined before.The quantity RðZ; WÞ are the radiative corrections described below.
In this formalism, the Fermi function is defined from the Coulomb amplitudes α k of the relativistic electron wave functions, These wave functions are numerical solutions of the Dirac equation for the Coulomb potential of a nucleus modeled as a uniformly charged sphere.Indeed, no analytical solution exists even for such a simple potential; however, the method from [149] allows for a precise, and fast, calculation of the Coulomb amplitudes.The method inherently accounts for the finite nucleus size, while other methods usually require an analytical correction [L 0 in Eq. (A4)].The total angular momentum change ΔJ ¼ jJ i − J f j and the parity change π i π f between the initial and final nuclear states are from the input ENSDF file and determine the nature of the transition.Given that L ¼ 1 if ΔJ ¼ 0 or 1 for an allowed transition, and L ¼ ΔJ for any (L − 1)th forbidden unique transition, the theoretical shape factor can be expressed as The λ k parameters are defined from the Coulomb amplitudes α k by In the case of forbidden nonunique transitions, the structures of the initial and final nuclear states must be taken into account, which greatly complicates the calculation.The usual approximation consists of treating such a transition as a forbidden unique transition of identical ΔJ.The validity of this approximation, minutely tested in [62], FIG. 15.Low-energy part of the β spectral shape of the ground state to ground-state transition in 214 Pb decay.This first forbidden nonunique transition was calculated as allowed in every case but with different levels of approximations, as described in the text.The four spectra are normalized by area over the full energy range.See text for details on the shape of each spectrum.can be demonstrated only for some first forbidden nonunique transitions, which are then calculated as allowed.Its generalization to every forbidden nonunique transition is implemented in BetaShape.
The radiative corrections are nonstatic Coulomb corrections from quantum electrodynamics.They can be split into two parts: the inner corrections, which are independent of the nucleus, and the outer corrections, which depend on the nucleus.Only the latter depend on the β particle energy.The outer radiative corrections RðZ; WÞ take into account the internal bremsstrahlung process, by which the β particles lose energy in the electromagnetic field of the nucleus.For allowed transitions, analytical corrections were derived in [150,151] and are implemented in the first version of BetaShape as described in [62].
Finally, the spectral shape is modified by applying the screening correction SðW; ZÞ.The BetaShape program includes an analytical correction based on the work of Bühring [152] that is more precise than Rose's correction.The most realistic, spatially varying screened potentials at the time were of Hulthén type (see [153] and references therein).Bühring first developed a version of the Dirac equation that correctly includes Hulthén's potentials but simplified the angular momentum dependency, allowing analytical solutions to be established [154].He then performed in [152] a radial expansion at the origin of both the wave functions and the Coulomb potential, including Hulthén's screened potential, and retained only the dominant term.This procedure allows the determination of screened-to-unscreened ratios for the Fermi function F 0 L 0 and the λ k parameters, which are then used to correct for screening in Eq. (A3).Therefore, the quantity SðZ; WÞ in Eq. (A3) is more a symbolic notation.In BetaShape, this approach is used with Salvat's screened potentials [155], which can be expanded at the origin as where β is determined from the parameters A i and α i given in [155], These potentials are widely used for their precision and completeness.It is noteworthy that Bühring's correction does not create any nonphysical discontinuity in the spectrum as in Rose's correction.However, it tends to greatly decrease the emission probability at low energy.

Improved calculations
When high precision at low energy is required, the modeling of β − decays must include the atomic screening and exchange effects.The two approximate screening corrections previously described are not sufficient.The exchange effect is even more significant and comes from the indistinguishability of the electrons.The regular, direct decay corresponds to the creation of the β electron in a continuum orbital of the daughter atom.In the exchange process, the β electron is created in an atomic orbital of the daughter atom and the atomic electron which was present in the same orbital in the parent atom is ejected to the continuum.This process leads to the same final state as the direct decay, i.e., one electron in the continuum, and is possible because the nuclear charge changes in the decay.
Precise relativistic electron wave functions are necessary to calculate such effects.The numerical procedure was described in detail in [156], with the nucleus modeled as a uniformly charged sphere.For the continuum states, the Coulomb potential includes the appropriate Salvat screened potential.The wave functions, and therefore the Fermi function F 0 L 0 and the λ k parameters, inherently take into account the screening effect.For the bound states, an exchange potential has to be added to this Coulomb potential and a specific procedure was implemented to ensure good convergence to precise atomic energies.In [156], the one-electron energies from [157] were considered while in the present work, the more accurate orbital energies from [158] that include electron correlations are used.
A precise description of the exchange effect was set out in detail in [60,159], but only for the allowed transitions.In such a case, β electrons are created in continuum states with quantum number κ ¼ AE1 and the selection rules imply that exchange can only occur with atomic electrons of identical κ, i.e., in s 1=2 (κ ¼ −1) and p 1=2 (κ ¼ þ1) orbitals.The influence of the exchange effect can then be taken into account through a correction factor on Eq. (A3), The total exchange correction is defined by All primed quantities refer to the daughter atom and to the parent atom otherwise.The large and small components of the relativistic electron wave functions, respectively, g c κ and f c κ for the continuum states and g b n;κ and f b n;κ for the bound states, respectively, are calculated at the nuclear radius R. The quantities T −1 and T þ1 depend on the overlaps between the bound states of the parent atom and the continuum states of the daughter atom with energy E, The sums are running over all occupied orbitals of the daughter atom of same quantum number κ.
It is noteworthy that in [156], only the s 1=2 orbitals were taken into account, following the prescription in [60].The "new screening correction" proposed in [156] was necessary to reproduce the experimental β spectra of 63 Ni and 241 Pu, but was later found to be incompatible with a rigorous derivation of the β spectrum starting from the decay Hamiltonian and the corresponding S matrix.If correct screening and exchange effect with s 1=2 and p 1=2 orbitals are considered, together with precise atomic orbital energies, excellent agreement over the entire energy range of the two spectra is obtained.
Finally, more precise radiative corrections have been considered compared with those previously described.They were developed using more recent mathematical techniques and a significant change in the correction terms was found [160].Describing the various changes is out of the scope of the present work; however, many details can be found in [161].The influence of these new radiative corrections on the integrated β spectrum is given for 20 superallowed transitions in [162], for which an excellent agreement is obtained with the present implementation.It appears that these corrections are significantly smaller than the previous ones, especially for high atomic numbers.

Application to the transitions of interest
These different models have been applied to the ground state to ground-state transitions in 212 Pb, 214 Pb, and 85 Kr decays.The resulting spectra are similar to each other in the major part of the energy range, except at low energy.
The differences are illustrated in Fig. 15 for the lowenergy region of the 214 Pb β spectrum.The yellow curve is the GEANT4 RDM model as described in the Appendix and the nonphysical discontinuity due to the screening correction is clearly visible at 12 keV.The red curve is from IAEA LiveChart, thus generated by the first version of the BetaShape program as described in the Appendix.One can see the effect of Bühring screening correction that tends to decrease the emission probability.The cyan and blue curves were determined as described in the Appendix, without and with the atomic exchange correction, respectively.The screening effect is found to have a much smaller influence on the spectral shape when determined using a full numerical procedure than when applying an analytical approximation.However, the atomic exchange effect has a strong influence, as expected from previous studies [156].
Both transitions in 212 Pb and 214 Pb ground state to ground-state decays were calculated as allowed, accordingly with the approximation described in the Appendix.It is important to keep in mind that formally, such first forbidden nonunique transitions should be determined including the structure of the initial and final nuclear states, a much more complicated calculation that is beyond the present scope.
The transition in 85 Kr decay is first forbidden unique and can thus be calculated accurately without nuclear structure.The description of the exchange effect from [60,159] used here is only valid for allowed transitions.For a first forbidden unique transition, one can expect a contribution of the κ ¼ AE2 atomic orbitals but the exact solutions have still to be derived.However, the spectral shape is derived from a multipole expansion of the nuclear and lepton currents, as shown in the shape factor in Eq. (A5).Therefore, one can expect that the allowed exchange correction should give the main contribution, and this was done to determine the 85 Kr β spectrum.
The tritium β spectrum used in this work was obtained from the IAEA LiveChart [61], thus calculated using the standard Fermi function without corrections.As 3 H decays via an allowed transition, this spectrum is sufficiently precise at energies above 0.5 keV, as confirmed experimentally in [82].

Uncertainties
The dominant contribution to the continuous XENON1T low-energy background comes from 214 Pb β decay.We thus focus the uncertainty discussion on the 214 Pb ground state to ground-state transition, calculated for the final model (blue curve in Fig. 15).
The transition energy is directly given by the Q value [163]: Q β ¼ 1018ð11Þ keV.This uncertainty can be propagated by calculating the spectrum at ð1018 AE 11Þ keV, namely, at 1σ.The result is an envelope centered on the spectrum calculated at the Q value, which provides an uncertainty on the emission probability for each energy bin.The relative uncertainty is 1.7% below 10 keV and 1.1% at 210 keV.However, most of this uncertainty is removed because the 214 Pb spectrum is left unconstrained in the fitting procedure.The remaining uncertainty component on the emission probability is ∼0.5% for each energy bin, in which the shape of the spectrum cannot vary steeply.
The atomic screening effect only slightly modifies the shape of the β spectrum.Its uncertainty contribution can thus be safely ignored.The atomic exchange effect strongly affects the spectral shape below 5 keV, and its accuracy depends on the atomic model used.For the β spectra of 63 Ni and 241 Pu, the residuals between their high-precision measurement and the improved calculation in the Appendix showed that the agreement is better than the statistical fluctuations due to the number of counts in each energy channel, from 0.5 keV to the endpoint energy.A conservative value of 1% for each energy bin is the maximum relative uncertainty and is the value adopted here.
The 214 Pb transition of interest is first forbidden nonunique.As explained in the Appendix, the nuclear structure should be taken into account for such a transition because it has an influence on the spectral shape.Treating it as an allowed transition induces an inaccuracy which cannot be estimated by comparison with a measured spectrum-no measurement has been reported so far.In the same mass region, the 210 Bi decay exhibits also a first forbidden nonunique, ground state to ground-state transition with a comparable Q value, and an experimental shape factor is available.As can be seen in [164], treating this transition as allowed leads to an important discrepancy with measurement.The question is then how this observation can be used for assessing an uncertainty to the 214 Pb spectral shape.
First, allowing the rate normalization to be free in the (1, 210) keV region absorbs the vast majority of any difference.Second, the nuclear structures of 210 Bi and 214 Pb are not identical.The 210 Bi decay can be seen as two nucleons in the valence space above the doubly magic 208 Pb core, with the initial configuration ðp; 1h 9=2 Þðn; 2g 9=2 Þ and two protons in the 1h 9=2 orbital in the final state.However, this picture is too simple to be accurate because the core is not really inert.Nucleons from the core can give contributions to the β-decay matrix elements, mainly through meson exchange effects and core polarization effects [165].In 214 Pb decay, a single proton in the 1h 9=2 orbital is present in the final state and in the initial state, six neutrons are spread over the orbitals of the valence space but tend to couple to each other through pairing and dominantly occupy the 2g 9=2 orbital.Contributions from the core nucleons can be expected to be relatively small compared to the main ðn; 2g 9=2 Þ → ðp; 1h 9=2 Þ transition.In addition, even though it is difficult to predict if the nuclear structure component shifts the spectrum to lower energies, as for 210 Bi, or to higher energies, a steep variation at low energy is not realistic.
To conclude, we conservatively estimate a relative uncertainty on the spectral shape of 5% due to the nuclear structure component and an additional 1% for the energy dependency of the relative uncertainty on the maximum energy.Thus, a 6% total uncertainty on the spectral shape is estimated for the 214 Pb β-decay model in this work.

FIG. 4 .
FIG.4.A zoomed-in and rebinned version of Fig.3 (top), where the data display an excess over the background model B 0 .In the following sections, this excess is interpreted under solar axion, neutrino magnetic moment, and tritium hypotheses.

FIG. 6 .
FIG.6.Fit to220 Rn calibration data with a theoretical β-decay model (see the Appendix) and the efficiency nuisance parameter, using the same unbinned profile likelihood framework described in Sec.III C.This fit suggests that the efficiency shown in Fig.2describes well the expected spectrum from 214 Pb, the dominant background at low energies.

3
Although geographical and temporal HT abundances in the atmosphere vary due to anthropogenic activities, HT that reaches the Earth's surface undergoes exchange to HTO within 5 hours[108,109].E. APRILE et al.PHYS.REV.D 102, 072004 (2020) 072004-12

FIG. 7 . 13 Figure 8 (
FIG.7.Fits to the data under various hypotheses.The null and alternative hypotheses in each scenario are denoted by gray (solid) and red (solid) lines, respectively.For the tritium (a), solar axion (b), and neutrino magnetic moment (c) searches, the null hypothesis is the background model B 0 and the alternative hypothesis is B 0 plus the respective signal.Contributions from selected components in each alternative hypothesis are illustrated by dashed lines.Panel (d) shows the best fits for an additional statistical test on the solar axion hypothesis, where an unconstrained tritium component is included in both null and alternative hypotheses.This tritium component contributes significantly to the null hypothesis, but its best-fit rate is negligible in the alternative hypothesis, which is illustrated by the orange dashed line in the same panel.

FIG. 14 .
FIG. 14.Comparison between DAMA expected signals and XENON1T data (signal plus background).Dotted lines represent the expected signal spectra of selected masses in the XENON1T detector if the DAMA modulated signals are interpreted as WIMPs scattering on electrons through axial-vector interactions.XENON1T data are indicated by black points and the background model B 0 is illustrated by the red line.The right bound of the shaded region shows the threshold in this analysis.

TABLE I . Summary of components in the background model B 0 with expected and fitted number of events in the 0.65 tonne- year exposure of SR1. Both numbers are within the (1, 210) keV ROI and before efficiency correction. See text for details on the various components.
The efficiency parameter, which is dominated by detection efficiency and does not change with time.(ii) The 214 Pb component, which was determined to have a constant rate in time using detailed studies of the α decays of the 222 Rn and 218 Po, as well as the coincidence signature of 214 Bi and 214 Po.