Constraints on dark matter scattering with long lived mediators from observations of the Sun with the Fermi Large Area Telescope

The Sun is among the potential targets for indirect dark matter (DM) searches. DM particles from the Milky Way halo can undergo elastic interactions with the solar nuclei, in which they progressively lose energy until they are trapped by the solar gravitational field and sunk into the core of the Sun. DM particles captured inside the Sun can then self-annihilate into standard model (SM) particles that, with the exception of high-energy neutrinos, are likely re-absorbed in the interior of the Sun and cannot reach the Earth. However, in recent works [1, 2, 3, 4, 5, 6, 7] an alternative scenario was proposed, in which pairs of DM particles $\chi$ annihilate into pairs of long-lived mediators $\phi$, which can escape from the Sun and then decay, yielding SM particles in the final states, such as gamma rays or electron-positron pairs, that can reach the Earth and be detected.

In our previous works [8, 9] we analyzed the data collected by the Large Area Telescope (LAT) onboard the Fermi Gamma-ray Space Telescope mission [10] during its first 10 years of operation to search for possible signals from solar DM in the gamma-ray and in the electron-positron channels, and we were able to derive constraints on the DM-nucleon scattering cross sections. These constraints were obtained under the assumption that gamma rays or electron-positron pairs are directly produced in the mediator decays, i.e. in the processes $\phi\rightarrow\gamma\gamma$ or $\phi\rightarrow e^{+}e^{-}$. Under these assumptions, the resulting spectra of gamma rays and of electrons-positrons from DM are expected to exhibit a box-like feature, with the upper edge of the box corresponding to the mass $m_{\chi}$ of the candidate DM particle [11].

In this work we have updated the dataset of Ref. [8], including the data collected by the LAT in its first 13.5 years of operation. In addition, we have implemented a more general analysis approach, which can be also applied to search for possible solar DM signals from WIMP annihilations into long-lived mediators decaying into several different channels with gamma rays in the final states. In all these scenarios, the possible solar DM gamma-ray signals extend over a wide energy range, and the evaluation of the energy spectra requires a dedicated simulation, which is discussed in Sec. 2.

We have studied the spectra of gamma rays from a small Region of Interest (RoI) centered on the Sun and we have evaluated the possible signal counts from DM following two different approaches discussed in Sec. 3. The constraints on the DM signal intensities have been evaluated by requiring that the expected counts from DM annihilations in any channel do not exceed the possible signal counts.

In this work we are considering a physics scenario with pairs of DM particles $\chi$ annihilating inside the Sun into pairs of long-lived mediators $\phi$, which can escape from the Sun and then decay into channels with gamma rays in the final states. These photons can be produced in the decay channel $\phi\rightarrow\gamma\gamma$ or as final state radiation (FSR) through the decay channels $\phi\rightarrow b\bar{b},\tau^{+}\tau^{-},\mu^{+}\mu^{-},\ldots$

In the case of the $\phi\rightarrow\gamma\gamma$ decays, in Ref. [8] we have shown that, for a light mediator (i.e. $m_{\phi}\ll m_{\chi}$), the gamma-ray spectrum at Earth can be well modeled as a simple box, with its upper edge at an energy corresponding to the mass of the DM particle $m_{\chi}$. Similary, in the case of FSR, the gamma-ray spectrum extends over an energy range up to $m_{\chi}$. In all these cases we have evaluated the expected gamma-ray fluxes at Earth using the WimpSim²²2http://wimpsim.astroparticle.se/ version 5.0 Monte Carlo simulation code [12].

The DM gamma-ray flux expected at Earth is given by:

where $E$ is the gamma-ray energy, $m_{\chi}$ and $m_{\phi}$ are the masses of the DM particle $\chi$ and of the mediator $\phi$ respectively, $\sigma$ is the elastic DM-nucleon scattering cross section and $L$ is the boosted mediator decay length. The term $\Gamma_{\text{cap}}(m_{\chi},\sigma)$ in the r.h.s. of the previous equation is the DM capture rate in the Sun, which depends on $m_{\chi}$ and $\sigma$, while the term $\varphi(E;m_{\chi},m_{\phi},L)$ is the flux at Earth of gamma rays per DM annihilation, which depends on the masses $m_{\chi}$ and $m_{\phi}$ and on the mediator decay length $L$. This term can be expressed as:

where $D=1\hbox{$\,{\rm AU}$}$ is the Sun-Earth distance and $Y(E;m_{\chi},m_{\phi},L)$ is the gamma-ray yield per DM annihilation, which depends only on the kinematics of the mediator decay³³3The mass $m_{\chi}$ of the parent DM particle sets the energy of the mediator in the lab frame, since the annihilation process $\chi\chi\rightarrow\phi\phi$ is assumed to happen with both DM particles nearly at rest.. Eq. 2.1 has been written under the assumption of equilibrium between the DM capture and annihilation processes [9].

The DM capture rate is evaluated using the DARKSUSY code version 6.1.0 [13, 14, 15], under the assumptions of the default settings, with a local DM density $\rho_{\odot}=0.3\hbox{$\,{\rm GeV/cm^{3}}$}$, a Maxwellian velocity distribution with average velocity $v_{\odot}=20\hbox{$\,{\rm km/s}$}$ and a velocity dispersion $v_{rms}=270\hbox{$\,{\rm km/s}$}$. The calculation has been performed for both the spin-dependent and spin-independent cases. We have evaluated the capture rate using for the reference cross section $\sigma_{0}=10^{-40}\hbox{$\,{\rm cm^{2}}$}$. The capture rates corresponding to different cross sections can be easily obtained, since $\Gamma_{\text{cap}}$ is proportional to $\sigma$. Full details on the evaluation of the capture rates are given in ref. [9].

As mentioned above, to evaluate the gamma-ray yields, we have used the WimpSim code, that can simulate the production of neutrinos, gamma rays and other particles from DM annihilations in the Sun and in the Earth. The DM annihilations are simulated with the event generator Pythia-6.4.26⁴⁴4http://home.thep.lu.se/~torbjorn/Pythia.html [16], which is also used to simulate the subsequent decay chains of the annihilation products. In our simulations we used the med_dec module of WimpSim, which was specifically developed for describing the propagation and the decays of long-lived mediators. The mediators propagate out from the solar core, travelling along straight lines. The distance travelled before the decay is sampled from an exponential distribution, with an average boosted decay length $L$ set by the user.

In this work we have studied several mediator decay channels, performing a dedicated simulation campaign, in which we have scanned a wide range of possible values of the masses $m_{\chi}$ and $m_{\phi}$. For each pair of values $(m_{\chi},m_{\phi})$, a set of $10^{5}$ events have been simulated. In our simulations we used med_dec with the assumption $L=R_{\odot}$; the gamma-ray flux for different values of $L$ can be evaluated by taking into account the survival probability of the mediator in its journey from the Sun to the Earth [9]. After each simulation run, summary tables are generated with the phenomenological description of the events. In addition, event files are also generated, which contain the relevant information about the gamma rays produced in each event. The information stored in these files is used to evaluate the differential flux at Earth of gamma rays.

In Fig. 1 the outputs of some simulation runs are shown. The plots show the expected gamma-ray fluxes per DM annihilation from $1\hbox{$\,{\rm TeV}$}$ DM particles annihilating into mediators of different masses, which in turn can decay into the $b\bar{b}$, $\tau^{+}\tau^{-}$ and $\mu^{+}\mu^{-}$ channels. The gamma-ray fluxes at Earth have been evaluated either without any selection on the arrival directions of gamma rays or selecting only those photons with arrival directions lying within $2\hbox{${}^{\circ}$}$ from the Sun. This choice reflects the data selection implemented in the analysis (see Sec. 3). As can be seen from the examples shown in Fig. 1, this selection reduces the intensity of the possible DM signal at low energies, i.e. in the region where $E/m_{\chi}<10^{-3}$. This is due to the fact that the angular dispersion of photons is $\sim 1/\gamma_{\phi}$, where $\gamma_{\phi}$ is the Lorentz factor of the mediator. Assuming that DM particles annihilate at rest, $\gamma_{\phi}=m_{\chi}/m_{\phi}$ and, therefore, the heavier the mediator, the higher will be the angular dispersion and the lower will be the fraction of photons within a cone of $2\hbox{${}^{\circ}$}$ aperture from the Sun.

As mentioned above, in this work we have analyzed a data sample collected by the Fermi LAT during its first 13.5 years of operation. We have implemented an event selection based on our previous work [8] and we have developed a novel analysis procedure, which allows us to constrain both scenarios with direct and indirect gamma-ray production in the mediator decays.

In Figs. 4 and 5 the upper limits at $95\%$ CL on the spin-dependent ($\sigma_{SD}$) and spin-independent ($\sigma_{SI}$) DM-nucleon cross sections are shown as functions of the mass $m_{\chi}$ of the DM particles and for the mediator decay channels $\phi\rightarrow b\bar{b}$, $\phi\rightarrow\tau^{+}\tau^{-}$ and $\phi\rightarrow\mu^{+}\mu^{-}$. These constraints have been evaluated for a set of mediator masses $m_{\phi}$ in the kinematically allowed range, under the assumptions of a mediator decay length $L=R_{\odot}$. In the plots we show the limits evaluated with the two analysis approaches discussed in sec. 3. As expected, we see that the constraints on the cross sections obtained with the optimistic approach are in general stronger than those obtained with the conservative approach. The differences between the limits obtained with the two analysis approaches are typically within two orders of magnitude.

For all the channels investigated in this analysis the limits on the cross sections $\sigma_{SD}$ and $\sigma_{SI}$ are in the ranges from $10^{-45}$ up to $10^{-39}$ $\,{\rm cm^{2}}$ and from $10^{-47}$ up to $10^{-42}$ $\,{\rm cm^{2}}$ respectively. The strongest constraints are obtained for the lowest allowed values of $m_{\phi}$, i.e. when the mediator is highly boosted ($m_{\phi}\ll m_{\chi}$). We remark here that the limits on the DM-nucleon scattering cross section depend on the calculation of the reference capture rate. As discussed in Sec. 2, for this calculation we have used the DARKSUSY code with the default settings for the parameters describing the local DM density and its velocity distribution. Any variations of these parameters would result in a change of the reference capture rate, and consequently of the limits on the scattering cross sections. The values of the default parameters of DARKSUSY are widely used in the literature, and a study of the variations of the capture rate with these parameters would go beyond the goals of the present work.

As mentioned above, in our analysis we have considered DM masses in the kinematically allowed range. However, for masses below a few $\,{\rm GeV}$, the process of DM evaporation in the Sun interior can be non-negligible [31], and should be taken into account in the evaluation of the limits.

The first three panels of Fig. 6 show the upper limits obtained in this work on the spin-dependent DM-nucleon scattering cross section $\sigma_{SD}$ in the scenarios in which the mediator can decay in the channels $b\bar{b}$, $\tau^{+}\tau^{-}$ and $\mu^{+}\mu^{-}$. The limits shown in Fig. 6 have been evaluated assuming a mediator decay length $L=R_{\odot}$ and choosing for $m_{\phi}$ the lowest kinematically allowed value in each channel. We have also repeated our analysis procedure using the gamma-ray fluxes evaluated using med_dec without requiring the maximum angular separation of $2\hbox{${}^{\circ}$}$ between the arrival direction of gamma rays at Earth and the direction of the Sun (dashed lines in the figure). We find that if the DM fluxes are evaluated without any cut on the arrival directions of photons at Earth, the constraints on the cross sections are stronger. This feature is clearly evident in the region of low DM masses ($m_{\chi}<1\hbox{$\,{\rm TeV}$}$) and is due to the fact that, as discussed in Sec. 2, the cut on the arrival directions reduces the intensity of the expected DM gamma-ray signal in the low-energy region. Since for low DM masses the limits on the DM signal intensity are basically set by the limits on the gamma-ray counts in the low-energy bins ($E<10\hbox{$\,{\rm GeV}$}$), a lower number of expected counts from DM gamma rays in those bins will result in an increased limit on the DM signal intensity and consequently on the cross section.

In Fig. 6 we also compare the limits on $\sigma_{SD}$ obtained in our analysis with those published in refs. [32, 33], obtained from other analyses of the Fermi LAT and of the HAWC data, and with the constraints at $90\%$ CL from the direct measurements performed by the PICO-60 experiment [34]. Our limits are comparable to other measurements for a DM mass range up to $\sim 1\hbox{$\,{\rm TeV}$}$.

Finally, we have also reconsidered the scenario which we studied in our previous work [8], with the mediator decaying directly into pairs of gamma rays. The last panel of Fig. 6 shows the upper limits on $\sigma_{SD}$ evaluated in this scenario compared with those published in ref. [8]. The limits obtained with the optimistic analysis approach are consistent with those of our previous work up to DM masses of $\sim 100\hbox{$\,{\rm GeV/c^{2}}$}$, while they become consistent with those obtained with the conservative approach for DM masses in the $\,{\rm TeV}$ range.

In Fig. 7 we summarize the constraints obtained in this work with the optimistic analysis approach for both the spin-dependent and spin-independent DM-nucleon scattering cross sections $\sigma_{SD}$ and $\sigma_{SI}$. The limits shown in the figure have been obtained by studying several possible mediator decay channels ($b\bar{b}$, $\tau^{+}\tau^{-}$, $\mu^{+}\mu^{-}$, $\gamma\gamma$, $gg$, $WW$, $ZZ$ and $t\bar{t}$). We have assumed a mediator decay length $L=R_{\odot}$ and for each channel we have plotted the limits obtained by selecting the lowest kinematically allowed value of $m_{\phi}$. Our results are compared with those obtained from the direct measurements performed by the XENON1T [35] and PICO-60 [34] experiments respectively. In the case of the spin-dependent cross section, the limits obtained in this work are a few orders of magnitude stronger than those obtained by PICO-60 in the whole DM mass range explored, while in the case of the spin-independent cross section the limits are consistent with those from XENON1T. We also remark that the constraints quoted in refs. [35] and [34] are upper bounds at $90\%$ CL, while here $95\%$ CL limits are presented.

The evaluation of the constraints on the DM-nucleon cross section is strongly model-dependent, since the mediator properties determine the final results. However, we point out that the model considered in this work is the same investigated by other authors in their recent works [32, 4, 9, 33]. In particular, in Ref. [32], the authors evaluated the limits on $\sigma_{SD}$ in a scenario with the mediator decaying into the $b\bar{b}$, $\tau^{+}\tau^{-}$ and $\gamma\gamma$ channels with a decay length $L=R_{\odot}$, combining Fermi and HAWC 3-years observations. The same scenario was investigated also in Ref. [33] for the $\mu^{+}\mu^{-}$ decay channel. The final limits quoted in these works are those obtained with the lowest allowed mediator mass, since, as already discussed, this condition provides the strongest limits.

Finally, we point out here that the limits presented in this work have been evaluated under the assumption of equilibrium between the DM capture and annihilation processes in the solar environment and assuming $100\%$ branching ratios for the processes $\chi\chi\rightarrow\phi\phi$ and $\phi\rightarrow XX$.

The results obtained in this work show the potentiality of solar gamma rays as a probe for indirect DM detection, since the limits obtained with this analysis are comparable or even stronger than those currently quoted in the literature. We also remark that searches from solar DM could benefit from an accurate model of the steady solar emission, that could improve the description of the background.

The Fermi LAT Collaboration acknowledges generous ongoing support from a number of agencies and institutes that have supported both the development and the operation of the LAT as well as scientific data analysis. These include the National Aeronautics and Space Administration and the Department of Energy in the United States, the Commissariat à l’Energie Atomique and the Centre National de la Recherche Scientifique / Institut National de Physique Nucléaire et de Physique des Particules in France, the Agenzia Spaziale Italiana and the Istituto Nazionale di Fisica Nucleare in Italy, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), High Energy Accelerator Research Organization (KEK) and Japan Aerospace Exploration Agency (JAXA) in Japan, and the K. A. Wallenberg Foundation, the Swedish Research Council and the Swedish National Space Board in Sweden.

Additional support for science analysis during the operations phase is gratefully acknowledged from the Istituto Nazionale di Astrofisica in Italy and the Centre National d’Études Spatiales in France. This work performed in part under DOE Contract DE-AC02-76SF00515.