Proton-philic spin-dependent inelastic Dark Matter (pSIDM) as a viable explanation of DAMA/LIBRA-phase2

We show that the Weakly Interacting Massive Particle scenario of proton-philic spin-dependent inelastic Dark Matter (pSIDM) can still provide a viable explanation of the observed DAMA modulation amplitude in compliance with the constraints from other experiments after the release of the DAMA/LIBRA-phase2 data and including the recent bound from COSINE-100, that uses the same $NaI$ target of DAMA. The pSIDM scenario provided a viable explanation of DAMA/LIBRA--phase1 both for a Maxwellian WIMP velocity distribution and in a halo-independent approach. At variance with DAMA/LIBRA-phase1, for which the modulation amplitudes showed an isolated maximum at low energy, the DAMA/LIBRA-phase2 spectrum is compatible to a monotonically decreasing one. Moreover, due to its lower threshold, it is sensitive to WIMP-iodine interactions at low WIMP masses. Due to the combination of these two effects pSIDM can now explain the yearly modulation observed by DAMA/LIBRA only when the WIMP velocity distribution departs from a standard Maxwellian. In this case the WIMP mass $m_{\chi}$ and mass splitting $\delta$ fall in the approximate ranges 7 GeV $\lesssim m_{\chi}\lesssim$ 17 GeV and 18 keV$\lesssim\delta\lesssim$29 keV. The recent COSINE-100 bound is naturally evaded in the pSDIM scenario due to its large expected modulation fractions.

We show that the Weakly Interacting Massive Particle scenario of proton-philic spin-dependent inelastic Dark Matter (pSIDM) can still provide a viable explanation of the observed DAMA modulation amplitude in compliance with the constraints from other experiments after the release of the DAMA/LIBRA-phase2 data. The pSIDM scenario provided a viable explanation of DAMA/LIBRA-phase1 both for a Maxwellian WIMP velocity distribution and in a haloindependent approach. At variance with DAMA/LIBRA-phase1, for which the modulation amplitudes showed an isolated maximum at low energy, the DAMA/LIBRA-phase2 spectrum is compatible to a monotonically decreasing one. Moreover, due to its lower threshold, it is sensitive to WIMP-iodine interactions at low WIMP masses. Due to the combination of these two effects pSIDM can now explain the yearly modulation observed by DAMA/LIBRA only when the WIMP velocity distribution departs from a standard Maxwellian. In this case the WIMP mass mχ and mass splitting δ fall in the approximate ranges 6 GeV < ∼ mχ < ∼ 17 GeV and 17 keV< ∼ δ < ∼ 29 keV.

I. INTRODUCTION
About one quarter of the total mass density of the Universe [1] and more than 90% of the halo of our Galaxy are believed to be constituted by Dark Matter (DM) and Weakly Interacting Massive Particles (WIMPs) are one of the most popular candidates to compose it. The scattering rate of DM WIMPs in a terrestrial detector is expected to present a modulation with a period of one year due to the Earth revolution around the Sun [2].
The DAMA collaboration [3][4][5] has been measuring for more than 15 years a yearly modulation effect in their sodium iodide target. Such effect has a statistical significance of more than 9σ and is consistent with what is expected from DM WIMPs. However, in the most popular WIMP scenarios the DAMA modulation appears incompatible with the results from many other DM experiments that have failed to observe any signal so far.
This has lead to extend the class of WIMP models. In particular, one of the few phenomenological scenarios that have been shown to explain the DAMA effect in agreement with the constraints from other experiments is proton-philic spin-dependent inelastic Dark Matter (pSIDM) [6,7] for WIMP masses 10 GeV < ∼ m χ < ∼ 30 GeV and a mass splitting 10 keV < ∼ δ < ∼ 30 keV. Recently the DAMA collaboration has released first result from the upgraded DAMA/LIBRA-phase2 experiment [8]. Compared to the previous data the two most important improvements are that now the exposure has almost doubled and that the energy threshold has been lowered from 2 keV electron-equivalent (keVee) to 1 keVee. Moreover, an important difference with the result of DAMA/LIBRA-phase1 is that the new DAMA/LIBRA-phase2 spectrum of modulation amplitudes no longer shows a maximum, but is rather monotonically decreasing with energy. In light of these differences in the present paper we wish to update the assessment of pSIDM with the new DAMA/LIBRA-phase2 data, both in a scenario where the WIMP speed distribution f (v) is given by a standard Maxwellian and using a halo-independent approach where f (v) is not fixed.
In the present paper we will show that pSIDM can still provide a viable explanation of the modulation effect after DAMA/LIBRA-phase2. In particular, while the pSIDM scenario was able to explain DAMA/LIBRA-phase1 both for a Maxwellian f (v) and in a haloindependent approach [6,7] in the present paper we will show that for a Maxwellian WIMP velocity distribution it provides a poor fit to the new DAMA data and for a range of the pSIDM parameters in tension with the null results of other DM searches. On the other hand in a halo-independent approach the pSIDM scenario is still viable.
The paper is organized as follows. In Section II we outline the main features of the pSIDM scenario; in Section III A we analyze the DAMA data adopting a standard Maxwellian for the WIMP velocity distribution; in Section III B we analyze the DAMA data in a haloindependent approach. We provide our conclusions in Section IV.

II. THE PSIDM SCENARIO
The most stringent bounds on an interpretation of the DAMA effect in terms of WIMP-nuclei scatterings are obtained by detectors using xenon (XENON1T [9], PANDA [10], LUX [11]) and germanium (CDMS [12][13][14][15]) whose spin is mostly originated by an unpaired neutron while, on the other hand, both sodium and iodine in DAMA have an unpaired proton. This implies that if the WIMP particle interacts with ordinary matter predominantly via a spin-dependent coupling which is suppressed for neutrons it can explain the DAMA effect in compliance with xenon and germanium bounds [16,17]. Actually, present limits from xenon detectors require to tune the neutron/proton coupling ratio c n /c p to a small but non-vanishing value [6]. In the following we will adopt the xenon-phobic combination c n /c p =-0.028, that minimizes the xenon spin-dependent response using the nuclear structure functions in [18]. This scenario is still constrained by droplet detectors and bubble chambers (COUPP [19], PICASSO [20], PICO-60 [21])) which all use nuclear targets with an unpaired proton ( 19 F and/or 127 I). As a consequence, this class of experiments rules out a DAMA explanation in terms of WIMPs with a spin-dependent coupling to protons [6,17,22].
In Ref. [6] Inelastic Dark Matter [23] (IDM) was proposed to reconcile the above scenario to fluorine detectors. In IDM a DM particle χ 1 of mass m χ1 = m χ interacts with atomic nuclei exclusively by up-scattering to a second heavier state χ 2 with mass m χ2 = m χ + δ. A peculiar feature of IDM is that there is a minimal WIMP incoming speed in the lab frame matching the kinematic threshold for inelastic upscatters and given by: with µ χN the WIMP-nucleus reduced mass. This quantity corresponds to the lower bound of the minimal velocity v min (also defined in the lab frame) required to deposit a given recoil energy E R in the detector: with m N the nuclear mass. In particular, indicating with v * N a min and v * F min the values of v * min for sodium and fluorine, and with v esc the WIMP escape velocity, constraints from WIMP-fluorine scattering events in droplet detectors and bubble chambers can be evaded when the WIMP mass m χ and the mass gap δ are chosen in such a way that the hierarchy: is achieved, since in such case WIMP scatterings off fluorine turn kinematically forbidden while those off sodium can still serve as an explanation to the DAMA effect. So the pSIDM mechanism rests on the trivial observation that the velocity v * min for fluorine is larger than that for sodium.

III. ANALYSIS
The expected rate in a given visible energy bin 2 of a direct detection experiment is given by: with ǫ(E ′ ) ≤ 1 the experimental efficiency/acceptance. In the equations above E R is the recoil energy deposited in the scattering process (indicated in keVnr), while E ee (indicated in keVee) is the fraction of E R that goes into the experimentally detected process (ionization, scintillation, heat) and q(E R ) is the quenching factor, is the probability that the visible energy E ′ is detected when a WIMP has scattered off an isotope T in the detector target with recoil energy E R , M is the fiducial mass of the detector and T exp the live-time exposure of the data taking. In Eq.(4) the differential recoil rate dR χT (t)/dE ee is given by: where ρ WIMP is the local WIMP mass density in the neighborhood of the Sun (in the following we will assume the standard value ρ WIMP =0.3 GeV/cm 3 ), f ( v T , t) is the WIMP velocity distribution (whose boost in the Earth rest frame induces a time-dependence), N T the number of the nuclear targets of species T in the detector (the sum over T applies in the case of more than one target), while: with m T the mass of the nuclear target, j χ =1/2 the spin of the WIMP, E max R = 2µ 2 χT /m T v 2 T and σ 0 the point-like WIMP-nucleon cross section. In the following, for the calculation of the squared amplitude |M T | 2 we will use the spin-dependent nuclear form factors from [18] 1 for all nuclei with the exception of caesium and tungsten, for which we follow the same procedure adopted in Appendix C of [24].
In particular, in each visible energy bin DAMA is sensitive to the yearly modulation amplitude S m , defined as the cosine transform of with T 0 =1 year and t 0 =2 nd June, while other experiments put upper bounds on the time average S 0 :

A. Maxwellian analysis
In this Section we assume that the WIMP velocity distribution in the Galactic rest frame is a standard isotropic Maxwellian at rest, truncated at the escape velocity v esc , Here u is the WIMP speed in the Galactic rest frame, v 0 the galactic rotational velocity at the Earth's position, Θ is the Heaviside step function, and with z = v esc /v 0 . The WIMP speed distribution in the laboratory frame can be obtained with a change of reference frame. It depends on the speed of the Earth with respect to the Galactic rest frame, which neglecting the ellipticity of the Earth orbit, is given by In this formula, v ⊙ is the speed of the Sun in the Galactic rest frame, v ⊕ is the speed of the Earth relative to the Sun, and γ is the ecliptic latitude of the Sun's motion in the Galaxy. We take cos γ ≃ 0.49, v ⊕ = 2π(1 AU)/(1 year) ≃ 29 km/s, v ⊙ = v 0 + 12 km/s, v 0 = 220 km/s [25], and v esc = 550 km/s [26]. The velocity integral in Eq. (5) for the truncated Maxwellian distribution is computed from the expression of the speed distribution. We have obtained S 0 and S m by expanding it to first order in v ⊕ /v ⊙ .
where we consider 14 energy bins of width 0.5 keVee from 1 keVee to 8 keVee. Here and in the next Section we fix the experimental input (exposure, energy resolution, quenching factors, efficiency, measured count rates, etc.) for both the DAMA/LIBRA experiment and for other DM searches as described in appendix B of [24] and appendix A of [28].
The 5-σ best-fit DAMA region in the (m χ -σ 0 ) plane for the pSIDM scenario is compared to the corresponding 90% C.L. upper bounds from other DM searches in Fig. 2 (both in this Figure and in Fig. 3 the legend shows the actual list of experiments that we have tested and that yield some bound, although for readability purposes not all of them fall within the plot boundaries). In the same plot the IDM mass splitting is fixed to the absolute minimum of the χ 2 , δ=18.3 keV. As can be seen from such figure the DAMA effect is in strong tension with the upper bounds from PICO60, KIMS and PICASSO. We have also performed a combined fit including the upper bounds from such experiments with the addition of COUPP and XENON1T, finding χ 2 min =41.1 with a pvalue 1.5×10 −3 and 18 dof. Including v 0 and u esc as nuisance parameters in the χ 2 (we assume v 0 =(220± 20) km/s [25] and u esc =(550± 30) km/s [26]) does not improve the fit (we find χ 2 min =40.965). This confirms that, at variance with the analyses of Ref. [6,7], after the release of the DAMA/LIBRA-phase2 data the pSIDM scenario in the Maxwellian case is ruled out. There are two reasons for this. The first reason is that while the DAMA/LIBRA-phase1 data where only sensitive to scattering events off sodium, the DAMA/LIBRA-phase2 data have a lower threshold and are now also sensitive to scattering events off iodine for E ′ <2 keVee at low WIMP masses. This makes it more difficult to fit the model to the data since in the pSIDM scenario the scaling between the cross sections off iodine and sodium is fixed (the parameter c n /c p , that would allow to change such scaling is locked to the combination that suppresses the response on xenon). Moreover, in the scenario described in Section II a minimal value of the mass splitting parameter δ is required in order to comply with the condition of Eq. (3), which, at the same time automatically implies that the recoil energy E * R =E R (v * N a min ), and so a single maximum of the modulation amplitude spectrum, falls inside the range of the DAMA signal [6] (the energy E * R maximizes the velocity integral in Eq. (15)). Indeed, the DAMA/LIBRA-phase1 data showed a single maximum in the 2.5 keVee< E ′ <3 keVee energy bin in the measured modulation amplitudes [3,4], implying an acceptable fit for the pSIDM model. On the other hand the DAMA/LIBRA-phase2 data show an energy spectrum of the modulation amplitudes more compatible to a monotonically decreasing one, closer to what expected for elastic scattering. As a consequence of this the DAMA/LIBRA-phase2 χ 2 pulls to low values of the δ mass splitting (indeed, the Maxwellian best-fit configuration m χ = 12.1 GeV, δ=18.3 keV falls below the halo-independent compatibility region discussed in the next Section and shown in Fig. 4), entering in conflict with the requirement of Eq.(3) 2 .

B. Halo-independent analysis
In the halo-independent approach [29] the expected rate in a direct detection experiment is recast in the form [30]: where the dependence on astrophysics is contained in the halo function:η and the WIMP velocity distribution is contained in the velocity integral: while the response function (v min ) is given by: Notice that for a standard spin-dependent interaction the scattering amplitude in Eq.(6) does not depend on v T so the term v 2 T in the equation above cancels out in the product v 2 T dσ T /dE R . Due to the revolution of the Earth around the Sun, the velocity integralη(v min , t) shows an annual modulation that can be approximated by the first terms of a harmonic series, (17) with the only requirement that |η 1 | ≤η 0 . In this approach measured rates R i (with i = 0, 1) are mapped into suitable averages of the two halo functions η i . Averagesη i [vmin,1,vmin,2] (i = 0, 1) using R(v min ) in Eq. (16) as a weight function can then be directly obtained from the experimental data R i [30]: The result of such procedure is shown in Fig. 3 (v min ) to the left of v min,1 and to the right of v min,2 are each separately 16% of the total area under the function. This gives the horizontal width of the crosses corresponding to the rate measurements in Fig. 3. On the other hand, the horizontal placement of the vertical bar in the crosses corresponds to the average of v min , i.e., v min (vertical bar) The extension of the vertical bar shows the 1σ interval around the central value of the measured rate.
To compute upper bounds onη 0 from upper limits R lim on the unmodulated rates, we follow the conservative procedure in Ref. [29]. Sinceη 0 (v min ) is by definition a non-decreasing function, the lowest possibleη 0 (v min ) function passing through a point (v 0 ,η 0 ) in v min space is the downward step functionη 0 θ(v 0 −v min ). The maximum value ofη 0 allowed by a null experiment at a certain confidence level, denoted byη lim (v 0 ), is then determined by the experimental limit on the rate R lim The corresponding upper limits at 90% C.L. are shown as continuous lines in Figs. 3 for the same experiments shown in Fig. 2.
For the specific benchmark m χ =11.4 GeV, δ=23.7 keV shown in Fig. 3 one can see that pSIDM cannot be ruled out as an explanation of the DAMA/LIBRA effect since in all the energy range of the signal one has |η 1 [vmin,1,vmin,2] | ≪η lim . The same benchmark is represented by a starred point in Fig 4 and lies inside the the closed contour where the compatibility factor defined as [6]: is less than unity. In the equation above [v i min,1 , v i min,2 ] and v i 0 represent intervals and averages of v min for each of the i=1...14 DAMA/LIBRA bins below E ′ =8 keVee, while σ i is the 1-σ fluctuation onη 1 . In particular, the requirement D(m DM , δ) <1 ensures that within the closed contour of Fig. 4 no 1-σ interval of the quantitiesη 1 [vmin,1,vmin,2] obtained from the DAMA/LIBRA data lies completely above any of the upper boundsη lim .
From Fig. 4 one can see that in a haloindependent approach the pSIDM scenario can explain the DAMA/LIBRA data for 6 GeV < ∼ m χ < ∼ 17 GeV and 17 keV< ∼ δ < ∼ 29 keV. obtained from the DAMA/LIBRA data lies completely above any of the upper boundsη lim . The starred point corresponds to the benchmark shown in Fig. 3.

IV. CONCLUSIONS
We have shown that the Weakly Interacting Massive Particle scenario of proton-philic spin-dependent inelas-tic Dark Matter (pSIDM) can still provide a viable explanation of the observed DAMA modulation amplitude in compliance with the constraints from other experiments after the release of the DAMA/LIBRA phase-2 data. The pSIDM scenario provided a viable explanation of DAMA/LIBRA-phase 1 both for a Maxwellian WIMP velocity distribution and in a halo-independent approach. At variance with DAMA/LIBRA-phase1, for which the modulation amplitudes showed an isolated maximum at low energy, the DAMA/LIBRA-phase2 spectrum is compatible to a monotonically decreasing one. Moreover, due to its lower threshold, it is sensitive to WIMP-iodine interactions at low WIMP masses. Due to the combination of these two effects pSIDM can now explain the modulation observed by DAMA/LIBRA only when the WIMP velocity distribution departs from a standard Maxwellian. In this case the WIMP mass m χ and mass splitting δ fall in the approximate ranges 6 GeV < ∼ m χ < ∼ 17 GeV and 17 keV< ∼ δ < ∼ 29 keV.