Search for Dark Matter Annihilation Signals from Unidentified Fermi-LAT Objects with H.E.S.S.

The presence of dark matter (DM) is suggested by a wealth of astrophysical and cosmological measurements; however, its underlying nature is yet unknown. Among the most promising candidates are weakly interacting massive particles (WIMPs): particles thermally produced in the early universe with mass and coupling strength at the electroweak scale predict a present relic density (Steigman et al. 2012) consistent with that observed today (Adam et al. 2016). WIMP self-annihilation would produce Standard Model particles including gamma-rays, which for a long time have been recognized as a prime messenger to indirectly detect DM annihilation or decay. Gamma-rays are not deflected by magnetic fields and therefore point back to their sources. Among the most promising DM targets observed by ground-based imaging atmospheric Cerenkov telescopes (IACTs) such as H.E.S.S. are the Galactic center (Abdallah et al. 2016, 2018a) and nearby dwarf galaxies (Aharonian et al. 2008a; Abramowski et al. 2011, 2014; Abdalla et al. 2018b).

Other compelling and complementary DM targets for IACTs are DM subhalos populating galactic halos (see, e.g., Kamionkowski et al. 2010). The observed universe today is believed to have formed hierarchically with the smallest structures first. DM particles first collapse into gravitationally bound systems that later merge to form the first subhalos, which subsequently form more massive ones. The merging history leads to DM halos massive enough to retain gas and trigger star formation and give rise to the galaxies we observe today. However, most of the subhalos remain completely dark. Assuming that DM is made of WIMPs, they could shine in gamma-rays. The annihilation process in subhalos could be frequent enough to be detectable by IACTs provided that the WIMPs are sufficiently massive. If DM subhalos are made of WIMPs, 10⁻⁴ to 10¹⁰ M_⊙ mass subhalos are expected to lie in DM halos of Milky Way– (MW-) sized galaxies (Diemand et al. 2008; Springel et al. 2008b), with the most massive of them (≳10⁸ M_⊙) hosting dwarf galaxies.

DM subhalos are predicted to be compact and concentrated and are not expected to harbor conventional astrophysical high-energy emitters. Provided they are close and/or massive enough, DM annihilations in these objects could produce gamma-ray fluxes detectable with satellite and ground-based experiments such as IACTs (Calcaneo-Roldan & Moore 2000; Tasitsiomi & Olinto 2002; Stoehr et al. 2003; Koushiappas et al. 2004). However, their actual location in the galaxy is not known. Their search can be performed using all-sky gamma-ray observations (see, for instance, Diemand et al. 2007) such as with the Large Area Telescope (LAT) instrument onboard the Fermi satellite (see, for instance, Berlin & Hooper 2014) or wide-field surveys carried out with IACTs (see, for instance, Aharonian et al. 2008b; Brun et al. 2011).

All-sky Fermi-LAT observations revealed a population of sources that lack associated signals from observations at other wavelengths (Ajello et al. 2017; Abdollahi et al. 2020). These sources are therefore classified as unidentified Fermi objects (UFOs). If the DM particle mass lies below 100 GeV, some UFOs detected by Fermi-LAT could be potentially described by DM models (Belikov et al. 2012; Zechlin et al. 2012; Bertoni et al. 2015, 2016; Calore et al. 2017; Coronado-Blazquez et al. 2019b). Identifying some of the UFOs as DM subhalos requires, however, higher (≳100 GeV) DM masses given their hard gamma-ray spectra in the energy range of a few tens to 100 GeV. Such objects are therefore excellent targets for IACTs to perform searches for tera-electron-volt DM subhalos. In 2018 and 2019, the H.E.S.S. collaboration carried out an observational campaign for a selection of the most promising UFOs in order to probe their potential tera-electron-volt-mass DM-induced emission.

The paper is organized as follows. Section 2 presents the expected DM-induced gamma-ray signals from Galactic subhalos. Section 3 describes the selection procedure of UFOs as DM subhalo candidates relevant for H.E.S.S. observations along with the Fermi-LAT data analysis of the selected UFOs as DM subhalo candidates. In Section 4, the H.E.S.S. observations and data analysis of the selected UFOs are presented. The constraints from H.E.S.S. observations on DM-induced emission models for the UFOs are derived in Section 5. Sections 6 and 7 are devoted to the discussion of the results obtained in this work and to a general summary, respectively.

2.1. Expected Gamma-Ray Flux from DM Annihilation

The energy-differential gamma-ray flux expected from the self-annihilation of Majorana DM particles of mass m_DM can be expressed as

$$\begin{eqnarray}\begin{array}{rcl}\displaystyle \frac{{\rm{d}}{{\rm{\Phi }}}_{\gamma }}{{\rm{d}}{E}_{\gamma }}({E}_{\gamma },{\rm{\Delta }}{\rm{\Omega }}) & = & \displaystyle \frac{\langle \sigma v\rangle }{8\pi {m}_{\mathrm{DM}}^{2}}\sum _{f}{\mathrm{BR}}_{f}\displaystyle \frac{{\rm{d}}{N}^{f}}{{\rm{d}}{E}_{\gamma }}\,J({\rm{\Delta }}{\rm{\Omega }}),\\ \mathrm{with}\quad J({\rm{\Delta }}{\rm{\Omega }}) & = & {\displaystyle \int }_{{\rm{\Delta }}{\rm{\Omega }}}{\displaystyle \int }_{\mathrm{los}}{\rho }^{2}(s(r,\theta )){ds}\,d{\rm{\Omega }}.\end{array}\end{eqnarray} \tag{ 1 }$$

where 〈σ v〉 is the thermally averaged velocity-weighted annihilation cross section and ${\sum }_{f}{\mathrm{BR}}_{f}{{dN}}^{f}/{{dE}}_{\gamma }$ is the sum of the annihilation spectra dN^f /dE_γ per annihilation in the final states f with associated branching ratios BR_f . The expected DM annihilation signal consists of a continuum spectrum of gamma-rays extending up to the DM mass and possibly a line-like feature close to the DM mass. The former contribution arises from the hadronization and/or decay of quarks, heavy leptons, and gauge bosons involved in the annihilation process. The latter comes from the direct annihilation into γ X with X =γ, h, Z or a non–Standard Model neutral particle, providing a spectral line at an energy ${E}_{\gamma }={m}_{\mathrm{DM}}[1\,-{({m}_{X}/2{m}_{\mathrm{DM}})}^{2}]$. When the DM particles self-annihilate into charged particles, gamma-rays are produced via processes involving virtual internal bremsstrahlung and final state radiation. These processes provide an additional bump-like feature that peaks at an energy close to the DM mass.

The term J(ΔΩ), hereafter referred to the J factor, corresponds to the integration of the square of the DM density over the line of sight s and solid angle ΔΩ. As opposed to objects with measured stellar dynamics like dwarf galaxies, UFOs have unknown distances to Earth and their J factors cannot be derived from stellar kinematics.

2.2. Expected Subhalo J-factor Distribution in the MW

Cosmological N-body simulations (see, for instance, Diemand et al. 2008; Springel et al. 2008b) predict MW-sized galaxy halos to harbor as of yet unmerged DM substructures called galactic subhalos. These simulations make robust predictions on the slope and normalization of the subhalo mass function defined as $\mathrm{dln}N/\mathrm{dln}M\propto {M}^{-{\alpha }_{{\rm{m}}}}$, with a slope α _m ≃ 1.9 for MW-like galaxies (see, for instance, Diemand et al. 2008; Springel et al. 2008a; Gao et al. 2012; Fiacconi et al. 2017). Galactic subhalos are not expected to host conventional high-energy astrophysical sources, and a large number density of subhalos in MW-like galaxies with high DM concentrations is expected. Only close by and massive enough subhalos are expected to provide detectable gamma-ray signals.

The cosmological simulations provide the abundance of the resolved subhalos, their radial distribution, and structural properties such as their mass and concentration and suggest a DM distribution in the subhalos following a cuspy density profile that can be well described by Navarro–Frenk–White (NFW) (Navarro et al. 1997) or Einasto (Springel et al. 2008c) parameterizations. The limited spatial resolution of the current simulations can only probe the slope of the radial DM density distribution in a subhalo for the most massive subhalos. Assuming a parameterization of the DM density distribution in subhalos, the distribution of J factors of the galactic subhalo population, dN/dJ, can be derived. The cumulative J-factor distribution, N(J) ≡ N(≥ J), is defined as the number of subhalos with a J factor higher than or equal to specified.

In order to derive the J-factor distribution of DM subhalos in the MW, the CLUMPY code v3.0.0 (Charbonnier et al. 2012; Bonnivard et al. 2016; Hütten et al. 2019a) is used. One thousand simulations of an MW-like galaxy are performed with a smooth NFW (Navarro et al. 1997) DM main halo profile with the parameters corresponding to the best-fit NFW parameters presented in a recent study of DM distribution in the MW (Cautun et al. 2020). For each simulation, the subhalo parameters were chosen similar to the ones used in Hütten et al. (2016) for the "HIGH" model. The power-law slope of the subhalo mass function is chosen to be α _m = 1.9 (Diemand et al. 2008); the number of objects between 10⁸ and 10¹⁰ M_⊙ is taken as N_calib = 300 following Springel et al. (2008b); and the subhalo mass-concentration relation follows the distance-dependent prescription of Moliné et al. (2017). From each simulation, the Galactic coordinates of all subhalos and their J factors integrated in circular regions with 0.1° radius around the centers of gravity of the subhalos are derived.

The cumulative J-factor distribution N(≥ J) is shown in the upper panel of Figure 1 for subhalos located at Galactic latitudes ∣b∣ > 5°. The solid red curve shows the averaged distribution computed from all of the realizations, and the shaded region shows the formal 1σ statistical dispersion calculated over all simulated MW-like galaxies. In the lower panel of Figure 1, the red dashed curve illustrates the probability of finding in any simulation at least one subhalo with a J factor higher than specified. The blue dotted curve corresponds to the probability of finding three or more subhalos. The horizontal black dashed line gives the 5% probability: the probability on average of finding one or more subhalos with J factor J ≥ 3 × 10²⁰ GeV² cm⁻⁵ (three or more subhalos with J factors J ≥ 1 × 10²⁰ GeV² cm⁻⁵) is only about 5%.

Zoom In Zoom Out Reset image size

A hypothetical gamma-ray emission from a DM subhalo would show up in the all-sky gamma-ray surveys (Kamionkowski et al. 2010) as unidentified objects, that is, detected with the Fermi satellite but without counterparts at any other wavelengths. A smoking-gun signature for DM detection is a very distinct energy cutoff close to the DM particle mass, assuming the two-body annihilation process taking place almost at rest. For sufficiently large DM particle masses, that is, above a few hundred giga-electron-volts, the energy cutoff could be too high in energy to be measurable by the LAT instrument onboard the Fermi satellite within a reasonable observation time. The combination of Fermi-LAT and IACT observations is therefore mandatory for DM subhalo searches with unidentified sources detected by the Fermi-LAT.

3.1. Selection Procedure for H.E.S.S. Observations

3.2. Fermi-LAT Data Analysis on the Selected Objects

3.3. DM Models for the Selected Sources

Figure 2 shows DM annihilation models together with Fermi-LAT flux measurements. Model predictions for DM masses of 1 and 10 TeV, respectively, are plotted separately for the W⁺ W⁻ and τ⁺ τ⁻ annihilation channels. Some DM models are able to qualitatively describe the observed gamma-ray flux from the selected UFOs. For instance, the predictions shown for m_DM = 1 TeV describe well the Fermi-LAT data for 3FHL J0929.2-4110 in the W⁺ W⁻ and τ⁺ τ⁻ annihilation channels. In contrast, for 3FHL J2030.2-5037, the predictions for m_DM =1 TeV in the τ⁺ τ⁻ annihilation channel provide a poor description of the emission, while a better one is obtained in the W⁺ W⁻ one.

Zoom In Zoom Out Reset image size

In order to assess quantitatively the above statements on viable DM-induced emission models, the spectrum of each UFO is explicitly modeled with a DM annihilation–induced spectral template. ⁴³ For a given m_DM and annihilation channel, the model is characterized only by the overall normalization of the spectra given by 〈σ v〉J. To identify the range of viable parameters for DM annihilation, a scan over a large range of 〈σ v〉J is performed.

A TS is defined as a difference between best-fit log-likelihood functions for models with no DM emission (${{ \mathcal L }}_{0}$, "background-only" hypothesis) and the model (${ \mathcal L }$), which includes the UFO source described by the corresponding parameter 〈σ v〉J: ${TS}=-2\,\mathrm{log}({ \mathcal L }/{{ \mathcal L }}_{0})$ (Mattox et al. 1996). ⁴⁴ Negative TS values correspond to the detection of the source, that is, adding a source with a corresponding parameter improves the fit compared with the background-only hypothesis.

The left panels of Figure 3 illustrate the results for a single UFO 3FHL J0929.2-4110 for the W⁺ W⁻ (top) and τ⁺ τ⁻ (bottom) annihilation channels, respectively. For each 〈σ v〉J, the color scale shows the TS values. Assuming that the TS follows a χ² distribution, a TS equal to −9, (resp. −25) would correspond to a 3σ (resp. 5σ) detection for 1 degree of freedom. The dashed cyan and orange lines show the detection region that corresponds to the improvement of TS by −9 and −25, respectively. The right panels in Figure 3 present the results for the analysis of the combined data sets of the three selected UFOs (3FHL J0929.2-4110, 3FHL J1915.2-1323, and 3FHL J2030.2-5037) obtained through the combination of the log-likelihood profiles from individual objects, for the W⁺ W⁻ (top) and τ⁺ τ⁻ (bottom) annihilation channels, respectively.

Zoom In Zoom Out Reset image size

4.1. Data Analysis

H.E.S.S. is an array of five IACTs located in the Khomas Highland in Namibia, at an altitude of 1800 m. The array is composed of four 12 m diameter telescopes (CT1-4) and a fifth 28 m diameter telescope (CT5) at the middle of the array. The observations presented here were performed in 2018 and 2019 with the full five-telescope array for the selection of unidentified Fermi-LAT objects presented in Table 2. They were carried out in the wobble mode, where the telescope pointing direction is offset from the nominal target position by an angle between 05 and 07. The observations for the data analysis are selected according to the standard run selection criteria (Aharonian et al. 2006). After the calibration of raw shower images recorded in the camera, the reconstruction of the direction and energy of the gamma-ray events is performed with a template-fitting technique (de Naurois & Rolland 2009) in which the recorded images are compared to precalculated showers computed from a semianalytical model. An energy resolution of 10% and an angular resolution of 006 at 68% containment radius for gamma-ray energies above 200 GeV are achieved. The results described here have been cross-checked with an independent calibration and analysis chain yielding compatible results (Parsons & Hinton 2014).

The selected UFOs are assumed to be pointlike sources according to the PSF of Fermi-LAT, which reaches ∼01 above 100 GeV. Given the H.E.S.S. PSF, the ROI, hereafter referred to as the ON source region, is therefore defined as for pointlike-emission searches for H.E.S.S., and the ROI is taken as a disk of 012 radius. The residual background is measured in OFF regions according to the MultipleOff technique (Aharonian et al. 2006). For each telescope pointing position, the OFF regions are defined at the same distance from the pointing position as for the ON region, which leads to identical acceptances in the ON and OFF regions. An excluded region defined as a disk of radius equal to twice the ON region radius is used in order to avoid any leakage of the searched signal into the OFF regions. The α parameter is defined as the ratio between the solid angle size of the OFF and ON regions by α = ΔΩ_OFF/ΔΩ_ON. The excess significance in the ROI is computed following the statistical approach of Li & Ma (1983). Table 3 summarizes for each UFO the live time, the mean zenith angle of the observations, the ON and OFF counts, the α parameter averaged over all the observations, and the excess significance in the ROI. No significant gamma-ray excess is found, neither in the ON source region nor anywhere else in the field of view.

Table 3. H.E.S.S. Data Analysis Results for Each UFO

Name	Live time	Mean zenith angle	NON	NOFF	$\bar{\alpha }$	Significance
	(hours)	(degrees)	(counts)	(counts)		(σ)
3FHL J0929.2-4110	27.4	29.0	424	5884	13.9	0.1
3FHL J1915.2-1323	3.6	19.4	87	1181	13.9	0.2
3FHL J2030.2-5037	9.8	31.3	160	2192	13.9	0.1
3FHL J2104.5+2117	6.8	46.7	73	853	13.9	1.1

Note. The second and third columns give the live time and mean zenith angle of the H.E.S.S. observations, respectively. Count numbers measured in the ON and OFF regions are given in the fourth and fifth columns, respectively, with the α parameter averaged over all observations, $\bar{\alpha }$, given in the sixth column. The seventh column provides the measured excess significance between the ON and OFF counts.

Download table as:  ASCIITypeset image

Since no significant excess is found from any of the selected UFOs, energy-differential observed flux upper limits are computed at 95% C.L. assuming the best-fit power-law spectral index derived from the Fermi-LAT data analysis. Expectations for the energy-differential flux upper limits are computed from 100 independent Poisson realizations of the measured OFF counts for the ON and OFF regions. From the realizations, mean and standard deviation values are extracted and used to compute the 68% and 95% containment bands. The differential flux upper limits are shown in Figure 2 with an energy binning of 0.5 dex. The expectations and the containment bands for the upper limits shown in the figure are computed including 25% systematic uncertainty. Assuming a power-law spectral index of 2 would change the differential flux upper limits by less than 6%.

4.2. Upper Limit Computation for the DM Emission

A binned Poisson maximum likelihood analysis is performed to search for spectral features expected from DM annihilation signals with respect to the background. For each UFO, the H.E.S.S. energy range is divided into 62 logarithmically spaced bins from 100 GeV up to 70 TeV. For a given DM mass and annihilation channel, the Poisson likelihood function in the energy bin i can be written as

$$\begin{eqnarray}&&\begin{array}{l}{{ \mathcal L }}_{{\rm{i}}}({N}^{{\rm{S}}},{N}^{{\rm{B}}}| {N}_{\mathrm{ON}},{N}_{\mathrm{OFF}},\alpha )=\displaystyle \frac{{\left({N}_{{\rm{i}}}^{{\rm{S}}}+{N}_{{\rm{i}}}^{{\rm{B}}}\right)}^{{N}_{\mathrm{ON},\ {\rm{i}}}}}{{N}_{\mathrm{ON},\ {\rm{i}}}!}{e}^{-({N}_{{\rm{i}}}^{{\rm{S}}}+{N}_{{\rm{i}}}^{{\rm{B}}})}\\ \quad \times \displaystyle \frac{{\left({N}_{{\rm{i}}}^{{{\rm{S}}}^{{\prime} }}+\alpha {N}_{{\rm{i}}}^{{\rm{B}}}\right)}^{{N}_{\mathrm{OFF},\ {\rm{i}}}}}{{N}_{\mathrm{OFF},\ {\rm{i}}}!}{e}^{-({N}_{{\rm{i}}}^{{{\rm{S}}}^{{\prime} }}+\alpha {N}_{{\rm{i}}}^{{\rm{B}}})}.\end{array}\end{eqnarray} \tag{ 2 }$$

The variables N_ON,i and N_OFF,i stand for the number of measured events in the ON and OFF regions, respectively; ${N}_{{\rm{i}}}^{{\rm{B}}}$ is the expected number of background events in the ON region; and ${N}_{{\rm{i}}}^{{\rm{S}}}$ and ${N}_{{\rm{i}}}^{{{\rm{S}}}^{{\prime} }}$ are the expected number of DM signal events in the ON and OFF regions, respectively. They are computed by folding the expected theoretical DM flux given in Equation (1) with the energy-dependent acceptance and energy resolution of H.E.S.S. for the considered data set. The term ${dN}/{{dE}}_{\gamma }^{f}$ in Equation (1) is extracted from Cirelli et al. (2011) for each assumed DM mass and annihilation channel. The energy resolution of H.E.S.S. is represented by a Gaussian function with a width of σ_E/E = 10% above 200 GeV. The UFOs are also pointlike sources for H.E.S.S.; therefore no leakage is expected in the background region and ${N}_{{\rm{i}}}^{{{\rm{S}}}^{{\prime} }}$ is taken to ${N}_{{\rm{i}}}^{{{\rm{S}}}^{{\prime} }}\equiv$ 0. The likelihood function for a given object ${ \mathcal L }$ is defined as ${ \mathcal L }={\prod }_{{\rm{i}}}{{ \mathcal L }}_{{\rm{i}}}$.

Since no significant excess is found in any of the selected UFOs by H.E.S.S., upper limits can be derived assuming UFOs are DM-induced gamma-ray emitters from a likelihood ratio TS given by

$$\begin{eqnarray}&&{TS}=-2\,\mathrm{log}\displaystyle \frac{{ \mathcal L }({N}^{{\rm{S}}}(\langle \sigma v\rangle J),\widehat{\widehat{{N}^{{\rm{B}}}}}(\langle \sigma v\rangle J)| {N}_{\mathrm{ON}},{N}_{\mathrm{OFF}},\alpha )}{{ \mathcal L }(\widehat{{N}^{{\rm{S}}}}(\langle \sigma v\rangle J),\widehat{{N}^{{\rm{B}}}}| {N}_{\mathrm{ON}},{N}_{\mathrm{OFF}},\alpha )}.\end{eqnarray} \tag{ 3 }$$

The variable $\widehat{\widehat{{N}_{{\rm{i}}}^{{\rm{B}}}}}$ is obtained through a conditional maximization, achieved by solving $d{ \mathcal L }/{{dN}}_{{\rm{i}}}^{{\rm{B}}}=0$; $\widehat{{N}_{{\rm{i}}}^{{\rm{S}}}}$ and $\widehat{{N}_{{\rm{i}}}^{{\rm{B}}}}$ are computed using an unconditional maximization. Following the procedure defined in Cowan et al. (2011), upper limits are computed assuming a positive signal, that is, $\widehat{{N}_{{\rm{i}}}^{{\rm{S}}}}\gt 0$. The 〈σ v〉J value for which the TS value is equal to 2.71 provides the one-sided 95% C.L. upper limit on the quantity 〈σ v〉J.

The hypothesis that all UFOs are indeed DM subhalos, but too faint to be detected as such in the tera-electron-volt energy range with the given exposure, can be tested when the individual data sets are combined. If no significant overall excess is found in the combined data set, combined upper limits on 〈σ v〉J can be derived versus the DM mass assuming J to be an averaged J factor. In Equation (3), the likelihood function is replaced by the combined likelihood expressed as ${{ \mathcal L }}_{\mathrm{comb}}={\prod }_{j=1}^{{{\rm{N}}}_{{targets}}}{{ \mathcal L }}_{{\rm{j}}}$, where ${{ \mathcal L }}_{{\rm{j}}}$ is the likelihood of the target j.

No significant excess is measured in any of the H.E.S.S. data sets of the selected UFOs. The 95% C.L. upper limits on 〈σ v〉J are derived versus the DM mass using Equation (3). Figure 4 shows for each UFO the upper limits as a function of the DM mass for the W⁺ W⁻ and τ⁺ τ⁻ annihilation channels, respectively. In most of the DM mass range , the strongest constraints are reached for 3FHL J0929.2-4110 observations. For a 1 TeV DM mass, the constraints 〈σ v〉J = 5.5 × 10⁻⁵ GeV² cm⁻² s⁻¹ and 1.9 × 10⁻⁵ GeV² cm⁻² s⁻¹ in the W⁺ W⁻ and τ⁺ τ⁻ annihilation channels, respectively, for 3FHL J0929.2-4110.

Zoom In Zoom Out Reset image size

The combined analysis of the four H.E.S.S. data sets does not show any significant excess. Therefore, the UFO data sets are combined and upper limits on 〈σ v〉J are computed. Given its possible association with an active galactic nucleus, the source 3FHL J2104.5+2117 is removed to provide conservative combined upper limits. The right panel of Figure 3 shows the combined 95% C.L. upper limits on 〈σ v〉J as a function of the DM mass for the W⁺ W⁻ and τ⁺ τ⁻ annihilation channels, respectively. The analysis of the combined data sets allows for an improvement of about 10% and 20% for 1 TeV DM mass in the W⁺ W⁻ and τ⁺ τ⁻ annihilation channels, respectively, with respect to the most constraining upper limit from the individual UFO data sets. The combined limits reach 3.7 × 10⁻⁵ and 8.1 × 10⁻⁶ GeV² cm⁻² s⁻¹ in the W⁺ W⁻ and τ⁺ τ⁻ channels, respectively, for a 1 TeV DM mass.

In order to derive the J-factor values required to explain the UFO emission in terms of DM models, the value of the annihilation cross section expected for thermal WIMPs (〈σ v〉_th ≃ 3 × 10⁻²⁶ cm³ s⁻¹) is used (Steigman et al. 2012). The J-factor upper limits for the DM models of the UFOs as a function of the DM mass are given at 95% C.L. in Figure 5. For a 1 TeV DM mass in the W⁺ W⁻ channel, the J-factor values are constrained to be between 2.4 × 10²⁰ and 1.3 × 10²¹ GeV² cm⁻⁵ for DM models with TS ≤ −25 (which corresponds to ≥5σ confidence interval assuming TS follows χ² distribution). For a DM mass of 10 TeV in the W⁺ W⁻ channel, all the J-factor values for DM models with TS ≤−25 are ruled out at 95% C.L. by the H.E.S.S. constraints. In the τ⁺ τ⁻ channels, the H.E.S.S. constraints are even stronger. For 300 GeV DM mass, the allowed J-factor values are between 1.4 × 10²⁰ and 5.9 ×10²⁰ GeV² cm⁻⁵ for TS ≤−25 DM models. The H.E.S.S. upper limits restrict the J factors to lie in the range 6.1 × 10¹⁹–2.0 × 10²¹ GeV²cm⁻⁵ and the masses to lie between 0.2 and 6 TeV in the W⁺W⁻ channel. For the τ⁺ τ⁻ channel, the J factors lie in the range 7.0 × 10¹⁹–7.1 × 10²⁰ GeV²cm⁻⁵ and the masses lie between 0.2 and 0.5 TeV.

Zoom In Zoom Out Reset image size

Using predictions of N-body cosmological simulations, the number of subhalos with a J factor higher than a given value for an MW-like galaxy can be extracted as displayed in Figure 1. The probability of having at least three subhalos with a J factor higher than 10²⁰ GeV² cm⁻⁵ is below 5% (Figure 1, blue dotted line). According to this prediction, the interpretation of the UFO emissions in terms of DM particle annihilations in Galactic DM subhalos can be further constrained from Figure 5 to m_DM ≲ 1 TeV for W⁺ W⁻ and m_DM ≲ 0.3 TeV for τ⁺ τ⁻ channels.

A substantial number of UFOs may emit gamma-rays from DM annihilations in subhalos. However, some of them could be active galactic nuclei or other types of galaxies that so far lack detection at other wavelengths. Alternative astrophysical interpretations of UFOs as pulsars or low-luminosity globular clusters hosting millisecond pulsars (Mirabal et al. 2016) may be less plausible since typical gamma-ray spectra for these types of objects are characterized by an energy cutoff at energies of a few giga-electron-volts.

The cumulative J-factor distribution is in very good agreement with the results of Hütten et al. (2016) for the "HIGH" model intended to predict the highest possible number of subhalos in a typical MW-like galaxy. The real number of DM subhalos can be an order of magnitude smaller as shown for the predictions in the "LOW" model of Hütten et al. (2016). The choice of the number of subhalos of masses between 10⁸ and 10¹⁰ M_⊙ of N_calib = 300 is motivated by the output of DM-only simulations Springel et al. (2008b). Baryon feedback can significantly reduce this value, up to a factor of two (Mollitor et al. 2015; Sawala et al. 2016). This would make the highest J-factor values even more unlikely. As discussed in Coronado-Blázquez et al. (2019a), subhalos with the highest J factors should appear as extended sources for Fermi-LAT given its PSF of about 01 above 10 GeV. However, even these brightest DM subhalos would produce faint gamma-ray sources whose spatial extension would be challenging to measure for Fermi-LAT. On the simulation front, further work is likely needed to use predictions for subhalo angular sizes in MW-like galaxies to definitely rule out pointlike UFOs as potential DM subhalos.

The DM density distribution of the Galactic halo is assumed here to follow an NFW parameterization. Incorporating hydrodynamics and baryon feedback in cosmological simulations tends to soften the inner cusp of the DM profiles in MW-like galaxies, leading to a flattening of order 1 kpc (see, for instance, Chan et al. 2015). However, these predictions for the expected DM distribution have large uncertainties due to the effects of baryonic physics. The resolution limit of the simulations at sufficiently small distances also becomes pertinent. Alternative Galactic mass models can be used to describe subhalo parameters for MW-like galaxies (Catena & Ullio 2010; McMillan 2011; Stref & Lavalle 2017; McMillan 2017). Using different Galactic mass models, the subhalo luminosity functions derived in Stref & Lavalle (2017) provide compatible results. Considering a core profile would make the high DM mass exclusion of the DM models for the UFO emission even stronger. DM subhalos could also have cored profiles with lower DM concentration, which would make them more prone to tidal disruption. This would therefore both decrease the normalization of the J-factor distribution function and shift it to lower values. The J-factor distribution prediction is based on DM-only Via Lactea-II simulations using WMAP cosmology. The simulations based on the most recent cosmological results from the Planck mission with added baryonic physics could somewhat change the predicted properties of the MW subhalos. The effect of baryon feedback and tidal effects induced by both DM and baryons is likely to alter the DM concentration of subhalos (Despali & Vegetti 2017; Stref & Lavalle 2017). Details of the modeling of tidal disruption of Galactic DM subhalos on the brightest subhalo can be found in Di Mauro et al. (2020). Therefore, the cumulative J-factor distribution would be shifted to lower values. In this case, the DM-induced interpretation of the UFO emission would be even more constrained given that the probability of finding high J-factor values would be even smaller. The presence of baryons affects both the amplitude and the slope of the DM halo and subhalo mass functions. It can reduce the slope by a few percent in the 10⁶–10⁹ M_⊙ subhalo mass range (Benson 2020), which can also alter the rate of subhalos with large J factors. No cut is applied to the maximal value of the subhalo mass when computing the cumulative J-factor distribution. While subhalos with masses of about 10⁷ M_⊙ or higher may be able to host star formation and actually be dwarf galaxies, simulations including hydrodynamics and feedback physics in addition to the gravitational effects for the expected DM distribution in both the main halo and its subhalos such as in Zhu et al. (2016) show that a significant fraction of subhalos with masses of about 10⁹ M_⊙ is found to host no stars. For subhalo masses larger than 10⁷–10⁸ M_⊙, a noticeable fraction of them can start triggering star formation and indeed form a faint dwarf galaxy. Arguably these masses critically depend on the baryonic physics implementation in the simulations and its associated feedback. Considering subhalos with masses lower than 10⁷ M_⊙ implies the probability of having at least one subhalo with J ≥ 3 × 10²⁰ GeV² cm⁻⁵ to be about 0.3%, therefore a factor of about 16 lower that in the case without mass cut.

The interpretation of UFOs as DM subhalos of tera-electron-volt-mass scale thermal WIMPs requires J factors in excess of a few 10²⁰ GeV² cm⁻⁵. Such J-factor values are only occasionally obtained in N-body simulations of MW-type galaxies. The highest subhalo J-factor values are usually subject to a large statistical variance. The precise value of the brightest subhalos can be subject to a large uncertainty, implying a factor of 10 uncertainty on the J-factor value for J ≳ 10²⁰ GeV² cm⁻⁵ in the "HIGH" model (Hütten et al. 2019b). ⁴⁵ Previous studies using Fermi-LAT data only searched for UFOs as promising DM subhalo candidates for DM masses below 100 GeV (see, for instance, Bertoni et al. 2015, 2016; Coronado-Blazquez et al. 2019b). For instance, in the study of Coronado-Blázquez et al. (2019a), masses of a few tens of giga-electron-volts are ruled out for canonical thermal WIMPs. The present analysis searches for UFOs as DM subhalos for the uncharted tera-electron-volt-mass thermal WIMP models.

The maximal value of subhalo J factors obtained in simulations is model dependent and can be increased even in comparison with the optimistic estimate of the J-factor distribution considered here as discussed below. The normalization of the subhalo mass function assumes usually about 10% of the total DM halo content to be in the form of subhalos, with the total DM halo density being normalized to the DM density at the location of the Sun ρ(r_⊙) = ρ_⊙ = 0.39 GeV cm⁻³. The precise value is subject to uncertainties of a factor of 1.5–2 (Read 2014). Varying the input parameters of the simulations in the relevant ranges of interest, such as ρ_⊙ = 0.6 GeV cm⁻³ and the scale radius of the main DM halo r_s = 25 kpc, would increase the highest J-factor values by a factor of two. The predicted cumulative J-factor distribution would therefore probe higher J-factor values. Considering substructures in galactic subhalos (see, for instance, Hiroshima et al. 2018) results in higher expected J-factor values for the Galactic subhalo population with typical increase factors of a few, therefore shifting the cumulative J-factor distribution to higher values.

The abovementioned large systematic uncertainties in the prediction of the J-factor distribution weaken significantly the constraints from cosmological simulations, which makes them comparable to or weaker than the H.E.S.S. constraints in, for example, the τ⁺ τ⁻ channel. This makes the model-independent H.E.S.S. constraints the only relevant ones for robust interpretation of UFOs as Galactic subhalos of annihilating DM.

In this work, a straightforward filtering of the unidentified sources in the 3FHL point-source catalog has been performed using selection cuts to identify the most promising DM subhalo candidates for DM masses above a few hundred giga-electron-volts. The data sets for the four UFOs were collected with the Fermi satellite in a 12 yr observation period. Using H.E.S.S. observations, no significant signal has been detected from any of the selected UFOs. From H.E.S.S. flux upper limits, DM models describing the UFO emissions with high significance are strongly constrained in the tera-electron-volt DM mass range for different annihilation channels. Assuming thermally produced DM particles, the DM models for the UFO emissions require high J-factor values. When the model-dependent predictions from N-body simulations of the MW-like subhalo population are taken into account, the required high J-factor values for the DM models explaining the UFOs as Galactic subhalos are unlikely. This can point toward interpretation of the UFOs as subhalos of relatively light WIMPs with masses m_DM ≲ 0.3 TeV, for which somewhat lower J factors are preferred. However, this could be in tension with constraints on thermal WIMPs from dwarf galaxy observations by Fermi-LAT (Ackermann et al. 2015; Albert et al. 2017).

The support of the Namibian authorities and of the University of Namibia in facilitating the construction and operation of H.E.S.S. is gratefully acknowledged, as is the support by the German Ministry for Education and Research (BMBF), the Max Planck Society, the German Research Foundation (DFG), the Helmholtz Association, the Alexander von Humboldt Foundation, the French Ministry of Higher Education, Research and Innovation, the Centre National de la Recherche Scientifique (CNRS/IN2P3 and CNRS/INSU), the Commissariat à l'énergie atomique et aux énergies alternatives (CEA), the U.K. Science and Technology Facilities Council (STFC), the Knut and Alice Wallenberg Foundation, the National Science Centre, Poland grant No. 2016/22/M/ST9/00382, the South African Department of Science and Technology and National Research Foundation, the University of Namibia, the National Commission on Research, Science & Technology of Namibia (NCRST), the Austrian Federal Ministry of Education, Science and Research and the Austrian Science Fund (FWF), the Australian Research Council (ARC), the Japan Society for the Promotion of Science, and the University of Amsterdam.

We appreciate the excellent work of the technical support staff in Berlin, Zeuthen, Heidelberg, Palaiseau, Paris, Saclay, Tübingen, and Namibia in the construction and operation of the equipment. This work benefitted from services provided by the H.E.S.S. Virtual Organisation, supported by the national resource providers of the EGI Federation.

In this section we present the results for the search of possible sources not accounted for by the considered model based on the 4FGL catalog of the UFOs vicinity. We built TS maps for 5° × 5° region centered at the position of each UFO. These maps illustrate the significance (∼$\sqrt{{TS}}$) above the background model of a test point-like source with a power-law spectrum characterized by a free normalization and slope fixed to −2, in each pixel.

For the background model, where the corresponding UFO source is removed, we consider the same spatial/spectral models as used for the analysis of Fermi-LAT data described in Section 3.2. Such a choice of the background model allows us to check the presence of unaccounted for sources close to UFOs' positions as well as verify the point-like spatial morphology of the emission associated with the UFO.

Figure 6 show the TS maps for all four UFO sources for energies above 10 GeV in Galactic coordinates with a pixel size of 005. No smoothing is applied to the maps. Cyan crosses show the positions of UFO sources. Green circles indicate the position of the nearby 4FGL sources included in the background model. The maps indicate the absence of significant residuals, which could affect the Fermi-LAT analysis of the UFOs.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In Figure 7 we show the contour of TS from Fermi-LAT datasets and the 95% C.L. upper limits on J as a function of the DM mass for 3FHL J1915.2-1323, 3FHL J2030.2-5037 and 3FHL J2104.5+2117. Contours and upper limits for the W⁺W⁻ and τ⁺ τ⁻ annihilation channels are shown in the left and right panels, respectively.