Dark Matter Interpretation of the Fermi–LAT Observation Toward the Galactic Center
The center of the Milky Way is predicted to be the brightest region of -rays generated by self-annihilating dark matter particles. Excess emission about the Galactic center above predictions made for standard astrophysical processes has been observed in -ray data collected by the Fermi Large Area Telescope. It is well described by the square of an NFW dark matter density distribution. Although other interpretations for the excess are plausible, the possibility that it arises from annihilating dark matter is valid. In this paper, we characterize the excess emission as annihilating dark matter in the framework of an effective field theory. We consider the possibility that the annihilation process is mediated by either pseudo-scalar or vector interactions and constrain the coupling strength of these interactions by fitting to the Fermi Large Area Telescope data for energies 1–100 GeV in the region about the Galactic center using self-consistently derived interstellar emission models and point source lists for the region. The excess persists and its spectral characteristics favor a dark matter particle with a mass in the range approximately from 50 to 190 (10 to 90) GeV and annihilation cross section approximately from 110 to 410 (610 to 210) cm/s for pseudo-scalar (vector) interactions. We map these intervals into the corresponding WIMP-neutron scattering cross sections and find that the allowed range lies well below current and projected direct detection constraints for pseudo-scalar interactions, but are typically ruled out for vector interactions.
Despite the overwhelming evidence from astrophysics and cosmology that roughly 80 of the matter in our Universe is in the form of dark, non-baryonic particles, how this so-called dark matter (DM) fits with the Standard Model (SM) of particle physics is currently unknown. Determining the nature of DM is one of the most pressing questions in the physical sciences, and a wide array of experiments are underway which hope to shed light on its identity by observing its interactions with the better understood particles of the SM.
Indirect detection is one of the promising avenues to elucidate the nature of DM. This method attempts to detect and discriminate the SM particles produced by DM particle annihilations (or decays) from those produced by conventional astrophysical processes. -rays of GeV energies are a particularly effective messenger because they propagate unhindered on galactic scales, and thus can be effectively traced back along the direction of their origin. In recent years, the Fermi Large Area Telescope (Fermi-LAT) has mapped out the -ray sky with the highest sensitivity of space-borne detectors to date, leading to the current best limits on the annihilation cross section for GeV DM annihilations that result in -rays.
Numerical simulations of galaxy formation offer clues as to where DM annihilation is expected to shine the most brightly. The simulations typically predict a large concentration of DM close to the Galactic center (GC), which smoothly falls off with Galactocentric radius. They also predict localized over-densities of DM, some of which correspond to dwarf spheroidal satellite galaxies. Both targets provide complementary regions of interest for DM searches. The DM related emission from the dwarf galaxies is expected to be of lower intensity, but to be relatively free of standard astrophysical backgrounds. Searches for -ray emission from dwarf satellites of the Milky Way have so far shown no convincing signal of DM annihilation Ackermann et al. (2015); Geringer-Sameth and Koushiappas (2011). In contrast, the GC is expected to produce a higher intensity annihilation signal. However, the region about the GC is strongly confused because of the intense interstellar emission and numerous discrete sources of -rays that are summed along and through the line-of-sight toward the GC. The estimation of these fore-/background contributions pose a significant challenge for detection of DM annihilation at the GC.
There seems to be an excess of -rays from the direction of the GC, above the expectations from astrophysics. This feature was first observed by Goodenough and Hooper Goodenough and Hooper (2009); Hooper and Goodenough (2011), and its general features, a spatial morphology remarkably consistent with predictions for a DM annihilation signal and a spectrum that peaks at a few GeV, persist in more recent analyses Hooper and Linden (2011); Abazajian and Kaplinghat (2012); Hooper and Slatyer (2013); Gordon and Macias (2013); Huang et al. (2013); Daylan et al. (2016); Abazajian et al. (2014); Zhou et al. (2015); Calore et al. (2015a); Abazajian et al. (2015); Calore et al. (2015b); Carlson et al. (2016). The Fermi-LAT collaboration has released its own analysis Ajello et al. (2016) of the -rays from the direction of the inner galaxy based on specialized interstellar emission models (IEMs) for estimating the fore-/background emissions, and enabling the analysis to make the first separation of the -ray emission of the kpc region about the GC from the rest of the Galaxy. Even with these IEMs, which represent the most sophisticated modeling to date, the excess persists. However, its spectral properties are strongly dependent on the assumed IEM, making it challenging to conclusively identify its origin. As a result, it remains unclear whether this signal arises from DM annihilation rather than from a currently unknown contribution from astrophysics such as a large population of milli-second pulsars, cosmic-ray (CR) proton or electron outbursts, additional cosmic ray sources, and/or emission from a stellar over-density in the Galactic bulge Hooper et al. (2013); Abazajian et al. (2014); Carlson and Profumo (2014); Petrovi? et al. (2014); Cholis et al. (2015a, b); Carlson et al. (2016); Macias et al. (2016). An interesting development is the use of statistical tools which indicate that GeV photons from the direction of the inner galaxy region show significantly more clustering than would be expected from Poisson noise from smooth components Lee et al. (2015); Bartels et al. (2016); Lee et al. (2016); McDermott et al. (2016). However, it remains difficult with the current models to disentangle whether this feature represents a property of the excess itself, or unmodelled variation in the background components Horiuchi et al. (2016).
While it is clearly premature to claim that the GeV excess represents a confirmed signal of DM annihilation, in this paper we extract the properties of the excess under the assumption that it does. We make simultaneous fits to the parameters of generic, realistic particle physics model of DM annihilation together with those defining the broad characterization of the possible fore-/backgrounds determined using the methodology of Ref Ajello et al. (2016). As a result, we can compare with the expectations for such models from direct searches for DM and colliders, finding that the null results of those searches play a significant role in shaping the allowed parameter space.
Our work is organized as follows. In Section II, we very briefly review the methodology of the Fermi-LAT analysis Ajello et al. (2016) to formulate realistic IEMs, which crucially define the fore- and backgrounds as well as the astrophysical contributions from the GC itself. This is followed in Section III by a revisitation of some of the most important morphological and spectral features of the signal: its centroid and whether there is evidence for two separate components with distinct morphologies and spectra. In Section IV, we define realistic flexible DM models described by effective field theories (EFTs), and perform a maximum likelihood (ML) fit to determine the ranges of their parameters capable of describing the excess together with the IEM parameters. We compare the ML regions of those models to direct and collider searches for DM in Section V. Section VI contains our conclusions and outlook.
Ii Interstellar Emission Model and Analysis
The analysis presented in this paper employs the same data as used by Ref Ajello et al. (2016): front converting events corresponding to the P7REP_CLEAN_V15 selection Ackermann et al. (2012), in the energy range 1-100 GeV, and with zenith angles less than 100. Exposure maps and the PSF for the pointing history of the observations were produced using the Fermi–LAT ScienceTools package (version 09-34-02)111Available at http://fermi.gsfc.nasa.gov/ssc/data/analysis. Events are selected from approximately 62 months of data, from 2008-08-11 until 2013-10-15. We note that for high statistics analyses such as the one presented here a notable difference is not expected in the results obtained with the P7REP_CLEAN_V15 data processing and those processed using Pass 8 Atwood et al. (2013); this is confirmed by several previous analyses Lee et al. (2016); Carlson et al. (2016); Weniger ().
ii.2 Interstellar Emission Models
The interstellar emission is the largest contribution to the -ray emission toward and through the line-of-sight toward the GC. To separate the contribution by the Galaxy between our location and the inner 1 kpc region about the GC, and that on the other side of the GC, specialized IEMs (four in total) were developed for the Ref Ajello et al. (2016) analysis. The methodology employed templates calculated using the well-known GALPROP CR propagation modeling code222A description of the GALPROP code is available at http://galprop.stanford.edu that were scaled to the data outside of the inner region about the GC. Under the assumption of Galactocentric azimuthal symmetry, these IEMs were used to estimate the fore-/background emission over the region, enabling the separation. Employing this prescriptive methodology ensures that minimal biases are introduced when fitting to the inner region. In addition, point source lists were developed for each IEM with the properties of the individual point sources obtained in a combined ML fit over the region. The construction of each IEM and its associated point-source list/model is a critical improvement over earlier works because the residual emission is strongly dependent on modeling both the over the region self-consistently.
The four distinct IEMs from Ref Ajello et al. (2016) are labeled:
OB stars, intensity-scaled
OB stars, index-scaled
The IEMs differ in the assumed distribution of the sources of CRs as tracing either the distributions of pulsars or OB stars; and in the procedure employed to scale the -ray intensity of the fore-/background components outside of the region to the data, either by scaling the normalization of the model templates for intensity-scaled IEMs, or scaling the normalization and spectral index (the latter only for gas-related templates interior to the solar circle) for the index-scaled IEMs. Notably, it was found that the data are compatible with a contribution from -rays from DM annihilation, and that the agreement between the data and the model significantly improves for all four IEMs when an additional component with a DM annihilation morphology is included in the fit.
ii.3 Analysis Procedure
We employ the procedure developed by the Fermi–LAT Collaboration in Ajello et al. (2016), which performs a ML fit of a model consisting of one of the four IEMs and its corresponding list of point sources to the data in the region. For each model, we include a DM annihilation contribution (described below) and perform the fit using the gtlike package of the Fermi–LAT ScienceTools. The results of the fit are the coefficients of the interstellar emission components from within the the innermost 1 kpc, as well as those describing the DM model under consideration. All point sources with a test statistic (defined as in Mattox et al. (1996)) are included in the model. Their fluxes and spectra are determined by iterative fits, with each iteration freeing the spectral parameters for a subset of point sources in order of decreasing TS.
Iii Morphology and Spectral Characteristics
The DM spatial distribution used in this paper is described in this section. Because Ajello et al. (2016) tested spatial templates fixed at the position of Sgr A* we investigate the possibility of an offset from this location by refitting the DM spatial distribution and scanning the ML grid about the GC. If a large offset is found, it might challenge a DM interpretation of the excess. For some IEMs the DM spectrum obtained by Ajello et al. (2016) extended beyond 10 GeV, but a dedicated study of the spatial distribution GeV was not made; this is also investigated in this section.
iii.1 Dark Matter Component
The results of numerical simulations for galaxy formation can broadly be described by the Navarro, Frenk, and White (NFW) profile Navarro et al. (1997):
For this analysis, we use a scale radius = 20 kpc, and corresponding to a local DM density = 0.4 . Two values for the inner slope of the DM distribution are considered, = 1, 1.2. The more cuspy distribution = 1.2 is motivated by the possibility of halo contraction due to the influence of baryons, which are typically not included in the simulations Diemand et al. (2005). The square of the NFW distribution is used as a template for DM annihilation, and we refer to it as the “NFW profile” (for ) or “NFW-c” (for = 1.2).
iii.2 NFW Centroid
The centroid of the Milky Way DM halo is conventionally centered at the location of Sgr A*. Because a large offset from this location might disfavor a DM interpretation, we verify that the centroid of the excess is sufficiently close. An offset between the centroid of the DM halo and Sgr A* as large as approximately 2 is consistent with numerical DM simulations, with the largest offsets tending to correlate with flatter central profiles Kuhlen et al. (2013); Lena et al. (2014). An offset in the centroid position was previously reported in Calore et al. (2015a); Linden et al. (2016), while other studies of the GC excess have found it to be consistent with Sgr A*.
We investigate the centroid position of the excess by scanning the ML for different locations near Sgr A*, for each of the four IEMs. A power-law with exponential cut-off is employed for the spectral model, following Ajello et al. (2016). The scan is performed by making the ML fit following Sec. II with the DM template centered at each point of a grid with spacing centered on Sag A*. The results of the scan are shown in Fig. 1, where the color scale shows the as a function of Galactic latitude and longitude. The intersections of the dotted grid lines correspond to the points where the likelihood is evaluated. The circle indicates the position of Sgr A*, and the triangle is the most likely position of the centroid for that IEM. We find that the centroid position is offset from Sgr A* for all four IEMs, with the Pulsars, index-scaled model displaying the largest offset, both in longitude (0.6) and latitude (0.2). The other three models prefer an offset only in longitude (within 0.4 up to the grid accuracy). Based on the scan, Sgr A* is not favored as the location of the NFW centroid for all four IEMs, however its position is roughly consistent with a DM interpretation for the GC excess and imperfections in the IEMs could plausibly introduce an offset. We therefore assume for the remainder of this paper that the DM distribution is centered at Sgr A*.
iii.3 Multiple Component Fit
|Fit components (1+2)||log L||2log L|
|NFW + NFW||-82870||844|
|NFW + Gas template||-82942||700|
|NFW + 1 Gauss||-82968||648|
|NFW + 2 Gauss||-82932||720|
|NFW + 5 Gauss||-82951||682|
|NFW + 10 Gauss||-82950||684|
|Fit components (1+2)||log L||2log L|
|NFW + NFW||-82972||914|
|NFW + Gas template||-83068||722|
|NFW + 1 Gauss||-83096||666|
|NFW + 2 Gauss||-83065||728|
|NFW + 5 Gauss||-83147||564|
|NFW + 10 Gauss||-83111||636|
Whether the high-energy tail ( GeV) of the GeV excess spectrum is related to that at lower energies remains an open issue. In Ajello et al. (2016), the excess emission above 10 GeV is most prominent in the intensity-scaled IEMs. For the index-scaled variants however, it is largely attributed to interstellar emission (see also Daylan et al. (2016)). The origin of the GeV excess has been previously investigated by several studies. In Horiuchi et al. (2016), the excess emission above 10 GeV is found to cut off in the innermost few degrees about the GC (unlike the excess at a few GeV) and therefore to have a different spatial morphology; secondary emission from unresolved millisecond pulsars is proposed as an interpretation. In Linden et al. (2016), the excess emission above 10 GeV is found to have a similar radial profile as the peak emission. Ref Linden et al. (2016) also discusses the interplay with the Fermi Bubbles, although the bubble morphology close to the Galactic plane is uncertain.
Here we investigate the morphology of the GeV excess emission present for the Pulsars and OB stars, intensity-scaled IEMs. We perform a ML fit over the 1-100 GeV energy range with two components to model the GC excess: an NFW template; and a second component that has either an NFW, gas, or a 2D gaussian (with half-width, half maximum of , , , or ) morphology. These are the same templates that were employed by Ajello et al. (2016). Six template combinations for the two intensity-scaled models are therefore tested. The spectrum for each template is modeled as a power law with an exponential cutoff function. The ML fit is performed iteratively, as described in section II, and the results are shown in Tables 1 and 2 for the Pulsars and OB stars, intensity-scaled IEM, respectively. The NFW + NFW combination is favored over all of the others considered, for both IEMs.
In Fig. 2 the differential fluxes integrated over the region for the two component fits, along with the fractional residuals, are shown for the Pulsars, intensity-scaled model. The contribution to the flux from each of the two spatial components and the IEM are shown, with the IEM broken down into the contributions from inverse Compton (IC), emission from the inner 1 kpc (“ring 1” in the legend), and from the point sources. For each of the six combinations we consider, the low energy excess is better described by an NFW morphology. The more peaked 2D gaussian templates ( and ) have spectra that peak in the few GeV energy range and cutoff at higher energies. Note that their contribution is always well below the contribution assigned to the NFW template. On the other hand, the spectra for the broader 2D gaussian templates ( and ) are more prominent at higher energies, suggesting that the high-energy tail of the GeV excess is consistent with an extended component in the region. The NFW morphology, which is peaked towards the GC and broadly extended in the region, is better suited to model the excess emission over the full energy range compared to the other options we have considered. However, due to the limitations of the IEMs together with the limited statistics at the higher energies, it is difficult to conclude decisively whether or not the high-energy tail is a true feature of the GC excess. Given the current preference for a single NFW morphology for both low and high energy components, we include the full energy range when comparing with the DM scenarios in Section IV below.
Iv Dark Matter Interpretation
In this section we fit the parameters of particle physics models of DM, together with the parameters describing the fore-/backgrounds, extracting a comprehensive DM interpretation of the GC excess. As described in more detail below, we employ a parameterization of the DM particle physics model which allows for distinct annihilation rates into up-type quarks, down-type quarks, and leptons. Our parametrization has more flexibility than the often-considered annihilation into a single channel of SM particles and, in this sense, is better able to capture a wider array of realistic particle physics models for DM annihilation than those typically used in indirect searches.
iv.1 EFT Description of Dark Matter Interactions
We consider two representative EFTs that describe the DM interactions with the SM fermions. These theories form part of a universal set of operators to which any theory of DM flows at low energies, well below the masses of the particles responsible for communicating between the SM and the dark matter Beltran et al. (2009); Cao et al. (2011); Beltran et al. (2010); Goodman et al. (2010, 2011); Kumar and Marfatia (2013). Such models have previously been considered to describe the GC excess Alves et al. (2014); Liem et al. (2016). More generalized constructions are employed here, and their parameters are fit together with the IEM parameters as described in Section II. Of course, models with light mediators are also interesting, and worthy of investigation in their own right Boehm et al. (2014); Abdullah et al. (2014); Martin et al. (2014); Berlin et al. (2014); Balzs and Li (2014); Ko and Tang (2015); Carpenter et al. (2016); Escudero et al. (2016). We leave exploration of such theories for future work.
Both of our considered EFTs are chosen such that they mediate -wave (velocity-unsuppressed) annihilation, because a -wave annihilation mechanism would require such strong interactions to overcome the innate suppression that it is likely to already be ruled out by direct and/or collider searches. We further restrict them to follow the principle of minimal flavor violation (MFV) D’Ambrosio et al. (2002), such that the most stringent constraints from flavor-violating observables are mitigated by small Yukawa interactions. We consider models containing either pseudo-scalar or vector Lorentz structures described by Lagrangians and (respectively, in the fermion mass basis),
where is the sum over fermion flavor with the indicated relative weighting of for the pseudo-scalar (vector) interaction types, as dictated by the leading terms consistent with MFV. The are parameters with dimensions of energy which specify the separate interaction strengths between the DM and up-type quarks, down-type quarks, and charged leptons. Together with the DM mass, , these coefficients specify the point in parameter space for the DM model. They represent generalizations (in that they allow the couplings of up-type and down-type quarks and leptons to vary independently) of the commonly considered interactions D4 and D5 used in DM searches via direct detection and at colliders Goodman et al. (2010).
iv.2 -ray Flux from Dark Matter Annihilation
where indicates averaging over the DM velocity profile, (1) for quarks (leptons) counts their color degrees of freedom, and is the appropriate for the fermion under consideration. The inclusive cross section for annihilation into up-type quarks, down-type quarks, and charged leptons is the sum of the individual cross sections for all three flavors of each fermion type, and the total cross section is the sum of the three inclusive cross sections. In presenting results, we typically trade the three parameters for and the fractional cross sections , , and (with ). It is easy to map these back into the parameters using the appropriate single channel cross section from Eqs. (4) and (5).
The -ray intensity and spectrum from DM annihilation is constructed by summing over all of the annihilation channels:
where is the number of rays per annihilation into the channel, generated from the PPPC 4 DM ID package Cirelli et al. (2011) based on fits to Pythia 8.1 Sjostrand et al. (2008), and = 2(4) for Majorana (Dirac) DM. The integral is the J-factor, obtained by integrating the DM density corresponding to either an NFW or NFW-c distribution, Eq. (1), over the line of sight () in direction .
To determine the preferred DM model parameters for each IEM, we fix the DM mass in the range from 10 – 250 GeV in 10 GeV increments. For each mass hypothesis the analysis procedure of Section II determines the fitted values of the DM model parameters , , and , along with the coefficients of the interstellar emission components from within the innermost 1 kpc and point sources, as usual. We repeat this scan for both NFW and NFW-c annihilation morphologies and for both the pseudo-scalar and vector models described above. We find that the DM component is detected with high statistical significance for all IEMs, and for pseudo-scalar as well as vector interactions. The likelihood values for pseudo-scalar interactions are summarized in Table 3.
|IEM||log L (null hypothesis)||log L (NFW)||log L (NFW-c)|
|OB stars, index-scaled||-82993||-82779||-82806|
|OB stars, intensity-scaled||-83429||-83081||-83117|
iv.3 Results for Pseudo-scalar Interactions
In Fig. 3, we display the likelihood profile as a function of the DM mass for each of the IEMs for the NFW-c annihilation morphology. The results for the NFW morphology are qualitatively similar. Each of the four IEMs shows a clear preference for particular DM masses, but there is considerable variation between them, with the index-scaled models favoring a mass around 50 GeV, while the intensity-scaled models favor higher masses 200 GeV. The results are consistent with the results obtained by Ajello et al. (2016), where the spectrum of the GC excess for the index-scaled IEMs displays a lower energy cutoff compared to the intensity-scaled IEMs. The spectra we consider here correspond to motivated DM scenarios, in contrast with the simpler assumptions made for the spectral model by Ajello et al. (2016).
In Fig. 4, we present the ML fractions into the three annihilation channels as a function of the DM mass, for each of the IEMs with the NFW-c annihilation morphology. These also vary considerably from one IEM to another, and are characterized by one channel or another typically dominating at any given DM mass hypothesis: charged leptons at lower masses GeV; down-type quarks in the range GeV; and up-type quarks above 180 GeV and at lower masses GeV. The lepton flux declines steeply above 20 GeV, and its contribution to the flux is smaller for the index-scaled models (Pulsars in particular) compared to the intensity-scaled ones. This reflects in part the lower energy cutoff of the GC excess spectrum for the index-scaled models and the harder -ray spectra produced by charged leptons compared to quarks. Also of note is the sharp transition from annihilation into down-type quarks to up-type quarks at the top mass threshold, GeV. This follows because the pseudo-scalar model annihilations are dominated by the heaviest quark kinematically accessible, and top quarks produced close to at rest decay into GeV bottom quarks, corresponding to the ML region at GeV.
The best-fit DM mass for the Pulsars (OB stars) index-scaled IEM is GeV ( GeV), and in both cases annihilation is predominantly into bottom quarks333The grid spacing is taken into account in the quoted uncertainties on the DM mass.. These results are compatible with the findings of previous studies Agrawal et al. (2015); Balzs and Li (2016) interpreting the spectrum of the excess as presented in Ref. Ajello et al. (2016). The intensity-scaled IEMs favor higher DM masses, GeV and GeV, for the Pulsars and OB stars variants, respectively, and primarily favor annihilation into top quarks. We note that the likelihood profile for the OB stars, intensity-scaled IEM is rather flat around the minimum, which yields a higher uncertainty in the best-fit DM mass, compared to the other IEMs. The uncertainties on the flux fractions into up-type and down-type quarks in this mass range are also somewhat larger.
The differential fluxes for the ML model (and the data points) are shown for each IEM in Fig. 5. Individual model components are displayed separately, including the contribution to the DM flux from each annihilation final state, as well as their sum. The contribution from each DM annihilation channel illustrates the fact that the integrated DM flux originates primarily from annihilations into quarks with the harder spectrum from annihilation into leptons becoming important at higher energies, particularly for the intensity-scaled IEMs. The -ray emission correlated with gas from the innermost 1 kpc is sub-dominant in the region. Fig. 5 also shows the fractional residuals as a function of energy. The agreement between data and model is at the level of a few % or better up to GeV for all IEMs, and is generally worse at higher energies for all but the Pulsars, index-scaled IEM. It is plausible that the energy cutoff at the DM mass in the annihilation spectrum limits its ability to describe the excess at the higher energies while simultaneously providing a good fit to the data in the few GeV range. We note that the fractional residuals based on realistic DM models including up-type, down-type, and lepton final states generally improve (for the same number of free parameters) over the results in Ajello et al. (2016) based on a power law with exponential cutoff spectrum.
Residual count (data-model) maps are shown in Fig. 6 for the energy bands , , and GeV, for each IEM. Structured excesses and deficits remain that may be attributed to imperfect modeling of the interstellar emission. Because of this, we do not rule out the DM models corresponding to IEMs with larger fractional residuals as these discrepancies might be explained by limitations in the IEMs. There is better agreement with the data when the DM spectrum is modeled with power law functions in 10 independent energy bins as done in Ajello et al. (2016); perhaps unsurprising given the larger number of free parameters for the spectral model.
The differential flux from the total DM annihilation component for both profiles (NFW, NFW-c) and all four IEMs are summarized in Fig. 7. The bands represent the 1 fit uncertainty on the flux summing the up-type, down-type, and lepton final states. For the index-scaled variants of the IEMs, the spectrum peaks at a few GeV, while for the intensity-scaled counterparts the peak shifts to higher energies. This is consistent with the requirement that the high energy tail in the spectrum for the intensity-scaled IEMs, predominantly from annihilations into leptons, has to cutoff at the same energy (corresponding to the DM mass) as the contribution to the flux from annihilations into up-type and down-type quarks, which dominate the DM flux at lower energies. Finally, we note that the flux for NFW-c profile is smaller compared to the NFW profile. As a consequence, a simple rescaling based on factors when comparing fit results obtained with different profiles is not accurate, as the flux assigned to the DM component has a dependence on the specific morphology.
We translate the DM template flux for each IEM into the inclusive annihilation cross section, with the results shown in Fig. 8. Also shown for comparison is the predicting saturation the measured DM relic density for a standard cosmology Steigman et al. (2012). The results for the index-scaled models are comparable to those found in most of the earlier studies of the GeV excess Hooper and Goodenough (2011); Abazajian and Kaplinghat (2012); Hooper and Slatyer (2013); Gordon and Macias (2013); Huang et al. (2013); Daylan et al. (2016); Abazajian et al. (2014); Zhou et al. (2015); Calore et al. (2015a); Abazajian et al. (2015); Calore et al. (2015b); Carlson et al. (2016). The intensity-scaled models however are consistent with larger DM masses and cross sections, as first discussed in Agrawal et al. (2015), based on the spectra from Ajello et al. (2016).
iv.4 Results for Vector Interactions
The analysis for the vector-type DM interactions proceeds very similarly to the analysis of the pseudo-scalar interactions described above. For each IEM and both NFW and NFW-c morphologies, the DM mass is scanned and the couplings to up-type quarks, down-type quarks, and charged leptons is fit. The results are presented in Figs. 9 and 10, respectively, for each IEM with the NFW-c profile (the results for the NFW profile are qualitatively similar.) Similarly to pseudo-scalar interactions, lower DM masses are favored by the index-scaled IEMs, compared to the intensity-scaled. However, in general, lower DM masses are favored for the vector interaction models than for the pseudo-scalar ones for the same IEM. In addition, because the coupling to SM fermions is assumed to be flavor-universal for the vector interaction model, there is no sharp transition in behavior at the top quark mass. For the Pulsars, index-scaled IEM, there are two close-to-degenerate minima in the likelihood profile, with the lower mass dominated by annihilations into leptons444For annihilations into leptons, secondary -ray emission via IC processes is neglected. Note that for DM masses GeV, IC photons are mainly produced at energies GeV Cirelli et al. (2013a); Lacroix et al. (2014). . The fitted values of and the DM mass for each of the IEMs and DM profiles are shown in Fig. 11.
V Comparison with Other Searches
As seen in sections IV.3 and IV.4, DM interpretations of the GC excess cover a broad range of masses ( GeV) and , depending on the IEM, DM profile, and interaction type. One crucial avenue toward exploring a DM hypothesis for the excess is to compare the regions of parameter space best describing the excess with the results from other searches for DM. Null results of such searches can sharpen the target parameter space or even exclude candidate explanations, whereas positive results could strengthen a DM interpretation of the excess and better define the characteristics of candidate models.
v.1 Indirect Searches
For masses in the range GeV, the strongest constraints from indirect detection are generally from Fermi–LAT observations of dwarf spheroidal galaxies Ackermann et al. (2015). These limits appear to constrain the region relevant for explanations of the GC excess, but are derived from less theoretically motivated DM annihilation models where the DM annihilates into one species of SM fermion at a time. As such, they do not precisely apply to the models considered here, although similar conclusions are likely. The bound based on the assumption of annihilation into , corrected to account for Dirac (rather than Majorana) DM particles, is shown on Figures 8 and 11 for reference. The dwarf spheroidal bounds for annihilations into leptons are not displayed in these figures. Although they would in principle be more pertinent to constrain our low mass, vector interaction results, they are still not adequate as the final state channel we consider here is an equal weight mixture of , , and therefore not directly comparable.
The limitations in the IEMs, modeling uncertainties in the dwarf halos Hayashi et al. (2016); Bernal et al. (2016); Ichikawa et al. (2016); Klop et al. (2016), modifications to the particle physics model for DM Kopp et al. (2016), and large uncertainties in the J-factor for the GC Abazajian and Keeley (2016), all widen the relative uncertainties when confronting the parameters describing the GC excess with the limits from observations of dwarf spheroidal galaxies. Because of this, care must be taken when contrasting these limits with a DM interpretation of the GC excess.
The particle physics models under consideration also lead to annihilations producing anti-matter, such as positrons or anti-protons. Positrons in particular show excess production compared to naive expectations Adriani et al. (2009); Aguilar et al. (2013), leading to limits which do not significantly constrain the parameters for the GC excess Bergstrom et al. (2013). Recently, Ref. Cuoco et al. (2016) performed a detailed analysis of the anti-proton spectrum measured by AMS-02 Aguilar et al. (2016), and also found an indication for an excess component roughly consistent with the parameter space describing a DM interpretation of the GC excess (see also Giesen et al. (2015) for a less optimistic view). The interpretation of CR anti-matter measurements is complicated by propagation, energy losses, and other modeling uncertainties related to particle fragmentation, as well as the spatial distribution of astrophysical sources. Consequently, the interpretation of these data in terms of DM is unclear.
v.2 Direct Searches
Coupling to quarks implies coupling to hadrons, and thus is bounded from direct searches for DM scattering with heavy nuclei. Models with pseudo-scalar interactions map onto a scattering cross section which is both suppressed by the small velocities of DM in the Galactic halo and are also spin-dependent. As a result, the expectation is that the constraints from direct searches result in mild constraints. In contrast, vector interactions lead to velocity-unsuppressed spin-independent scattering and are strongly constrained by direct searches. For the vector models, which contribute to the spin-independent cross section , we follow the usual convention mapping onto this quantity defined at zero relative velocity. For pseudo-scalar interactions, we compute the integrated cross section for DM scattering with a nucleon by integrating over the recoil energy of the nucleus and the velocity of the DM, which we assume follows a Maxwellian distribution, using techniques developed in Fitzpatrick et al. (2013); Anand et al. (2014); Cirelli et al. (2013b); Gresham and Zurek (2014), (specifically using the code presented in Ref. Cirelli et al. (2013b)). This integrated cross section should be distinguished from usual spin-dependent cross section , defined at zero velocity scattering, and is a more appropriate measure of scattering which is strongly velocity-dependent.
In Figs. 12 and 13, we show the ML points for the pseudo-scalar and vector models mapped into the WIMP-neutron spin-dependent integrated cross section, respectively, for each IEM and both NFW and NFW-c. For comparison, the limits from the LUX search for DM scattering with Xenon are presented Akerib et al. (2016), also mapped into or the integrated cross section for spin-dependent scattering with neutrons. For the vector models, the limits from LUX easily exclude all of the ML points except for the point with dark matter masses around GeV which annihilates predominantly into leptons for the Pulsars, index-scaled IEM with NFW-c profile, which has sufficiently small coupling to quarks that the scattering with nuclei is highly suppressed. For the pseudo-scalar models, the predictions for the ML points lie well below below the LUX bounds, with the lower mass points potentially probed long-term by Darwin Aalbers et al. (2016), while the higher mass points are slightly above the neutrino floor Billard et al. (2014) and out of the reach of these experiments. These results illustrate the importance of the IEM modeling and its influence on characterization of the putative signal, which can lead to drastic differences in the expectations from complementary searches.
v.3 Collider Searches
Searches at the Large Hadron Collider (LHC) are more model-dependent, and can be classified based on the masses and couplings of the particles mediating the interaction. When such particles are heavy compared to the typical collider energies, they can be described by the same EFTs employed in this paper. The results of searches in this regime are typically not competitive with direct searches except at masses far below those of interest to describe the GC excess Aad et al. (2013); Khachatryan et al. (2015). For lighter mediating particles, the limits depend sensitively on the specific couplings to the DM as well as to the SM fermions. In particular, for values of the cross sections similar to what has been found in past characterizations of the GeV excess, cases where a pseudo-scalar mediator’s coupling to DM is significantly weaker than the coupling to quarks are mildly constrained by LHC data, and the opposite limit is essentially unconstrained Fan et al. (2016). Given the wide range of parameter space (which is even larger for the specialized IEM analysis considered here), it seems possible that the LHC could eventually hope to observe an excess consistent with a pseudo-scalar mediator interpretation if parameters are favorable. Similar remarks apply to the vector mediator models, although all but the Pulsars, index-scaled IEM with NFW-c profile are already excluded by direct detection experiments. This latter model is consistent with vanishing coupling to quarks, and thus is unlikely to be excluded by searches at the LHC.
The excess of GeV -rays from the direction of the GC is an indication that there is something in the -ray sky beyond our current knowledge. Whether this source ultimately proves to originate from DM annihilation or from a more conventional astrophysical source still remains to be determined, and is likely to require further experimental input. As part of this process, we have examined key aspects of the putative signal using the specialized IEMs, developed by the Fermi–LAT Collaboration Ajello et al. (2016). Our goal in characterizing potential DM explanations is to explore the implications from complementary searches, which can rule out or favor a DM interpretation.
Our results illustrate the impact of interstellar emission modeling on the extracted characteristics of the excess and highlight the need for improved modeling to capture a more realistic range of possibilities. As far as the gross characteristics of the excess are concerned, we find an offset of of the excess centroid from Sgr A* for all four IEMs considered. We further find no significant evidence that the tail of the excess has a different spatial morphology than the few GeV bump, with both high energy and low energy components favoring an NFW morphology compared to the other morphologies we have considered.
We also consider flexible and realistic particle physics models for DM interacting with up-type quarks, down-type quarks, and charged leptons, for two separate interaction types (pseudo-scalar and vector) leading to -wave annihilation. These theories are described by EFTs, valid when the momentum transfer is small compared to the masses of the particles mediating the interactions – to describe annihilation, this implies the mediators are heavier than the DM itself. We find that the choice of IEM has a large impact on the preferred DM mass, annihilation cross section, and primary annihilation channel. In particular, we identify regions with higher masses and annihilation predominantly into top quarks. Comparing the ML points in parameter space with direct and collider searches, we find that all of the vector models aside from one at DM mass GeV and annihilating into leptons are ruled out by null results from the LUX experiment. The pseudo-scalar models predict spin-dependent and velocity-dependent scattering with nuclei at a rate far below the current sensitivity, but in some cases within the grasp of future planned experiments. It would be interesting, but beyond the scope of this work, to extend our analysis beyond the EFT limit to the case of models where the DM can annihilate directly into the mediator particles themselves.
The GeV excess is a compelling hint that there is more to learn about the Galaxy. It is likely to take a combined effort of observation and interpretation to unravel its nature.
The authors are pleased to acknowledge conversations with D. Finkbeiner, D. Hooper, M. Kaplinghat, T. Slatyer, and C. Weniger. The work of CK and SM is supported in part by Department of Energy grant DE-SC0014431. The work of TMPT and PT is supported in part by National Science Foundation grants PHY-1316792 and PHY-1620638. GALPROP development is partially funded via NASA grants NNX09AC15G, NNX10AE78G, and NNX13AC47G.
- Ackermann et al. (2015) M. Ackermann et al. (Fermi-LAT), Phys. Rev. Lett. 115, 231301 (2015), arXiv:1503.02641 [astro-ph.HE] .
- Geringer-Sameth and Koushiappas (2011) A. Geringer-Sameth and S. M. Koushiappas, Phys. Rev. Lett. 107, 241303 (2011), arXiv:1108.2914 [astro-ph.CO] .
- Goodenough and Hooper (2009) L. Goodenough and D. Hooper, (2009), arXiv:0910.2998 [hep-ph] .
- Hooper and Goodenough (2011) D. Hooper and L. Goodenough, Phys. Lett. B697, 412 (2011), arXiv:1010.2752 [hep-ph] .
- Hooper and Linden (2011) D. Hooper and T. Linden, Phys. Rev. D84, 123005 (2011), arXiv:1110.0006 [astro-ph.HE] .
- Abazajian and Kaplinghat (2012) K. N. Abazajian and M. Kaplinghat, Phys. Rev. D86, 083511 (2012), [Erratum: Phys. Rev.D87,129902(2013)], arXiv:1207.6047 [astro-ph.HE] .
- Hooper and Slatyer (2013) D. Hooper and T. R. Slatyer, Phys. Dark Univ. 2, 118 (2013), arXiv:1302.6589 [astro-ph.HE] .
- Gordon and Macias (2013) C. Gordon and O. Macias, Phys. Rev. D88, 083521 (2013), [Erratum: Phys. Rev.D89,no.4,049901(2014)], arXiv:1306.5725 [astro-ph.HE] .
- Huang et al. (2013) W.-C. Huang, A. Urbano, and W. Xue, (2013), arXiv:1307.6862 [hep-ph] .
- Daylan et al. (2016) T. Daylan, D. P. Finkbeiner, D. Hooper, T. Linden, S. K. N. Portillo, N. L. Rodd, and T. R. Slatyer, Phys. Dark Univ. 12, 1 (2016), arXiv:1402.6703 [astro-ph.HE] .
- Abazajian et al. (2014) K. N. Abazajian, N. Canac, S. Horiuchi, and M. Kaplinghat, Phys. Rev. D90, 023526 (2014), arXiv:1402.4090 [astro-ph.HE] .
- Zhou et al. (2015) B. Zhou, Y.-F. Liang, X. Huang, X. Li, Y.-Z. Fan, L. Feng, and J. Chang, Phys. Rev. D91, 123010 (2015), arXiv:1406.6948 [astro-ph.HE] .
- Calore et al. (2015a) F. Calore, I. Cholis, and C. Weniger, JCAP 1503, 038 (2015a), arXiv:1409.0042 [astro-ph.CO] .
- Abazajian et al. (2015) K. N. Abazajian, N. Canac, S. Horiuchi, M. Kaplinghat, and A. Kwa, JCAP 1507, 013 (2015), arXiv:1410.6168 [astro-ph.HE] .
- Calore et al. (2015b) F. Calore, I. Cholis, C. McCabe, and C. Weniger, Phys. Rev. D91, 063003 (2015b), arXiv:1411.4647 [hep-ph] .
- Carlson et al. (2016) E. Carlson, T. Linden, and S. Profumo, Phys. Rev. D94, 063504 (2016), arXiv:1603.06584 [astro-ph.HE] .
- Ajello et al. (2016) M. Ajello et al. (Fermi-LAT), Astrophys. J. 819, 44 (2016), arXiv:1511.02938 [astro-ph.HE] .
- Hooper et al. (2013) D. Hooper, I. Cholis, T. Linden, J. Siegal-Gaskins, and T. Slatyer, Phys. Rev. D88, 083009 (2013), arXiv:1305.0830 [astro-ph.HE] .
- Carlson and Profumo (2014) E. Carlson and S. Profumo, Phys. Rev. D90, 023015 (2014), arXiv:1405.7685 [astro-ph.HE] .
- Petrovi? et al. (2014) J. Petrovi?, P. D. Serpico, and G. Zaharija?, JCAP 1410, 052 (2014), arXiv:1405.7928 [astro-ph.HE] .
- Cholis et al. (2015a) I. Cholis, D. Hooper, and T. Linden, JCAP 1506, 043 (2015a), arXiv:1407.5625 [astro-ph.HE] .
- Cholis et al. (2015b) I. Cholis, C. Evoli, F. Calore, T. Linden, C. Weniger, and D. Hooper, JCAP 1512, 005 (2015b), arXiv:1506.05119 [astro-ph.HE] .
- Macias et al. (2016) O. Macias, C. Gordon, R. M. Crocker, B. Coleman, D. Paterson, S. Horiuchi, and M. Pohl, (2016), arXiv:1611.06644 [astro-ph.HE] .
- Lee et al. (2015) S. K. Lee, M. Lisanti, and B. R. Safdi, JCAP 1505, 056 (2015), arXiv:1412.6099 [astro-ph.CO] .
- Bartels et al. (2016) R. Bartels, S. Krishnamurthy, and C. Weniger, Phys. Rev. Lett. 116, 051102 (2016), arXiv:1506.05104 [astro-ph.HE] .
- Lee et al. (2016) S. K. Lee, M. Lisanti, B. R. Safdi, T. R. Slatyer, and W. Xue, Phys. Rev. Lett. 116, 051103 (2016), arXiv:1506.05124 [astro-ph.HE] .
- McDermott et al. (2016) S. D. McDermott, P. J. Fox, I. Cholis, and S. K. Lee, JCAP 1607, 045 (2016), arXiv:1512.00012 [astro-ph.HE] .
- Horiuchi et al. (2016) S. Horiuchi, M. Kaplinghat, and A. Kwa, (2016), arXiv:1604.01402 [astro-ph.HE] .
- Ackermann et al. (2012) M. Ackermann et al. (Fermi-LAT), Astrophys. J. Suppl. 203, 4 (2012), arXiv:1206.1896 [astro-ph.IM] .
- Atwood et al. (2013) W. Atwood et al. (Fermi-LAT) (2013) arXiv:1303.3514 [astro-ph.IM] .
- (31) C. Weniger, private communication .
- Mattox et al. (1996) J. R. Mattox et al., Astrophys. J. 461, 396 (1996).
- Navarro et al. (1997) J. F. Navarro, C. S. Frenk, and S. D. M. White, Astrophys. J. 490, 493 (1997), arXiv:astro-ph/9611107 [astro-ph] .
- Diemand et al. (2005) J. Diemand, M. Zemp, B. Moore, J. Stadel, and M. Carollo, Mon. Not. Roy. Astron. Soc. 364, 665 (2005), arXiv:astro-ph/0504215 [astro-ph] .
- Kuhlen et al. (2013) M. Kuhlen, J. Guedes, A. Pillepich, P. Madau, and L. Mayer, Astrophys. J. 765, 10 (2013), arXiv:1208.4844 [astro-ph.GA] .
- Lena et al. (2014) D. Lena, A. Robinson, A. Marconi, D. J. Axon, A. Capetti, D. Merritt, and D. Batcheldor, Astrophys. J. 795, 146 (2014), arXiv:1409.3976 [astro-ph.GA] .
- Linden et al. (2016) T. Linden, N. L. Rodd, B. R. Safdi, and T. R. Slatyer, (2016), arXiv:1604.01026 [astro-ph.HE] .
- Beltran et al. (2009) M. Beltran, D. Hooper, E. W. Kolb, and Z. C. Krusberg, Phys. Rev. D80, 043509 (2009), arXiv:0808.3384 [hep-ph] .
- Cao et al. (2011) Q.-H. Cao, C.-R. Chen, C. S. Li, and H. Zhang, JHEP 08, 018 (2011), arXiv:0912.4511 [hep-ph] .
- Beltran et al. (2010) M. Beltran, D. Hooper, E. W. Kolb, Z. A. C. Krusberg, and T. M. P. Tait, JHEP 09, 037 (2010), arXiv:1002.4137 [hep-ph] .
- Goodman et al. (2010) J. Goodman, M. Ibe, A. Rajaraman, W. Shepherd, T. M. P. Tait, and H.-B. Yu, Phys. Rev. D82, 116010 (2010), arXiv:1008.1783 [hep-ph] .
- Goodman et al. (2011) J. Goodman, M. Ibe, A. Rajaraman, W. Shepherd, T. M. P. Tait, and H.-B. Yu, Nucl. Phys. B844, 55 (2011), arXiv:1009.0008 [hep-ph] .
- Kumar and Marfatia (2013) J. Kumar and D. Marfatia, Phys. Rev. D88, 014035 (2013), arXiv:1305.1611 [hep-ph] .
- Alves et al. (2014) A. Alves, S. Profumo, F. S. Queiroz, and W. Shepherd, Phys. Rev. D90, 115003 (2014), arXiv:1403.5027 [hep-ph] .
- Liem et al. (2016) S. Liem, G. Bertone, F. Calore, R. Ruiz de Austri, T. M. P. Tait, R. Trotta, and C. Weniger, JHEP 09, 077 (2016), arXiv:1603.05994 [hep-ph] .
- Boehm et al. (2014) C. Boehm, M. J. Dolan, C. McCabe, M. Spannowsky, and C. J. Wallace, JCAP 1405, 009 (2014), arXiv:1401.6458 [hep-ph] .
- Abdullah et al. (2014) M. Abdullah, A. DiFranzo, A. Rajaraman, T. M. P. Tait, P. Tanedo, and A. M. Wijangco, Phys. Rev. D90, 035004 (2014), arXiv:1404.6528 [hep-ph] .
- Martin et al. (2014) A. Martin, J. Shelton, and J. Unwin, Phys. Rev. D90, 103513 (2014), arXiv:1405.0272 [hep-ph] .
- Berlin et al. (2014) A. Berlin, P. Gratia, D. Hooper, and S. D. McDermott, Phys. Rev. D90, 015032 (2014), arXiv:1405.5204 [hep-ph] .
- Balzs and Li (2014) C. Balzs and T. Li, Phys. Rev. D90, 055026 (2014), arXiv:1407.0174 [hep-ph] .
- Ko and Tang (2015) P. Ko and Y. Tang, JCAP 1501, 023 (2015), arXiv:1407.5492 [hep-ph] .
- Carpenter et al. (2016) L. M. Carpenter, R. Colburn, J. Goodman, and T. Linden, Phys. Rev. D94, 055027 (2016), arXiv:1606.04138 [hep-ph] .
- Escudero et al. (2016) M. Escudero, A. Berlin, D. Hooper, and M.-X. Lin, (2016), arXiv:1609.09079 [hep-ph] .
- D’Ambrosio et al. (2002) G. D’Ambrosio, G. F. Giudice, G. Isidori, and A. Strumia, Nucl. Phys. B645, 155 (2002), arXiv:hep-ph/0207036 [hep-ph] .
- Cirelli et al. (2011) M. Cirelli, G. Corcella, A. Hektor, G. Hutsi, M. Kadastik, P. Panci, M. Raidal, F. Sala, and A. Strumia, JCAP 1103, 051 (2011), [Erratum: JCAP1210,E01(2012)], arXiv:1012.4515 [hep-ph] .
- Sjostrand et al. (2008) T. Sjostrand, S. Mrenna, and P. Z. Skands, Comput. Phys. Commun. 178, 852 (2008), arXiv:0710.3820 [hep-ph] .
- Agrawal et al. (2015) P. Agrawal, B. Batell, P. J. Fox, and R. Harnik, JCAP 1505, 011 (2015), arXiv:1411.2592 [hep-ph] .
- Balzs and Li (2016) C. Balzs and T. Li, (2016), arXiv:1609.00113 [hep-ph] .
- Steigman et al. (2012) G. Steigman, B. Dasgupta, and J. F. Beacom, Phys. Rev. D86, 023506 (2012), arXiv:1204.3622 [hep-ph] .
- Cirelli et al. (2013a) M. Cirelli, P. D. Serpico, and G. Zaharijas, JCAP 1311, 035 (2013a), arXiv:1307.7152 [astro-ph.HE] .
- Lacroix et al. (2014) T. Lacroix, C. Boehm, and J. Silk, Phys. Rev. D90, 043508 (2014), arXiv:1403.1987 [astro-ph.HE] .
- Hayashi et al. (2016) K. Hayashi, K. Ichikawa, S. Matsumoto, M. Ibe, M. N. Ishigaki, and H. Sugai, Mon. Not. Roy. Astron. Soc. 461, 2914 (2016), arXiv:1603.08046 [astro-ph.GA] .
- Bernal et al. (2016) N. Bernal, L. Necib, and T. R. Slatyer, (2016), arXiv:1606.00433 [astro-ph.CO] .
- Ichikawa et al. (2016) K. Ichikawa, M. N. Ishigaki, S. Matsumoto, M. Ibe, H. Sugai, and K. Hayashi, (2016), arXiv:1608.01749 [astro-ph.GA] .
- Klop et al. (2016) N. Klop, F. Zandanel, K. Hayashi, and S. Ando, (2016), arXiv:1609.03509 [astro-ph.CO] .
- Kopp et al. (2016) J. Kopp, J. Liu, T. R. Slatyer, X.-P. Wang, and W. Xue, (2016), arXiv:1609.02147 [hep-ph] .
- Abazajian and Keeley (2016) K. N. Abazajian and R. E. Keeley, Phys. Rev. D93, 083514 (2016), arXiv:1510.06424 [hep-ph] .
- Adriani et al. (2009) O. Adriani et al. (PAMELA), Nature 458, 607 (2009), arXiv:0810.4995 [astro-ph] .
- Aguilar et al. (2013) M. Aguilar et al. (AMS), Phys. Rev. Lett. 110, 141102 (2013).
- Bergstrom et al. (2013) L. Bergstrom, T. Bringmann, I. Cholis, D. Hooper, and C. Weniger, Phys. Rev. Lett. 111, 171101 (2013), arXiv:1306.3983 [astro-ph.HE] .
- Cuoco et al. (2016) A. Cuoco, M. Krmer, and M. Korsmeier, (2016), arXiv:1610.03071 [astro-ph.HE] .
- Aguilar et al. (2016) M. Aguilar et al. (AMS), Phys. Rev. Lett. 117, 091103 (2016).
- Giesen et al. (2015) G. Giesen, M. Boudaud, Y. Gnolini, V. Poulin, M. Cirelli, P. Salati, and P. D. Serpico, JCAP 1509, 023 (2015), arXiv:1504.04276 [astro-ph.HE] .
- Fitzpatrick et al. (2013) A. L. Fitzpatrick, W. Haxton, E. Katz, N. Lubbers, and Y. Xu, JCAP 1302, 004 (2013), arXiv:1203.3542 [hep-ph] .
- Anand et al. (2014) N. Anand, A. L. Fitzpatrick, and W. C. Haxton, Phys. Rev. C89, 065501 (2014), arXiv:1308.6288 [hep-ph] .
- Cirelli et al. (2013b) M. Cirelli, E. Del Nobile, and P. Panci, JCAP 1310, 019 (2013b), arXiv:1307.5955 [hep-ph] .
- Gresham and Zurek (2014) M. I. Gresham and K. M. Zurek, Phys. Rev. D89, 123521 (2014), arXiv:1401.3739 [hep-ph] .
- Akerib et al. (2016) D. S. Akerib et al., (2016), arXiv:1608.07648 [astro-ph.CO] .
- Aalbers et al. (2016) J. Aalbers et al. (DARWIN), JCAP 1611, 017 (2016), arXiv:1606.07001 [astro-ph.IM] .
- Billard et al. (2014) J. Billard, L. Strigari, and E. Figueroa-Feliciano, Phys. Rev. D89, 023524 (2014), arXiv:1307.5458 [hep-ph] .
- Aad et al. (2013) G. Aad et al. (ATLAS), JHEP 04, 075 (2013), arXiv:1210.4491 [hep-ex] .
- Khachatryan et al. (2015) V. Khachatryan et al. (CMS), Eur. Phys. J. C75, 235 (2015), arXiv:1408.3583 [hep-ex] .
- Fan et al. (2016) J. Fan, S. M. Koushiappas, and G. Landsberg, JHEP 01, 111 (2016), arXiv:1507.06993 [hep-ph] .