Searching for Light in the Darkness: Bounds on ALP Dark Matter with the optical MUSE-Faint survey

We use MUSE spectroscopic observations of the dwarf spheroidal galaxy Leo T between 470 and 935 nm to search for radiative decays of axion like particles (ALPs). Under the assumption that ALPs constitute the dark matter component of the Leo T halo, we derive bounds on the effective ALP-two-photon coupling. We improve existing limits by more than one order of magnitude in the ALP mass range 2.7-5.3 eV.


Introduction
Axion-like particles (ALPs) are compelling cold dark matter candidates. They are a generalization of the QCD axion, which was originally introduced to solve the strong charge-parity (CP) problem [1,2,3,4]. Numerous experimental strategies have been envisaged to detect ALPs (see Ref. [5] for a review), many of them exploiting their coupling to photons given by the operator L = − 1 4 g aγγ a F µνFµν , where a is the ALP field, F µν is the electromagnetic field strength,F µν its dual, and g aγγ the coupling constant. In astrophysical environments, this operator leads to photon signals from the radiative decay of ALP dark matter [6,7,8,9,10,11,12] and the conversion of ALPs into photons in the presence of magnetic fields [13,14,15,16,17,18,19]. For ALP masses in the eV range, the monochromatic line emission from ALP decays falls in the optical and near infrared bands. Upper limits on this signal have been obtained from observations of galaxy clusters in Ref. [20], see also Refs. [21,22] for previous analysis. Interestingly, ALPs of such masses have been proposed to explain an excess of the measured cosmic near infrared background and its angular anisotropies [23,24]. Recently the XENON1T experiment has reported an excess of electron recoil events over expected backgrounds with a significance around 3.5σ [25]. This anomaly can be interpreted in terms of ALPs with masses < 100 eV produced in the interior of the Sun and then detected by the Xenon1T experiment through the axioelectric effect, exploiting the couplings of ALPs with electrons 3 . While this result calls for further examination, it certainly spurs interest in ALPs.
Dwarf spheroidal galaxies are ideal targets to search for indirect dark matter signals, given their proximity and their high dark matter content. In this work, assuming that all dark matter is in the form of ALPs, we constrain its coupling to photons by searching for the optical line from ALP decays in the Leo T dwarf spheroidal galaxy. We exploit data from a Guaranteed Time Observing (GTO) programme targeting ultra-faint dwarf galaxies using the Multi Unit Spectroscopic Explorer (MUSE) at the VLT [26]. These spectroscopic data cover the range 470-935 nm with a medium spectra resolution (R = λ/∆λ > 10 3 ) and excellent sensitivity, covering the area around the centre of Leo T, extending up to approximately its half-light radius.
Previous astrophysical studies [20] placed an upper limit on g αγγ of ≈ 5 × 10 −12 GeV −1 for ALP masses between 4.5 and 5.5 eV. In the following we will demonstrate that the MUSE data on Leo T improves current constraints by more of an order of magnitude for ALP masses between 2.7 and 5.3 eV.
The paper is structured as follows. The data from MUSE observations are presented in Sec. 2. The computation of the ALP signal is detailed in Sec. 3. The statistical analysis and results are discussed in Sec. 4. We conclude in Sec. 5.

Observations and data reduction
The central region of Leo T was observed as part of MUSE-Faint [27], a GTO survey of ultra-faint dwarf galaxies (PI Brinchmann) with MUSEww, a large-field medium-resolution Integral Field Spectrograph on the Very Large Telescope (VLT). The data described here use 15 exposures of 900 seconds adding up to 3.75 hours of exposure time 4 , see also Vaz et al (in preparation). The data were taken in the Wide Field Mode with adaptive optics (WFM-AO), which provides a 1 × 1 arcmin 2 field of view with a spatial sampling of 0.2 arcsec pixel −1 and with a spatial resolution of 0.61 arcsec (full width-half maximum) at a wavelength of 7000Å. The data covers a wavelength range of 4700 − 9350Å with a wavelength sampling of 1.25Å. In order to avoid the light of the sodium laser of the adaptive optics system, a blocking filter removes light in the wavelength range 5820-5970Å (2.13-2.08 eV) which shows up as a gap in the constraints below.
Standard data reduction using the MUSE Data Reduction Software (DRS; version 2.8.1; Ref. [28]) was used, with the most salient points for this paper being that the data were flux calibrated using flux standards observed during the night, atmospheric emission lines, caused by Raman scattering of the laser light of the adaptive optics system, were removed, and we performed a subtraction of emission lines from the night sky. These latter have well-known wavelengths and lead to increased noise at these wavelengths.
Since a key aim of this paper is to place limits on the presence of emission lines across the field of view of Leo T, it is important that we have a good estimate of the noise in the data cube. The standard noise estimate is based on the uncertainty propagated through the DRS, but it is known [e.g. 29] that the MUSE DRS underestimates the uncertainties in the final data-cube since covariance terms are neglected. To address this point we re-estimated the pixel-to-pixel variance directly from each individual datacube following the approach described in Ref. [29], using SExtractor [30] to define a mask image. Finally, all single exposure datacubes were combined using MPDAF [31], creating the datacube that has been used in the following part of the analysis. We discuss the impact of different noise estimates have on the results further below.
The data contain a large number of stellar sources within the field of view, many from Leo T (to be discussed in Vaz et al in prep), but also some likely foreground stars from the Milky way. Moreover also some galaxies are present. To reduce the impact of these sources on the final results we mask the brightest ones. To construct the mask we took two steps: first we created a white-light image by summing over the wavelength axis, which increases the signal-to-noise in each source ensuring that even sources that are undetected in a single wavelength layer are masked out. We then ran SExtractor on the white-light image, with a detection threshold of 3σ (the results are not very sensitive to this particular choice, see also below). The resulting segmentation map is used to mask sources, while no attempt to subtract them has been made. Thus we will mainly consider data where no sources are detected in the white-light image.
Let us point out two caveats to this analysis. (a) Emission from the night sky cannot be perfectly subtracted, which has two consequences: one is that there are residuals, both positive and negative, around the strongest sky emission lines and the other is that the mean flux where no objects are found may differ from zero after sky subtraction. (b) There are a number of galaxies in the data cube that are mostly seen as line emitters and not identified as sources with the above procedure. Note that faint stars have a smooth spectrum with energy interspersed with absorption features and should not confound our results.
For (a) we could apply a sky subtraction correction code [e.g. ZAP 32], but this procedure might add subtle effects in the data that are undesirable for our analysis. Instead we include an arbitrary spatially flat term in each wavelength layer (and set it through the fit) to account for any zero-point offset. To address b) we used the ORIGIN emission line detection software [33] to blindly detect emission line sources which we then mask out. These emission line sources turn out to be rare and their masking does not change the results significantly.
In the left panel of Fig. 1, we show the aforementioned white-light image of Leo T with the 3σ segmentation mask overlaid. In the right panel we report with a blue line the flux density averaged over all the unmasked pixels of our region of interest and in the different wavelength channels of MUSE observations.

ALP signal
The flux density at wavelength λ produced by decays of ALPs inside Leo T from a given direction identified by θ can be computed as: where we have assumed spherical symmetry. Deviations from spherical symmetry are constrained to small values for LeoT [34], and, in any case, their possible presence would have a minor impact on the final bounds. The decay rate Γ a depends on the the ALP mass m a and the effective ALP-two-photon coupling g aγγ . In natural units, it reads Γ a = m 3 a g 2 aγγ /(64π). The dark matter spatial density distribution ρ a (r) is written as a function of the distance r from the center of the dwarf. This distance can be expressed in terms of the coordinate along the line of sight , the angle θ and the distance D of Leo T from us by means of r 2 = 2 + D 2 − 2 D cos θ . The angle θ combines the direction of observation θ, i.e. the angular off-set with respect to the Leo T centre, and the position Ω inside the observing angular beam. The latter is given by the function B(Ω). We assume a Gaussian function for both the energy and angular response of the detector. Their FWHM as a function of wavelength are taken from Ref. [35]. In Eq. (1) σ λ denotes the spectral resolution.
The wavelength of emission can be computed as λ em = c/ν em with ν em = m a /(4π). We neglect the velocity dispersion of ALPs in the Leo T halo, since it is significantly smaller than the spectral resolution of MUSE. Indeed, the velocity dispersion is σ v 10 km/s [36], which means σ v /c 3 × 10 −5 , while the spectral resolution is ∆λ/λ 5 × 10 −4 .
In our analysis we scan over g aγγ and m a , while we take a model for the ALP spatial distribution. Analyses of the velocity dispersion in dSph often provides the so-called D-factor. It is defined as: where dΩ = 2π sinθdθ and the angular region of integration is identified by the angleθ max . For our computation, we are interested in the function φ = d ρ DM , which enters in Eq. 1. To obtain it we consider the function D(θ) determined in Ref. [36] from observational data, and then we invert Eq. 2. The uncertainties we will quote in Sec. 4 are derived from the uncertainty on the D-factor from Ref. [36]. The latter is about a factor of four (at 68% C.L.) at ∼ 1 , which translates in approximately a factor of two for what concerns the bound on g aγγ . The flux density in Eq.1 grows ∝ g 2 aγγ and depends on the ALP mass through the ALP decay rate (Γ a ∝ m 3 a ), the spectral resolution and the observing angular beam at the wavelength of emission. In Fig. 1 (right), we show the expected signal in Leo T, taking two examples of ALPs with couplings significantly smaller than any existing bound: g aγγ = 2 × 10 −12 GeV −1 for m a = 3 eV (green) and g aγγ = 10 −12 GeV −1 for m a = 5 eV (red). Comparing them with the measured average emission (blue line), it is clear that such models are in the ball-park accessible by MUSE observations.

Statistical Analysis and Results
We assume the likelihood for the ALP diffuse emission to be described by a Gaussian likelihood: where S i th is the theoretical estimate for the flux density in the pixel i, S i obs is the observed flux density and σ i rms is the r.m.s. error, both described in Sec. 2. The theoretical estimate is provided by Eq. 1 plus a spatially flat term S λ,f lat that we include in the fit to the map at each wavelength (to account for incomplete sky subtraction) and we treat as a nuisance parameter. N pix is the total number of pixels in the area under investigations, that is chosen to be a circle of 30 of radius. N F W HM pix is the number of pixels within the MUSE angular FWHM. Bounds on the parameter g aγγ are computed at any given mass m a through a profile likelihood technique, namely "profiling out" the nuisance parameter S λ,f lat . We assume that λ c (g aγγ ) = −2 ln[L(g aγγ , S b.f. λ,f lat )/L(g b.f. aγγ , S b.f. λ,f lat )] follows a χ 2 -distribution with one d.o.f. and with one-sided probability given by P = ∞ √ λc dχ e −χ 2 /2 / √ 2 π, where g b.f. aγγ denotes the best-fit value for the coupling at that specific ALP mass. In other words, the 95% C.L. upper limit on g aγγ at mass m a is obtained by increasing the signal from its best-fit value until λ c = 2.71, keeping S λ,f lat fixed to its best-fit value.
Results are shown in Fig. 2. The solid blue curve refers to the bounds derived assuming the best-fit D-factor from Ref. [36], while the shaded area considers their 68% C.L. interval for the D-factor. The curves show rapid variation at lower energies/longer wavelengths reflecting the presence of strong OH emission lines from the night sky which increases the noise here, as mentioned in Section 2. The limit on g aγγ becomes more stringent as the mass m a increases. This is due to the scaling of the decay rate as Γ a ∝ m 3 a (mitigated by the fact that also the background increases as the wavelength decreases).
In Fig. 2, we also show the bound derived in Ref. [20] from observation of clusters, in Ref. [37] from the ratio of horizontal branch (HB) to red giants stars in globular clusters and, for reference, the preferred region for the QCD axion [38] 5 . Our results improve existing bounds by more than one order of magnitude. They exclude the possible interpretation of near infrared background anisotropies in terms of ALP dark matter [23] in the wavelength/mass range covered by our analysis.
To test the robustness of our results against different masking and error estimates, we conduct a few sanity checks. Excluded at 95% C.L. Figure 2: The solid blue curve shows the 95% C.L. upper limits on the effective ALP-two-photon coupling gaγγ as a function of the ALP mass, derived in this work. The violet area includes the uncertainties on the D-factor, taken from Ref. [36]. We also show the bound derived in Ref. [20] from observation of clusters, in Ref. [37] from the ratio of horizontal branch (HB) to red giants stars in globular clusters, and the preferred region for the QCD axion [38].
First, we perform the same analysis mentioned above but on maps where we discard the last step of data manipulation described in Section 2. Namely, we do not mask the faint sources (which are derived through segmentation of the white-light image by SExtractor and by running ORIGIN for emission line sources). The spatially averaged spectrum of these maps is shown in green in Fig. 3 (left). The spectrum is well above zero at nearly any wavelength, meaning that there is a significant residual emission, not related to the ALP signal. On the other hand, the resulting bounds (green line in the right panel) are only very mildly less constraining than in the reference analysis (blue line). At few wavelengths, the bounds becomes slightly more constraining in the unmasked case because in the masked map the fit shows a preference for the DM component over the flat term.
As a second test, we consider a different derivation of the measurement uncertainties, by computing the standard deviation in a region of 2 × 2 around the pixel i, instead of using the procedure described in Section 2. We find this alternative derivation to provide, on average, a slightly more optimistic estimate of the errors with respect to the reference analysis. However, the bounds are only marginally different, as can be seen by comparing red and blue lines in the right panel of Fig. 3, supporting the reliability of our analysis.
Finally, we have also tested that varying the spectral resolution within its uncertainties has a negligible impact on the bounds.

Conclusions
Dark matter in the form of ALPs can be searched for exploiting their coupling to photons. One strategy is to look for the monochromatic photon flux generated by ALP decays inside astrophysical structures. Dwarf spheroidal galaxies are particularly suitable targets since they are dark matter dominated objects and they are relatively close to us. Excellent spectral resolution and sensitivity are required to search for the faint lines from ALPs. For this purpose, we have exploited the spectroscopic observations obtained with the MUSE instrument at the VLT. In this work we have searched for ALPs signals from the Leo T dwarf spheroidal galaxy. The data cover wavelengths in the range 470-935 nm and therefore allow to test ALP masses between 2.7 and 5.3 eV.
We have derived new bounds on the effective ALP-two-photon couplings, which improve existing constraints from observations of clusters [20] and stars in globular clusters [37] by more than one order of magnitude. We have investigated the impact of different sources of systematical uncertainties in the analysis, namely the masking of faint sources present in our image, the uncertainties on the measurements, and those on the spectral resolution. In this respect, we have shown that our constraints are rather robust.
In the future we plan to extend our analysis including MUSE observations of additional dwarf spheroidal galaxies. A joint analysis of several targets will allow not only to improve current exclusion limits, but it will also provide a better handle to unveil a possible ALP signal, which should obviously show up at the same wavelength in all the targets.