Prospects for detection of target-dependent annual modulation in direct dark matter searches

Earth's rotation about the Sun produces an annual modulation in the expected scattering rate at direct dark matter detection experiments. The annual modulation as a function of the recoil energy $E_\text{R}$ imparted by the dark matter particle to a target nucleus is expected to vary depending on the detector material. However, for most interactions a change of variables from $E_\text{R}$ to $v_\text{min}$, the minimum speed a dark matter particle must have to impart a fixed $E_\text{R}$ to a target nucleus, produces an annual modulation independent of the target element. We recently showed that if the dark matter-nucleus cross section contains a non-factorizable target and dark matter velocity dependence, the annual modulation as a function of $v_\text{min}$ can be target dependent. Here we examine more extensively the necessary conditions for target-dependent modulation, its observability in present-day experiments, and the extent to which putative signals could identify a dark matter-nucleus differential cross section with a non-factorizable dependence on the dark matter velocity.


Introduction
Despite being the dominant form of matter in the Universe, the exact nature of the dark matter (DM) is still unknown. One of the most well-motivated candidates for DM is a particle with few GeV to hundreds of TeV mass and weak-scale interactions, referred to as a weakly interacting massive particle (WIMP). Efforts to shed light on the non-gravitational interactions of WIMP DM primarily focus on either detecting the byproducts of DM annihilation or decay (indirect detection), producing DM in the laboratory through collisions of standard model particles, or detecting interactions between DM in the galactic halo and terrestrial nuclei (direct detection).
Direct DM detection experiments attempt to gain insight into both the particle physics properties of DM and the local DM velocity distribution by observing the energy deposited by DM particles interacting with nuclei as they pass through detectors. A key feature of any convincing direct detection signal would be the annual modulation of the scattering rate caused by Earth's rotation around the Sun [1]. For DM velocity distributions that are locally smooth and isotropic in the galactic frame, it is usually expected that the differential rate for dark matter scattering off a target nuclide T is nearly sinusoidal and can be well represented by where E R is the nuclear recoil energy. Allowing the modulation amplitude S m (E R ) to assume both positive and negative values, the phase t 0 is independent of E R . Taking instead S m (E R ) to be non-negative, as we do in this paper, t 0 changes from early June at large E R to early December at small E R , with the transition occurring sharply at a single E R value. Accounting for the presence of anisotropy in the DM halo modifies this picture, most notably by modifying the E R dependence of the modulation phase. The extent to which various forms of anisotropy, including DM substructure, the gravitational focusing (GF) of DM particles by the Sun, and triaxial halo models, modify Eq. (1.1) has been investigated e.g. in [2][3][4][5][6][7][8][9][10][11][12][13]. At fixed recoil energies, experiments employing different target elements are not necessarily expected to measure the same modulation of the rate. However, for most interactions, some observables associated with the annual modulation like the modulation fraction or the time of maximum and minimum signal, t max and t min , do not depend on the target nuclide when expressed as functions of v min . This is the minimum speed a DM particle must have in Earth's frame to impart a recoil energy E R on a target nucleus. This definition naturally treats v min as an E R -dependent function. Alternatively, it is possible to think of E R as a v min -dependent function. In this context, E R is interpreted as the extremum energy (corresponding to a maximum energy if the scattering is elastic, and either a maximum or minimum energy if the scattering is inelastic) that can be imparted to a nucleus by an incoming DM particle traveling with speed v = v min in Earth's frame. For each nuclide there exists a bijective relation between E R and v min dictated by the scattering kinematics, and the choice of one or the other as the independent variable may lead to different insights. As commented above, for most interactions (e.g. the standard spin-independent (SI) and spin-dependent (SD) interactions) observables like t max and t min are nuclide-independent functions of v min (this is no longer true when expressed as functions of E R , since the E R -v min relation is target dependent). Therefore for studying the signal modulation for single-element targets it is convenient to adopt v min as the independent variable (averaging over different isotopes). For targets consisting of multiple elements, one must choose whether to treat E R or v min as the independent variable (see e.g. [14][15][16]). When we consider multiple targets in Sec. 3.1 we choose to return to using E R as the independent variable.
We pointed out in [17] that when the DM-nucleus differential cross section has a nonfactorizable velocity dependence, as for DM interacting through a magnetic dipole or an anapole moment, t max and t min are no longer target-independent functions of v min . Here, we reconsider the analysis performed in [17] and examine more extensively how target-dependent modulation arises, how various experiments can actually observe such a signal, and the extent to which putative signals could identify DM with a non-factorizable velocity dependence in its differential scattering cross section. Specifically, we consider how (i) integrating the scattering rate over a finite energy range, (ii) the presence of multiple target elements with non-negligible contributions to the rate, and (iii) different DM-nucleus scattering kinematics affect the potential observability of target-dependent modulation.
This paper is organized as follows. In Section 2 we introduce the formalism and discuss what conditions must be present for target-dependent modulation. In Section 3 we take the particular example of DM interacting with nucleons through an anomalous magnetic dipole moment and discuss how observables associated with the annual modulation of the rate depend on v min for specific targets employed in currents experiments. Additionally, we examine how experiments would view a signal arising from magnetic dipole DM a function of the observed energy E and the extent to which the expected signal would be distinguishable from a signal arising from a standard SI or SD contact interaction, for both elastic and inelastic scattering. We conclude in Section 4.

Direct detection rate
Direct DM detection experiments try to measure the recoil energy E R a nucleus initially at rest in the detector acquires after scattering with a DM particle with initial velocity v in the detector's rest frame. The differential scattering rate on a nuclide T per unit detector mass is where ρ = 0.3 GeV/cm 3 is the local DM density, m is the DM particle mass, C T is the nuclide mass fraction in the detector, m T is the target nuclide mass, and f (v, t) is the DM velocity distribution in Earth's frame. The energy dependence of v min (E R ), is dictated by the scattering kinematics, for instance for elastic scattering Experiments do not measure directly the recoil energy, but a proxy for it denoted here with E . This detected energy can e.g. be measured in keVee (keV electron-equivalent energy) or photoelectrons. For experiments that bin their data, the energy-integrated scattering rate between detected energies E 1 and E 2 is where (E ) is the counting efficiency and G T (E R , E ) describes the probability that an event detected with energy E resulted from a nuclear recoil having energy E R . G T (E R , E ) is frequently taken to be a Gaussian with mean value E = Q T E R , where Q T (E R ) is an element-dependent quenching factor. Typically one assumes the DM is on average at rest with respect to the galaxy, and the velocity distribution in the galactic frame f G (v) is smooth and isotropic. The DM velocity distribution in Earth's frame is then obtained via the Galilean transformation is the velocity of Earth with respect to the Sun and v the velocity of the Sun with respect to the galaxy. In this paper we choose to model the velocity of Earth with respect to the Sun following the procedure of Ref. [18], and take the velocity of the Sun with respect to the Galaxy to be v = (11, 232, 7) km/s in galactic coordinates. Furthermore, for concrete applications we assume the Standard Halo Model (SHM), in which the velocity distribution of the dark halo is a truncated Maxwellian, with galactic escape velocity v esc = 533 km/s [19] and velocity dispersion v 0 = 220 km/s [20]. The normalization, DM that is on average at rest with respect to the Galaxy has a preferred direction of motion in the Sun's reference frame. For this reason, DM particles viewed in the Sun's reference frame appear as a constant "wind", with velocities preferentially opposed to v . The gravitational potential of the Sun bends the trajectories of DM particles as they pass by, resulting in a focusing effect that is maximized at Earth's location when Earth is on the leeward side with respect to the Sun, occurring on March 1 st . This effect, referred to as GF, implies the DM density and velocity distribution at Earth's location depend on Earth's position relative to the Sun. The influence of GF is larger on slower moving particles as they spend more time in the Sun's gravitational potential, and is negligible on WIMPs traveling faster than a few hundred km/s in the Solar reference frame. The effect of GF is taken into account by replacing is the velocity a DM particle had asymptotically far away from the Sun's gravitational potential, such that its velocity when arriving at Earth is v [21]. Here u esc = 2GM /r ≈ 40 km/s is the escape velocity of the Solar System at Earth's location, r is the Sun-Earth distance,r is the unit vector pointing from the Sun to Earth, and v 2 ∞ = v 2 − u 2 esc .

Time dependence of the rate
For the commonly considered SI and SD contact interactions, the differential scattering cross section for DM-nucleus elastic scattering has the form where µ T is the DM-nucleus reduced mass, σ T is the total cross section for a point-like nucleus, and F T (E R ) is the appropriate nuclear form factor normalized as F T (0) = 1. This general form arises every time the scattering amplitude for a point-like nucleus is (at least approximately) independent of the scattering angle, i.e. of the recoil energy. In this case, T v 2 /m T is the maximum recoil energy a nucleus can get from scattering elastically with a DM particle with speed v. The effect of the finite size of the nucleus is then taken into account with the appropriate form factor. The differential rate for cross sections of the form in Eq. (2.7) then reads The modulation of the differential rate is determined solely by the time dependence in the velocity integral η 0 (v min , t), which is a target-independent function of v min , and therefore common to all experiments. Even though what enters the rate is the function η 0 (v min (E R ), t), which depends on the target through v min (E R ), one can express E R as a function of v min and study dR T /dE R (E R (v min ), t), which is proportional to the target-independent quantity η 0 (v min , t) (see e.g. [15,16]). The target-independent nature of the time dependence of the differential rate for the standard SI and SD contact interactions is a consequence of the fact that velocity and target dependence can be factored in the differential cross section shown in Eq. (2.7). One may then ask, in general, under what circumstances observables associated with the modulation of the rate, such as t max and t min , are target-dependent functions of v min . Following our preliminary study [17], we find that this can only happen when the following conditions are met: 1. the velocity and target dependence in the differential cross section cannot be factored, and 2. the scattering events that can be recorded by an experiment probe portions of the DM velocity distribution that are locally anisotropic in the galactic frame.
As shown in Ref [17], it is possible to meet both requirements and thus have a targetdependent modulation. Regarding point 2, anisotropy in the local DM velocity distribution can arise from an anisotropy in the smooth component of the halo, DM substructure, and gravitational interactions of DM with nearby massive objects such as the Sun. In this paper we choose to introduce anisotropy by only including the effect of GF of DM particles by the Sun because this anisotropy necessarily exists and is well understood [2,3,21].
Regarding point 1, the factorizable velocity and target dependence of the differential cross section, despite being very common, is not a completely general feature. The differential scattering cross section for DM interacting through a magnetic dipole [16, or an anapole moment [22,[43][44][45][46][47][48][49][50] actually contains two terms with unique velocity dependences and energy-dependent coefficients. These types of differential cross sections also appear with the interactions described by some of the effective operators studied e.g. in [51][52][53][54][55][56] (see [41,[57][58][59][60] for explicit formulas of scattering amplitudes). In all these examples, velocity dependences other than the dσ T /dE R ∝ 1/v 2 in Eq. (2.7) are present. This happens e.g. when higher order terms in the non-relativistic (small v) expansion of the scattering amplitude become important. To be concrete, we can take for example the scattering rate to be where we generalized the definition of the velocity integral in Eq. (2.10) to The interesting case for us is when r 0 and r 1 have similar magnitudes. The proportionality factor between r i and η i in Eq. (2.12) is in general E R dependent, and this dependence must balance the suppression provided by the extra powers of v in η 1 with respect to η 0 in order for r 0 and r 1 to be comparable. We will see below that the scattering rate of a DM particle interacting through an anomalous magnetic moment has exactly this form. As is clear from Eq. (2.11), the time dependence of the rate does not coincide in general with that of a single velocity integral, as it happened instead in the simple case of Eq. (2.9). It is therefore useful to denote with τ max (τ min ) the time of maximum (minimum) of each velocity integral, to distinguish it from the time of maximum (minimum) of the rate denoted t max (t min ). To understand the time-dependent behavior of η n we begin by considering the behavior of η 0 . The left panel of Fig. 1 shows η 0 evaluated at the first day of the month for the first six months of the year as a function of v min . Since the behavior of the curves is difficult to discern, we plot in the right panel of Fig. 1 the difference between each of the curves in the left panel and η 0 evaluated at March 1 st . Here, τ max , the time of maximum of the velocity integral, can be seen to transition from early January to early June as v min increases from ≈ 140 km/s to ≈ 260 km/s (actually, τ max occurs before January 1 st , during the month of December at low values of v min ). The inset in the right panel of Fig. 1 shows how η 0 changes with time should GF be neglected. Without GF, τ max still transitions from January 1 st to June 1 st , but this transition occurs very rapidly over a very narrow range of v min values. Fig. 2 shows η 1 as a function of v min for various fixed times. Unlike η 0 , there appears to be a fixed separation between the various fixed time curves across nearly all values of v min . This occurs because the additional factor of v 2 entering the velocity integral of η 1 weights the high velocity part of the spectrum, where the fixed time curves of η 0 are visibly separated. The inset of Fig. 2 zooms in on the low v min region to emphasize that τ max of η 1 does have a small v min dependence, transitioning from late May at small values of v min , to early June for v min 300.
For n > 1, one would expect the high end of the velocity distribution to become increasingly weighted, which within the SHM should result in a time dependence similar to that of η 1 , but even more independent of v min . This is shown in Fig. 3, where τ max and τ min are plotted for η 0 , η 1 , and η 2 . Instead of plotting τ min , we plot τ min −τ min , whereτ min ≡ τ max + 6 months. Fig. 3 shows the effect of including (solid) and neglecting (dashed) GF. 1 For η 2 , τ max is hardly affected by GF and thus only a single solid line is plotted. The results for τ min −τ min without GF are not shown as in this case τ min is nearly indistinguishable from  Figure 3. Time of maximum τ max (left) and minimum τ min (right) of η 0 , η 1 , and η 2 assuming the SHM, with (solid) and without (dashed) GF. In the right panel we plot τ min −τ min , whereτ min is τ max + 6 months. Neglecting GF, τ min is nearly indistinguishable fromτ min , and thus is not shown.
τ min . Fig. 3 shows that, within the SHM, η 0 is the only η n whose time-dependent behavior differs markedly from η n 1 . Thus, for the target-dependent features of the modulation to appear, assuming no other forms of anisotropy are present within the dark halo, the differential cross section must not only contain a non-factorizable velocity dependence, but one of the terms in the differential cross section must be proportional to η 0 . Should τ max and τ min of η 0 become v min dependent above 300 km/s, e.g. due to the presence of DM substructure [4], the approximate degeneracy of η n 1 (and near exact degeneracy of η n 2 ) would break and the previous requirement would no longer be necessary.
We would like to note that any η n can actually be rewritten in terms of, and thus computed from, which implies where we used the fact that η 0 (∞, t) = 0. With a similar set of manipulations, any arbitrary η n can be written in terms of any other arbitrary η n . Therefore, in principle, one may choose to express the rate in terms of any of the η n (or even in terms of f (v, t) itself, as shown in Eq. (18) of [16]). Some of the η n may have good properties for specific calculations, for example the normalization condition f (v, t) d 3 v = 1 can be written either as ∞ 0 η 0 (v min , t) dv min = 1 (see e.g. [61,62]) or η 1 2 (0) = 1. Moreover, whenever the velocity integrals need to be computed numerically (e.g. for complicated halo models, or when computing the effect of GF), Eq. (2.16) can be used to straightforwardly determine η n =0 once η 0 is known.
The different time dependence of the the various η n can be understood by looking at Eq. (2.16). Were it only for the first term on the right-hand side, η n =0 and η 0 would obviously have the same time dependence at fixed v min . Because of the second term, however, η n (v min , t) is a function of time that depends in a nontrivial way on η 0 (v, t) for all v v min .

Elastic scattering
We study here in detail the case of a Dirac fermion DM candidate χ elastically scattering with nuclei through an anomalous magnetic dipole moment λ χ , with interaction Lagrangian L = (λ χ /2)χσ µν χF µν . The differential cross section for elastic scattering off a target nuclide T with Z T protons and spin S T is with α = e 2 /4π the electromagnetic fine structure constant, m p the proton mass,λ T the nuclear magnetic moment in units of the nuclear magneton e/(2m p ) = 0.16 GeV −1 , and The differential cross section contains two terms, one arising from the charge-dipole interaction and the other arising from the dipole-dipole interaction. The former thus depends on the nuclear charge and contains a spin-independent form factor while the latter depends on the nuclear spin and contains a magnetic form factor. Both form factors are normalized to 1 at zero recoil energy. We compute the cross section with the formalism and form factors provided in [57,58].
Since the magnetic DM differential cross section contains terms proportional to η 0 (v min , t) and η 1 (v min , t), the modulation of the differential rate is a direct consequence of the interplay of these two functions and their respective coefficients. The relative importance of each of these functions is determined by the target and DM mass-dependent coefficients. We define r 0 (E R , t) and r 1 (E R , t) as the terms of the differential rate containing η 0 and η 1 respectively, as in Eqs. (2.11) and (2.12), andr 0 (E R ) andr 1 (E R ) to be their time average. The timeaveraged differential rate reads then dR T (E R )/dE R =r 0 +r 1 . Fig. 4 depicts the absolute value of the time-averaged rate fractions, as functions of v min , for six elements (fluorine, iodine, sodium, xenon, germanium, and argon) employed by current DM direct detection experiments. When more than one isotope is present, i.e. for germanium and xenon, r 0 and r 1 are understood to be summed over isotopes. Solid (dashed) lines correspond to a 100 GeV (1 TeV) DM particle.  The target dependence of t max and t min can be understood by combining the information on the time dependence of η 0 and η 1 in Fig. 3 with the information on the rate fraction of the corresponding element shown in Fig. 4. t max and t min as functions of v min are shown in Figs. 5-8 for magnetic DM scattering off fluorine, sodium, iodine, and xenon. We have chosen not to plot t max and t min for germanium and argon because the results for all DM masses  Figure 5.
Time of maximum t max (left) and minimum t min (right) of the differential rate for magnetic DM scattering off fluorine, plotted for various DM masses as a function of v min . The current low energy threshold for PICO has been mapped onto v min for each DM mass and is shown as a small solid dot. below 10 TeV are identical due to their small (germanium) or zero (argon) nuclear magnetic moment (see Ref. [17] for details). For each element, t max (left panels) and t min (right panels) are plotted for various DM masses ranging from 10 GeV to 10 TeV. Also shown, depicted as dots on the t max and t min curves, are the E R thresholds for LUX [63]  Figs. 5-8 show that t max and t min become target and DM mass independent for v min 300 km/s. This is due to the fact that the difference in the time-dependent behavior between η 0 and η 1 , which are central to the target-dependent features, vanish above v min ≈ 300 km/s (see Fig. 3), if the only source of anisotropy in the local halo is GF. Fig. 4 confirms that at sufficiently small values of v min the contribution to the differential rate from the term proportional to η 0 can be neglected. This is because the r 1 term contains the factor 1/v 2 min , which dominates the v min dependence of the rate at small v min values. Thus, in the small v min limit, t max occurs in late May and t min occurs in late November, regardless of the target element and DM mass. This behavior is a feature of elastic magnetic DM and other DM models could behave in a qualitatively different way.
For target elements that have a nonzero average nuclear magnetic moment (i.e., all elements considered here except argon), at large enough values of v min the dipole-dipole interaction inevitably becomes dominant, and thus r 0 > r 1 . This is because the spin-independent form factor in Eq. (3.1) decreases significantly faster than the magnetic form factor. Fig. 4 confirms that for all elements considered except argon, there exists a value of v min below which r 1 is the dominant contribution to the rate, and above which r 0 is the dominant contribution to the rate. The location in v min of this transition and how fast or gradual it is determine the unique element-dependent features of t max and t min in Figs. 5-8. The mass of the DM particle can have a large influence on the appearance of target-  Figure 6. Time of maximum t max (left) and minimum t min (right) of the differential rate for magnetic DM scattering off sodium, plotted for various DM masses as functions of v min . The current low energy threshold for DAMA has been mapped onto v min , assuming a quenching factor Q Na = 0.3, for each DM mass and is shown as a small solid dot. dependent features. Consider for instance the difference between a 100 GeV and 1 TeV DM particle scattering off xenon. For a 100 GeV DM particle, Fig. 4 shows that the v min point at which r 0 becomes dominant is around v min ≈ 400 km/s. Since this value of v min lies in the target-independent region, t max is effectively determined solely by the time dependence of η 1 . As the DM mass increases, the point at which r 0 becomes dominant with respect to r 1 shifts to lower values of v min . This is partly due to the fact that the v min value corresponding to a given E R decreases, but also because the terms 1/µ 2 T and µ 2 T /m 2 multiplying the SI component of Eq. (3.1) decrease. Consequently, for a 1 TeV DM particle scattering off xenon, the v min value at which r 0 becomes dominant appears in a v min region where the time dependence of η 1 and η 0 differ, leading to the appearance of a unique target-dependent feature in the t max and t min curves.
Up to this point we have only discussed how target-dependent modulation arises and how, under the assumption of magnetic DM, observables associated with the modulation of the rate in v min can potentially change. We have not yet discussed how these effects would manifest in present day experiments. To determine if experiments are capable of observing these target-dependent features, one must take into account the experimental threshold, the efficiency, the energy resolution, and the binning method.
The obvious requirement for these target-dependent effects to be observable, is that the experimental threshold in v min must be below 300 km/s. The threshold in v min depends on the threshold in E , the DM particle mass, and the scattering kinematics. Figs. 5 and 8 show that present thresholds are already low enough to give rise to a four month difference in t max for a 50 GeV DM particle scattering elastically off fluorine and xenon (while the 50 GeV curve is not shown for xenon, it directly overlaps with the 100 GeV curve), should the differential scattering rate be measured with perfect energy resolution, which is not possible for actual experiments.
Since we would like to see how observable this target dependence could be, we choose to consider experiments employing elements with large nuclear magnetic moments. For this Time of maximum t max (left) and minimum t min (right) of the differential rate for magnetic DM scattering off iodine, plotted for various DM masses as functions of v min . The current low energy threshold for DAMA has been mapped onto v min , assuming a quenching factor Q I = 0.09, for each DM mass and is shown as a small solid dot. reason we begin by considering the fluorine-based experiment PICO. PICO measures the energy-integrated rate as a function of threshold energy E th , and has an energy-dependent efficiency function that reduces the contribution of the scattering events near threshold. Figs. 5 and 9 can be used to understand how much the modulation features in the differential rate are erased in the energy-integrated rate. Fig. 9 depicts the time-averaged differential rate (summed over isotopic composition) for a 100 GeV DM particle scattering off fluorine, sodium, iodine, argon, germanium, and xenon, for magnetic DM as a function of E R . Fig. 9 includes both log-linear (left) and log-log (right) plots to show the different features of the spectra. If the differential rate were very steep, the integrated rate would be dominated by the differential rate at threshold, and thus have a similar annual modulation. As the differential rate flattens, an increasingly unweighted averaging occurs for all energies above threshold. The flattening of fluorine's differential rate occurs below PICO's 3.2 keV threshold, and thus the pronounced features appearing in t max of the differential rate should be strongly suppressed in the integrated rate. Fig. 10 depicts how PICO would realistically observe the time of maximum of the energyintegrated rate as a function of the threshold energy for a 100 GeV DM particle interacting through a magnetic dipole (solid blue line) or with the standard SI/SD contact interaction (dashed red line). As PICO does not provide an analytic form of their efficiency, we take the parametrization used by PICASSO, with α = 5 for fluorine [66]. We also assume a perfect energy resolution, G T (E R , E ) = δ(E R − E ). We have checked that the contribution from carbon is negligible for all energies so we consider only fluorine. Fig. 10 shows that the time of maximum of the rate as would be measured by PICO is nearly identical for the magnetic dipole interaction (dashed red line) and the standard SI/SD contact interactions (solid blue line), for all threshold energies we examined (larger than 0.1 keV). To determine if the two interactions could be differentiated by binning the data, we also consider a fluorine-based experiment capable of measuring the rate in 1 keV bins. For this hypothetical experiment we take the same efficiency function we used for PICO, and plot the result as horizontal bars in Fig. 10 for the magnetic dipole interaction (blue) and standard SI/SD contact interaction (red). The difference in the time of maximum of the energy-integrated rate for the two interactions in this hypothetical experiment ranges from 7 days to 20 days for threshold energies between 1 and 10 keV.
There are a number of reasons for the unique target-dependent features shown in Fig. 5 to be strongly suppressed when calculating the energy-integrated rate. First, the features in t max for the magnetic dipole interaction differ the most from the standard SI/SD contact interactions in the v min region where the r 0 and r 1 terms in Eq. (2.11) cross over. For fluorine, this occurs at very small v min values, v min 70 km/s. The top axis of Fig. 10 shows that this region of v min corresponds to very low energies, far below PICO's current threshold. Additionally, for elastic scattering E R ∝ v 2 min , and since the integration of the differential rate is over E R , the Jacobian's dependence on v min must be included in the integrand when performing the integral in v min instead. This additional factor increases the weight of the large v min region in the integration where the modulation is target independent. Finally, as previously mentioned, the differential rate decreases rather slowly as a function of E R , smearing the target-dependent features.
Let us see if other experiments could better preserve the target-dependent features. Let us consider DAMA/LIBRA, henceforth referred to as DAMA (or any of the upcoming DAMA-like experiment as KIMS-NaI, ANAIS, DM-Ice17, and SABRE, see e.g. [67,68] and references therein). DAMA is an interesting experiment to consider as both sodium and iodine have reasonably large nuclear magnetic moments and bin their data in small, 0.5 keVee, intervals. In the left panel of Fig. 11 we plot the time of maximum of the DAMA binned rate as a function of E for both the magnetic dipole interaction (blue) and the standard SI/SD contact interaction (red), assuming elastic scattering with a 100 GeV DM particle. Also depicted with a vertical dashed line is DAMA's current low energy threshold of  Figure 9. The time-averaged differential rate (summed over isotopes) in units of counts/(kg day keV) for a 100 GeV magnetic DM particle scattering off various elements as a function of recoil energy, shown on a semi-log (left) and log-log (right) plot. λ χ has been set to 10 −20 e cm.
2 keVee for the analysis of the modulated signal. The results for DAMA are calculated using quenching factors Q Na = 0.3 and Q I = 0.09, and a gaussian energy resolution function with standard deviation 0.448 √ E + 0.0091E [69]. The results for the two interactions are nearly indistinguishable above 4 keVee, and only differ by about a month in the lowest observable energy bin. It is worth mentioning that DAMA will soon extend their low-energy threshold down to 1 keVee which should result in a further observable difference between modulation arising from the standard SI/SD contact interactions and magnetic DM.
Like PICO, DAMA also sees a strong suppression in the target element dependent features of the modulation. The reason for the suppression in DAMA, however, is not primarily due to integrating over the differential rate, but rather due to the fact that DAMA has two non-negligible target elements. The independent contribution to the time-averaged differential rate from sodium (yellow) and iodine (green) as a function of detected energy is shown in the right panel of Fig. 11. Since each element has a different v min to E (average) mapping, and neither element dominates the differential rate in the 2-6 keVee range, the target-dependent region of t max for sodium is partially averaged with the target-independent region of iodine, leading to a large suppression of the target-dependent features. Furthermore, the small quenching factor of iodine pushes the most pronounced differences of the t max curve below threshold. The horizontal dashed lines in the left panel of Fig. 11 show how the v min values for sodium (yellow) and iodine (green) independently map to E , in average.
Since experiments do not know the DM particle mass or the scattering kinematics a priori, it is nontrivial to obtain t max as a function of v min from the data. For this reason, and because t max as a function of E R is necessarily known to be target element dependent, it is logical to ask how t max for magnetic DM differs from t max for the standard SI/SD contact interactions as functions of E R . This comparison is made in Fig. 12, where the left panel shows t max for SI/SD interactions while the right panel shows t max for magnetic DM, both as functions of E R . In both cases we assume a 100 GeV DM particle scattering elastically with various target elements (note that the curves for argon, germanium, and xenon in the right panel overlap almost entirely). For the standard contact interaction with only r 0 in the rate (see Eqs. (2.11) and (2.12)), as the SI/SD interaction, the differences in the curves is determined solely by the mass of the target nuclide. The largest difference in t max therefore occurs between fluorine and xenon and is around three months for recoil energies between 15 and 20 keV. While this is a rather large discrepancy, the shape of the t max curves for the standard SI/SD interactions are all stretched and compressed images of each other. In fact, all curves are obtained from the curve for η 0 in Fig. 3 with the E R -v min relation for elastic scattering in Eq. (2.2), which of course only differ in each case for the choice of m T . In this sense, t max and other observables associated with the modulation are not truly target dependent for interactions with only r 0 in the rate. The same cannot be said for magnetic DM. The right panel of Fig. 12 shows that the difference between various t max curves is more pronounced than when the standard interactions are considered, and furthermore, the curves have a more individualized shapes. The only exception are the curves for germanium, argon, and xenon, which completely overlap for a 100 GeV DM particle, a consequence of having a small or zero (for argon) average nuclear magnetic moment.

Inelastic scattering
Prior to this point we have only considered DM-nuclei elastic scattering. It has been shown that inelastic scattering, which can occur when there exist at least two DM particles with nearly degenerate masses m and m + δ with δ m, has the potential to significantly alter the scattering kinematics and the observed annual modulation [70,71]. The time-averaged differential event rate for a 100 GeV magnetic DM particle scattering off sodium (yellow) and iodine (green) as a function of detected energy.
Inelastic endothermic scattering occurs when the light DM state scatters into the heavy DM state, δ > 0. Since this process requires additional energy, only DM particles traveling at speeds greater than or equal to v T δ ≡ 2δ/µ T can scatter off a particular target T . If GF is the sole source of anisotropy, target-dependent modulation can only occur when speeds of about 200 km/s are probed. This implies that for a fixed DM mass, there exists a maximum mass splitting δ max for which target-dependent modulation can occur. For a 100 GeV DM particle scattering off fluorine, sodium, and iodine, this corresponds to values of δ max ≈ 3.3 keV, 4 keV, and 12 keV, respectively. These values of δ are quite small with respect to the typical momentum transfer in the interaction, and thus we expect the scattering kinematics to be almost elastic. Without an additional form of anisotropy, endothermic scattering is therefore ineffective in probing values of v min which can lead to target-dependent modulation.
Inelastic exothermic scattering, occurring when the heavier DM particle down-scatters into the lighter DM state (δ < 0), can be potentially more interesting for target-dependent modulation. To illustrate how exothermic scattering can alter the observed modulation, we plot in the left panel of Fig. 13 t max for DM interacting with various elements through the standard SI/SD contact interaction, assuming m = 100 GeV and δ = −10 keV, as a function or E R . This result is obtained by mapping the τ max (v min ) line corresponding to η 0 shown in Fig. 3 onto E R by using the E R -v min relation for inelastic scattering, (remember that for the SI/SD interaction t max coincides with τ max ). We have chosen not to plot t max for magnetic exothermic dark matter because, for all elements considered, the results mirror what would be expected should the differential cross section either be independent of velocity, or proportional to v −2 . That is to say for a given element, only the term proportional to η 0 or the term proportional to η 1 is relevant, never both. To understand why this is the case, it is necessary to first consider the differential cross section [25]: There are two additional terms with respect to the elastic case in Eq. (3.1), both contributing to the charge-dipole term for inelastic magnetic dark matter, one of which is proportional to E −1 R and the other to E −2 R . Both of these terms are contained within f 0 (see Sec. 3.1), and since the target dependence relies on the interplay between f 0 and f 1 , it is important to understand how these two new terms contribute to the relative rate fractions.
In Sec. 3.1, we showed that for elastic scattering f 1 is always the dominant contribution to the rate at low v min . This is a consequence of having a term proportional to v −2 min ∝ E −1 R . For inelastic magnetic DM, f 0 now has a term proportional to E −2 R , thus at very low energies r 0 is always the dominant contribution to the rate. This might be avoided, however, because there may exist a lower limit on E R which depends on v esc , and this may be above the region where E −2 R is the dominant factor (see Fig. 1 of [72]). At large energies, both of the new terms will be suppressed, and as for elastic scattering, the rate should be controlled by the term containing the magnetic form factor, r 0 (assuming the target element has a non-negligible nuclear magnetic moment). Whether r 0 or r 1 dominates the rate at intermediate energies depends strongly on the target element, the DM mass, and δ.
To illustrate how these variables affect the potential appearance of target-dependent modulation, we plot in the right panel of Fig. 13  show the terms proportional to f 0 and f 1 , respectively. The green region highlights values of E R where target-dependent modulation could potentially be observed (i.e. v min 200 km/s, assuming GF is the sole source of anisotropy), and the dot-dashed orange line depicts the energy corresponding to v min = 0 km/s. To compute the rate we again use the form factors provided in [57,58]. While these only apply to elastic scattering, [41] showed that they can be adapted to inelastic scattering by properly taking into account the modification to v ⊥ = v + q/2µ N , the component of v orthogonal to the momentum transfer q, due to inelastic kinematics (µ N being here the DM-nucleon reduced mass). Therefore one simply needs to replace the variable v ⊥ in the form factors of [57,58] with the true orthogonal component of the DM velocity for inelastic scattering, v ⊥ inel = v ⊥ + δq/|q| 2 . Two comments are in order. We previously stated that f 0 should be the dominant term at low values of E R due to the E −2 R term in the differential cross section. While this may not appear to be the case in Fig. 13, this is simply because we have not plotted the low E R regime, as it is not relevant for target-dependent modulation (low E R corresponds to large WIMP velocities where GF is unimportant). Next, for the current choice of parameters, f 0 is the only relevant term in the E R range where the effect GF is important, and thus the t max curve is identical to the fluorine curve shown in the left panel of Fig. 13. We stress that the unique target-dependent features seen in the t max and t min curves of Figs. 5-8 only arise if both f 0 and f 1 contribute in a non-negligible way within the region capable of probing low DM speeds.
It is interesting to see how changing m, δ, and the target element alter the results of Fig. 13. Changing the DM particle mass results in two distinct effects. Contrary to elastic scattering, lower values of m increase the importance of f 0 relative to f 1 at fixed E R , and thus the point at which f 0 becomes dominant relative to f 1 shifts to lower values of E R . The second and more import effect arises from changing the value of m in Eq. (3.4), which causes the E R range where the effect of GF is relevant in the right panel of Fig. 13 to shift. Using Eq. (3.4), one can see that decreasing the DM mass shifts the influence of GF to lower values of E R . We have checked that for δ 10 keV, lowering the DM particle mass to 10 GeV does not bring the point at which f 0 and f 1 cross into the region where target dependent modulation could occur.
Increasing the magnitude of δ (i.e. making δ more negative) also has two effects. First, it shifts the point at which f 0 and f 1 cross to higher values of E R . This effect is completely negligible, however, when compared with how this change in δ shifts the E R range where the effect of GF is relevant (see Eq. (3.4)).
The negligible nuclear magnetic moments of germanium, xenon, and argon lead to a complete dominance of f 1 over f 0 for essentially all values of E R , regardless of the DM mass and δ. This implies that inelastic magnetic DM scattering with these elements will always lead to an observation of t max between late May and early June, and the annual modulation will be consistent with inelastic scattering through differential cross sections that are independent of velocity. For iodine and sodium we have checked that the crossover from f 1 to f 0 as the dominant contribution to the rate, either always occurs far below threshold, or does not occur in the region where target-dependent modulation would arise. Identifying this type of scattering would then necessitate at least one experiment employing germanium, xenon, or argon, and another experiment employing fluorine, sodium, or iodine, to observe the annual modulation.

Identification of non-factorizable cross sections
The target-dependent effects described thus far have relied on two assumptions: experiments probe anisotropy in the DM halo and velocity and target dependence cannot be factored in the DM-nucleus differential scattering cross section. The question remains how a differential cross section of this form could be identified. A single experiment can never uniquely determine the underlying particle physics and astrophysics; it is only possible for a single experiment to say that their findings are consistent with some set of assumptions on the distribution of DM, the DM mass, a particular DM-nucleus interaction, etc. The most model-independent information is likely to come from a comparison of the outcomes of different experiments. We believe the most effective way to confirm the existence of a DM-nucleus cross section with a non-factorizable target and velocity dependence is to show that there exists no E Rv min relation capable of mapping observables associated with the modulation of the rate from experiments employing different target elements onto a unique function of v min . We emphasize however that finding unique functions of v min capable of reconciling the results of multiple experiments does not preclude the existence of non-factorizable differential cross sections. In the case of inelastic magnetic DM, elements with small average nuclear magnetic moments, e.g. germanium, xenon, argon, and carbon, will all yield similar results because the differential cross section is dominated by a single term, at least for the v min region where the local DM distribution is made anisotropic by GF.

Conclusions
It is typically assumed that observables associated with the annual modulation of the rate in direct detection experiments, when expressed as functions of v min (the minimum DM speed necessary to impart a given recoil energy to a target nucleus), are unique targetindependent functions. We have shown that this is not necessarily the case, and in fact the existence of a DM-nucleus differential cross section with a non-factorizable target and velocity dependence naturally leads to target-dependent modulation. The identification of this type of differential cross section is not straightforward and must be done through a process of elimination. In the event that multiple experiments with putative signals cannot find an E R -v min relation that can reconcile the differences between the observed modulations, one may then infer the potential existence of a non-factorizable differential cross section. We emphasize, however, that the reverse is not true. That is to say, finding an E R -v min relation that maps observables associated with the modulation from multiple experiments onto unique v min -dependent functions does not necessarily ensure that the modulation is target independent.
As a specific example, we have shown how t max (t min ), the time of maximum (minimum) of the differential rate, depends on the target nuclide for magnetic dipole DM elastically scattering with fluorine, germanium, iodine, sodium, argon, and xenon. We have also discussed how the annual modulation would appear should DM scatter inelastically with these elements. In our calculations we assume the SHM and included the effect of GF. We have shown that in an idealized experiment, the observed difference in t max for DM scattering off fluorine and xenon at a fixed value of v min could differ by as much as four months for DM masses above 50 GeV, however, accounting for the limitations of a realistic detector and integrating the differential rate can significantly suppress these differences. The plausible presence of DM substructure or forms of anisotropy other than GF could nevertheless enhance the target dependence of the modulation.