Upper Limits on the VHE gamma-ray Emission from the Willman 1 Satellite Galaxy with the MAGIC telescope

We present the result of the observation of the ultrafaint dwarf galaxy Willman 1 performed with the 17 m MAGIC telescope during 15.5 hr between March and May 2008. No significant $\gamma$-ray emission was found. We derived upper limits of the order of $10^{-12}$ ph cm$^{-2}$ s$^{-1}$ on the integral flux above 100 GeV, which we compare with predictions from several of the established neutralino benchmark models in the mSUGRA parameter space. The neutralino annihilation spectra are defined after including the recently quantified contribution of internal bremsstrahlung from the virtual sparticles that mediate the annihilation. Flux boost factors of three orders of magnitude are required even in the most optimistic scenario to match our upper limits. However, uncertainties in the dark matter intrinsic and extrinsic properties (e.g., presence of substructures, Sommerfeld effect) may significantly reduce this gap.


INTRODUCTION
Dwarf spheroidal galaxies (dSphs) are believed to be the smallest (size 1 ∼ kpc), faintest (luminosities 10 2 − 10 8 L ⊙ ) astronomical objects whose dynamics are dominated by dark matter (DM; Gilmore et al. 2008, and references therein). They are found as satellites orbiting in the gravitational field of a larger host galaxy (e.g., the Milky Way (MW)). Their member stars show large circular velocities and velocity dispersions that, combined with their modest spatial extent, can be interpreted with the presence of a large DM halo of the order of 10 5 − 10 8 M ⊙ . In recent years, the Sloan Digital Sky Survey (SDSS; York et al. 2000) led to the discovery of a new population of MW satellites, comprising about as many new objects as were previously known (Belokurov et al. 2004;Willman et al. 2005a;Zucker et al. 2006;Irwin et al. 2007;Walsh et al. 2007). This population of extremely low-luminosity galaxies is very interesting for DM searches and to study the galaxy formation at the lowest mass scale. The existence of a new class of ultrafaint MW satellites is also relevant because it provides a partial solution for the so-called missing satellite problem (Klypin et al. 1999;Simon & Geha 2007;Strigari et al. 2007;Madau et al 2008) by partially filling the gap between the predicted and the measured number of galactic subhalos. In this family is Willman 1, discovered by Willman et al. (2005a) and soon established as potentially the most DM-dominated dSph satellite of the MW (Martin et al. 2007;Sánchez-Salcedo & Hernandez 2007;Strigari et al. 2008). The physics of DM has gathered much interest in recent years, in particular after Wilkinson Microwave Anisotropy Probe (Spergel et al. 2007;Komatsu et al. 2008) measured its relic density with great precision. It is generally believed that DM manifests itself as a general class of weakly interacting massive particles that includes several candidates which satisfy both experimental constraints and theoretical requirements (see Bertone et al. 2005a, and references therein). Among them, one of the best theoretically motivated, for whom the relic density is calculated without fine tuning from its nature, is the neutralino, arising in SuperSymmetric (SUSY) theories beyond the Standard Model (Wess & Zumino 1974;Haber & Kane. 1985) and in particular in the mSUGRA extension (Chamseddine et al. 1982). The mSUGRA neutralino annihilations can be observed through the production of γ-rays. The main emission comes from secondary products of hadronization processes and from final state radiation. In addition, line emissions are found through direct processes such as χχ → γγ, and χχ → Z 0 γ, which provide γ-rays of energies E = m χ and E = m χ − m 2 Z 0 /m 2 χ respectively, even if those processes are at the loop level and therefore strongly suppressed.
Recently Bringmann et al. (2008a), following an earlier idea from Bergstrom (1989), showed that in some regions of the mSUGRA parameter space, a hitherto neglected contribution to γ-ray emission comes directly from charged sparticles mediating the annihilation into leptons, in processes like χχ → l + l − γ. They defined this intermediate state radiation as internal bremsstrahlung (IB). The IB mechanism permits to restore the helicity balance in processes that would otherwise be strongly forbidden, and supplies γ-rays toward high energies (E > 0.6 m χ ) of up to several orders of magnitude. This boosted emission is particularly interesting for groundbased γ-ray observation, as performed by Imaging Atmospheric Cherenkov Telescopes (IACTs), which are sensitive above ∼ 100 GeV where normally the IB boost takes place.
As the γ-ray flux is proportional to the square of the DM density, IACTs focus on concentrated DM objects. Dwarf galaxy satellites of the MW represent very good candidates, since they are the most DM dominated systems known in the universe, with very high mass-to-light ratios (M/L), close distance and reduced γ-ray background from unresolved conventional Galactic sources (i.e. stellar evolutionary remnants). Some dSphs have already been studied in γ-rays: Draco by MAGIC (Albert et al. 2008a) and Whipple (Wood et al. 2008); UMi by Whipple (Wood et al. 2008) and Sagittarius by HESS (Aharonian et al. 2007) without any sig-nificant observation of DM annihilations and only flux upper limits were estimated.
In this paper, we report results of the observation of the sky region around Willman 1 performed by the MAGIC telescope for a total of 15.5 hr between March and May 2008. After a brief description of Willman 1 in Section 2, in Section 3, we estimate the flux using benchmark models for the neutralino and a typical DM density profile. In Section 4, we briefly describe the MAGIC telescope and the Willman 1 data sample. In Section 5, we present and discuss the results of the observation and we set upper limits for the flux. In Section 6, we report our conclusions.  Willman et al. 2005a) as a faint overdensity of red, resolved stars, which was observed again the next year (Willman et al. 2005b) and named Willman 1. At that moment, this object represented the 10th dSph of the MW and, the first one discovered in 10 years. Further observations performed with the Keck/DEIMOS telescope confirmed the SDSS results (Martin et al. 2007), while a more recent observation is reported by Siegel et al. (2008).
Willman 1 is located at a distance of 38 ± 7 kpc in the Ursa Major constellation. It is characterized by a very low number of resolved stars, the total luminosity being L = 855 L ⊙ , and a very small half-light radius of r 1/2 = 21 ± 7 pc, almost 2 orders of magnitude smaller than other known dSphs. The source was defined by previous authors as an "extreme" dwarf galaxy, because some of its characteristics lie between those typical for a globular cluster (GC) and those expected in an extremely faint dSph. The large spread in metallicity of its stars favors the dSph interpretation rather than that of a GC, which would contain stars of a similar age and metallicity (Martin et al. 2007), even if this evidence was recently put under discussion by Siegel et al. (2008). Willman 1 is the least massive satellite galaxy known to date, with a total mass (M ∼ 5 × 10 5 M ⊙ ) about an order of magnitude smaller than those of the least massive satellite galaxies previously known. However, the corresponding mass-to-light ratio, M/L∼ 500 − 700 M ⊙ /L ⊙ , is one of the highest of all dSphs. This makes Willman 1 one of the most attractive dSph galaxies to look for DM at present (see, e.g., Bringmann et al. 2008b), its predicted DM annihilation flux being at least a factor of 3 larger than the second best DM candidate, according to recent work (Strigari et al. 2007).

THEORETICAL MODELING OF THE GAMMA-RAY EMISSION FROM WILLMAN 1
The γ-ray flux originating from DM particle annihilations can be factorized into a contribution called the astrophysical factor J(Ψ) related to the morphology of the emission region and a contribution called the particle physics factor Φ P P depending on the candidate particle characteristics: where E 0 is the energy threshold of the detector and Ψ the angle under which the observation is performed.
3.1. Astrophysical Factor At present, a concise and exact characterization of the DM density profile of Willman 1 is a delicate task, since observational data are still scarce. This paper is based on the studies from Strigari et al. (2008), who modeled the profile by using only 47 stars after removing those with unclear kinematics. Furthermore, to avoid membership problems, only the observational data related to the inner half of the galaxy were taken into account. It is important to note that a null/insignificant tidal stripping was assumed in order to carry out the modeling of the DM distribution as a system in dynamical equilibrium, a fact which is still under debate. For example, Willman et al. (2006) claim the existence of strong tidal debris, the evolution of the dwarf being strongly affected by tidal interactions with the MW still now, although DM constitutes 90% of its total mass. In the same line, Martin et al. (2007), following deep observations in the r band, infer that Willman 1 may probably be surrounded by tidal tails. The authors give as a plausible scenario that the dwarf could have been significantly tidally stripped but only in the past, when the object was more luminous and massive. At that age, Willman 1 could have lost most of its outskirts, only the innermost regions surviving intact. This picture would then allow the two contradictory arguments to coexist, since at least a correct modeling for the core of the dwarf may be possible assuming this region to be in dynamical equilibrium at present. A well established model for the DM distribution, and the only one so far applied to Willman 1 (Strigari et al. 2008), is the Navarro-Frenk-White (NFW) density profile (Navarro et al. 1997): where ρ s and r s are a typical scale density and radius respectively. The astrophysical factor can be written as: where Ψ 0 denoted the direction of the target. The first integral is performed over the spatial extension of the source, the second is performed over the line of sight variable λ. The density is convoluted with the Gaussian function B ϑr (θ) describing the telescope angular resolution where θ = Ψ − Ψ 0 is the angular distance with respect to the center of the object. We remark that the integration of Equation 3 involves foreground (MW halo) and extragalactic background whose contributions can be substantial (Elsaesser & Mannheim 2005). At a distance of 38 kpc, the scale radius corresponds to an extension of 0.54 • in the sky, which is well inside the MAGIC field of view (∼ 3.5 • ), but is rather extended compared with the telescope angular resolution of 0.1 • . This evidence is mitigated by the fact that the main emission still comes from the very core of the source, due to the very steep NFW profile at the center and the square density dependence. For this reason, we performed an analysis adapted for slightly extended sources (Sitarek & Mirzoyan 2008). To compute the astrophysical factor for Willman 1, we substitute r s = 0.18 kpc and ρ s = 4 × 10 8 M ⊙ /kpc 3 taken from Strigari et al. (2008) into Equation 2, and by computing Equation 3 we obtain J(Ψ 0 ) ∼ 3.5 × 10 17 GeV 2 /cm 5 . This value is possibly the largest over the rest of dSphs.

Particle Physics Factor
In many mSUGRA models, the lightest SUSY particle is one of the four neutralinos (χ 0 1...4 ) which are linear combinations of the BinoB 0 , the WinoW 0 and the two neutral higgsinosH 0 u ,H 0 d . Only five free parameters fully characterize mSUGRA: the scalar mass m 0 and the gaugino mass m 1/2 defined at the unification scale, the trilinear scalar coupling A 0 and the ratio tan β of the Higgs vacuum expectation values. In addition, one needs to define the sign of the Higgs mass parameter sign(µ). The typical DM γ-ray annihilation spectrum is a continuum with a sharp cutoff at the DM candidate mass and possibly with a bump at energies larger than 0.6 m χ in case the IB is present (see Figure 1).
The particle physics factor can be written as a product of two terms, the first depending only on the DM candidate mass and cross section, and a second term, depending on the γ-ray spectra of annihilation, which must be integrated above the energy threshold E 0 of the telescope: where σv χχ is the total averaged thermal cross section times the relative velocity of particles, m χ is the DM particle mass, and the factor 2 takes into account that the neutralino annihilates with itself. The γ-ray annihilation spectrum is composed of different contributions: The mSUGRA parameter space is conventionally described in a m 0 ⊕ m 1/2 plane, after having fixed the other free parameters. Usually, four zones are identified: the bulk region, with low m 0 and m 1/2 and neutralino masses at around 100 GeV, the focus point where m 0 and the neutralino are more massive, the funnel region where both m 0 and m 1/2 take large values, and the co-annihilation tail characterized by large m 1/2 . A neutralino shows different annihilation modes depending on the location in this plane. A representative set of benchmarks was defined by Battaglia et al. (2001Battaglia et al. ( , 2004. Hereafter, we use a subset of four slightly modified Battaglia models, as defined by Bringmann et al. (2008b) using DarkSUSY 4.01 (Gondolo et al. 2004), which include the contribution of IB in the computation of the cross sections and spectra: models I', J', K', F* for the bulk, co-annihilation, funnel and focus point regions respectively. All the defining parameters, as well as the derived Φ P P (> E 0 ), are listed in Table 1.
The total estimated flux due to DM annihilation computed according to Equation 1 is given in Table 2. These estimations must be taken with carefulness since there are some sources of uncertainty that may largely affect the values for the predicted flux: (1) our lack of knowledge of the DM density profile may change the astrophysical factor by more than 1 order of magnitude; (2) the presence of DM substructures in the dwarf may enhance the γ-ray flux at least by a factor of 2 − 3 according to N −body simulations (Diemand et al.  Bringmann et al. (2008b) and computation of the particle physics factor: m 1/2 and m 0 [GeV] are the gaugino and scalar mass respectively defined at the unification scale; tan β is the ratio of the Higgs expectation values; A 0 [GeV] the trilinear coupling constant and sign(µ) the sign of the Higgs mass; mχ [GeV] is the neutralino mass; σvχχ [cm 3 s −1 ] is the cross section times the relative velocity of DM particles and Φ P P (> 100) [cm 3 GeV −2 s −1 ] is the particle physics factor above 100 GeV.  2007a,b; Kuhlen et al. 2008) or even up to 20 according to Martinez et al. (2009). Substructures are in fact expected to be present not only in Willman 1 but also in any other DM halo, since cold dark matter (CDM) halos are approximately self-similar until a cutoff scale mass which lies in the range of 10 −4 -10 −12 M ⊙ (Profumo et al. 2006); (3) the exclusion of the baryons in the modelization of the total density profile. However, the effect of the adiabatic compression (Prada et al. 2004;Gnedin et al. 2004), although important for larger DM halos, will probably play a marginal role in the case of Willman 1, given its relatively low amount of baryons even in the central regions, where the effect is expected to be more important; (4) the inclusion of the Sommerfeld effect, recently investigated by Lattanzi & Silk (2008). This effect predicts a general enhancement S of the flux, in case of very low DM particle velocities v, proportional to the inverse of the particle velocity: S ∝ 1/v. In addition, for particular values of the DM mass and very small velocities, due to the presence of bound states, the Sommerfeld effect would give rise to resonances, which enhance the flux of a factor S ∝ 1/v 2 . Finally, we underline that our choice of benchmarks does not scan the complete parameter space and different neutralinos could have a larger expected flux. A deeper study of the parameter space is therefore very desirable.

MAGIC DATA
The MAGIC (Major Atmospheric Gamma Imaging Cherenkov) telescope is located on the Canary Island La Palma (2200 m asl, 28.45 • N, 17.54 • W). MAGIC is currently the largest IACT, having a 17 m diameter tessellated reflector dish. The faint Cherenkov light flashes produced by atmospheric showers initiated by VHE γ-rays in the top atmosphere are reflected into an ellipsoidal image in the focal plane of the telescope, where a camera consisting of 577 photomultipliers (pixels) records the image. More details can be found in Cortina et al. (2005).
MAGIC observed Willman 1 between March and May 2008. The source was surveyed at zenith angles between 22 • and 30 • , which guarantees the lowest energy threshold. The source was tracked for 16.8 hr plus another 9.3 hr in OFF observation mode, i.e. pointing to a dark patch in the sky close to Willman 1 where no γ-ray emission is expected, for background estimation. The main background Cherenkov telescopes have to deal with is produced by cosmic hadronic particles impinging onto the top atmosphere and generating electromagnetic subshowers that can mimic pure γ-rays showers, and by the night sky background. Background events are partly rejected at the trigger level and in the off-line analysis event selection, following a procedure called gamma/hadron (g/h) separation.
The analysis proceeds as follows (for a detailed description, see Albert et al. 2008b): data are calibrated and the number of photoelectrons per pixel extracted (Albert et al. 2008c), then an image cleaning selects pixels with at least six photoelectrons (three photoelectrons in the boundary of the image). Additional suppression of pixels containing noise is achieved by requesting a narrow time coincidence between adjacent pixels (∼ 7 ns). Based on the Hillas parameterization algorithm (Hillas 1985), the shower parameters are reconstructed. The hadronic background is suppressed with a multivariate method, the Random Forest (Breiman 2001;Albert et al. 2008d), that uses the Hillas parameters to define an estimator called hadronness by comparison with Monte Carlo (MC) γ-ray simulations. The hadronness expresses the likeness of an event to be a hadron and runs from 0 for gammas to 1 for hadrons. The Random Forest method is also used to estimate the energy of a reconstructed shower and the energy threshold is defined by the peak of the distribution of reconstructed MC gamma events. The g/h separation is optimized on a real data sample from the Crab Nebula, a supernova remnant and one of the brightest and stable γ-ray emitters, which is taken as standard candle in very high energy γ-ray astronomy. The optimization yields a best set of cuts in the Hillas parameters which defines the gamma and hadron acceptance of the analysis. In our case, the optimal set of cuts is obtained for an energy threshold of 100 GeV and a hadronness cut of 0.15. The overall data quality is very high, with only 7% data rejection, resulting in 15.5 hours effective observation time. Independent cross-checks were performed on the data giving compatible results.

RESULTS AND DISCUSSION
No significant γ-ray excess beyond 100 GeV above the background was observed in 15.5 hr of observation of the sky region around Willman 1, according to Figure 2, where the "α-plot" is reported. The α-parameter is the angular distance between the shower image main axis and the line connecting the image barycenter and the camera center. Due to their isotropic origin, hadronic events, in case they survive the analysis, are randomly oriented in the camera both in the ON-data and OFF-data sample. This is reflected into a rather smooth distribution of events in the α-plot, the nonperfect flatness being due to an increased camera acceptance for showers with small α. On the other hand, γ-rays trace back the source, and thus the orientation of the shower image is toward the center of the camera. Therefore, in case of positive detection, an excess of events in the ON-data above the OFF-data sample is expected at small α. A fiducial region α < 12 • is chosen where the signal is assumed with a cut slightly larger than for a pointlike source to take into account the moderate source extension. The OFFdata are normalized to the ON-data in the region where no signal is expected, i.e., between α = 30 • and α = 80 • . The significance is calculated using Equation 17 of Li & Ma (1983). The number of excess events N exc (> 100 GeV) = −223 ± 223 is calculated as the difference between the number of ON-events and the number of OFF-events in the fiducial α−region. The method from Rolke et al. (2005) was applied to estimate the upper limit in the number of excess events with a 90% confidence level and including 30% of systematic errors, giving as a result N ul exc (> 100 GeV) ∼ 191.4 events. This value is used to reconstruct the corresponding photon flux for a general γ-ray spectrum S(E) using: where E 0 GeV is the energy threshold, A cuts ef f (E) is the effective telescope area and ∆T is the effective observation time.
We applied Equation 5 for the four neutralino benchmarks defined. Results are reported in Table 2, where we also compare Φ u.l.
>E0 with the estimated flux of Table 1. We also report the upper limits on the boost factors that are required to match the two fluxes, again calculated separately for each neutralino model. The boost factor is defined as the ratio between the upper limit and the theoretical flux, and defines the minimal boost that the theoretical flux should be subject to in order to allow for a positive detection of the source. In order to provide results less dependent of the particular benchmark spectrum, we also calculated flux upper limits in four different energy bins [100-170, 170-350, 350-1000, 1000-20 000] GeV for a generic annihilation spectrum without cutoff and spectral index −1.5. Respectively, the resulting upper limits are [ 9.94, 4.75, 0.68, 0.35 ]×10 −12 photons cm −2 s −1 . Table 2 reveals that although we derived upper limits of the same order of magnitude for the four models considered, there are evident differences in the prospects of detection for different neutralinos. The boost factor largely depends on the benchmarks, but the main differences are connected to the particle physics factor, which varies by orders of magnitude among the different benchmarks, as shown in Table 1. The best prospects are for neutralinos in the f unnel region (model K ′ ) of the parameter space, for which the mass is large enough to place the cutoff well within the MAGIC energy threshold but still small enough not to reduce the particle physics factor Φ P P of Equation 4 too much. Next, with similar boost requirements follow the I ′ and F * models. In the former case, the effect of the IB plays an important role at energies close to the cutoff even if the neutralino mass is very close to the MAGIC energy threshold. In the latter case, although the IB effect is negligible, the signal is very extended in the energy region suitable for MAGIC, whereas the large mass makes the flux suppression too large. The worst case scenario is the co-annihilation neutralino. In this case, even if the IB contribution is large, the intrinsic total cross section is very low, which makes the flux very low compared with the others. The IB effect cannot counteract this intrinsic deficit.
Although the results of Table 2 are at first glance not looking overly promising for indirect detection of the discussed models in the immediate future, one should keep in mind that the estimation of the flux is model dependent and there are intrinsic effects that could boost up the flux for up to 1 order of magnitude if the effect of substructures, astrophysical uncertainties, and/or the inclusion of baryons is considered (as already discussed in Section 3.2) and for several orders of magnitude, if the DM particle velocity is small enough to give rise to the Sommerfeld effect. This could change completely the outlook, as also discussed in Pieri et al (2009).

CONCLUSION
In the context of DM searches, we have observed the Willman 1 dwarf galaxy with the MAGIC telescope for a total of 15.5 hr between March and May 2008. Willman 1 represents one of the best DM-dominated systems known in the universe to search for DM at present, ac -TABLE 2 Comparison of estimated integral flux above 100 GeV using Equation 1 for the benchmarks models defined in Table 1 and the upper limit in the integral flux Φ u.l. above 100 GeV coming from MAGIC data in units of photons cm −2 s −1 . On the rightmost column, the corresponding upper limit on the boost factor B u.l. required to match the two fluxes is calculated.
No γ-ray signal was detected above an energy threshold of 100 GeV. We have obtained different flux upper limits of the order of 10 −12 ph cm −2 s −1 separately for four benchmark models considered in the framework of mSUGRA.
Using the latest estimations of its structural parameters to build the DM density profile as well as the inclusion of the recently proposed IB mechanism, we calculated the boost factors needed to match the expected flux values from DM annihilation in Willman 1 with the upper limits obtained from the data. We can see that boosts in flux in the order of 10 3 are required in the most optimistic scenario considered. However, uncertainties in the DM distribution, the role of DM substructure or/and the contribution from the Sommerfeld enhancement may reduce this required boost significantly. It is expected that deeper observations of the Willman 1 dSph with the upcoming MAGIC II telescope will allow us to improve the flux limits presented here by a factor 2 − 10.