Constraints on dark matter annihilation by radio observations of M31

18 Pages·2013·5.67 MB·English

Constraints on dark matter annihilation by radio observations of M31
A.E.Egorov∗ and E.Pierpaoli†
Department of Physics and Astronomy, University of Southern California,
3620 McClintock Ave., SGM 408, Los Angeles, CA 90089, USA
(Dated: April3, 2013)
WeusedradioobservationsoftheneighbourgalaxyM31inordertoputconstraintsondarkmatter
particlemassandannihilationcrosssection. Darkmatterannihilation inM31haloproduceshighly
energetic leptons, which emit synchrotron radiation on radio frequencies in the galactic magnetic
field. We predicted expected radio fluxes for the two annihilation channels: χχ → b¯b and χχ →
τ+τ−. WethencomparedthemwithavailabledataonthecentralradioemissionofM31asobserved
3 by four radio surveys: VLSS (74 MHz), WENSS (325 MHz), NVSS (1400 MHz) and GB6 (4850
1 MHz). Assuming a standard NFW dark matter density profile and a conservative magnetic field
0 distributioninsidetheAndromedagalaxy,wefindthatthethermalrelicannihilationcrosssectionor
2 higherhσvi≥3·10−26 cm3/sareonlyallowed forWIMPmassesgreaterthan≈100GeVand≈55
r GeV for annihilation into b¯b and τ+τ− respectively. Taking into account potential uncertainties in
p thedistributionsofDMdensityandmagneticfield,thementionedWIMPlimitingmassescanbeas
A lowas23GeVforbothchannels,andashighas280and130GeVforannihilationintob¯bandτ+τ−
respectively. These mass values exceed the best up-to-day known constraints from Fermi gamma
2
observations: 40 GeV and 19 GeV respectively [A.Geringer-Sameth and S.M.Koushiappas, Phys.

Rev. Lett. 107, 241303 (2011)]. Precise measurements of the magnetic field in the relevant region
]
O and betterreconstruction of theDM densityprofile of M31 will beable toreducetheuncertainties
of our exclusion limits.
C
. PACSnumbers: 95.35.+d,98.56.Ne,95.85.Bh
h
p
-
o I. INTRODUCTION review). Conventional ΛCDM model of cosmology sug-
r gests WIMP velocity averagedannihilation cross section
st Thephysicalnatureofthedarkmatter(DM)iscontin- σv 3 10−26 cm3/s, which yields the correct cur-
h i ≈ ·
a uingtobenotunderstood. Themostprobablecandidate rent DM abundance in the Universe in case the DM is
[
for the role of DM is a Weakly Interacting Massive Par- thermally produced (see e.g. [1]). Indirect searches may

1 ticle (WIMP): a supersymmetric partner of a Standard either lead to the discovery of DM with given particle
v Model(SM)particle. Therearethreemainapproachesin annihilation cross section and mass 1, or constrain the
7 attempts to detect WIMPs: accelerator searches, direct crosssectionto be below the levelofa thermalrelic over
1
searches and indirect searches. Although great efforts all plausible range of WIMP masses (from several GeV
5
0 have been dedicated already to all of these directions, to several TeV). The latter case may imply that in fact
. the DM puzzle still seems to be very far from its final wedonotunderstandtheDMphenomenonincosmology
4 and clear solution. all that well, and a more sophisticated DM production
0
This article is related to the last mentioned search mechanism is needed.
3
1 strategy: indirect detection. The idea of indirect detec- Currently the most promising direction in indirect de-
: tion is based on the opportunity of pair annihilation of tection seems to be related to Fermi searches of the pri-
v WIMPs causing production of different highly energetic mary gamma emission produced directly by WIMP an-
i
X SM particles, which eventually decay into stable parti- nihilation(see[3]). Atthemoment,theseconstraintsare
r cles like leptons, protons and others. These yields may the strongest: they exclude WIMP masses smaller than
a emit electromagnetic radiation through various mecha- about 40 GeV and 19 GeV for b¯b and τ+τ− annihila-
nisms like synchrotron emission, Inverse Compton Scat- tion channels respectively for the standard thermal relic
tering (ICS), etc. in different astrophysical objects. By annihilation cross section. While the allowed region of
detecting such radiation we can infer WIMP properties. WIMP parameters space is gradually shrinking, WIMP
Since there are evidences of DM existence in any ob- massesof hundreds andthousands GeV are still allowed.
ject from dwarf galaxies to largest galaxy clusters, any Indirect searches may allow further insights in the DM
of these objects is potentially a good target for indi- phenomenon. Wavelength bands ranging from the radio
rect DM searches. Indirect searches of dark matter have tothegammaraysandvariousastrophysicalobjectsmay
been already extensively exploited considering all wave- beexploitedfurthertothisaim. Inthisworkwestudied,
length ranges from radio to gamma rays (see [1] for a
1 Note that annihilation cross section at the current epoch can
∗ [email protected] be larger in more general models, than the standard thermal
† [email protected];CaliforniaInstituteofTechnology,MC249-17, value even in the case of thermal production due to effects like
1200EastCaliforniaBlvd,Pasadena,CA91125, USA Sommerfeldenhancement -seee.g. [2].
2
inparticular,radioobservationsofthe closestbiggalaxy computation of the expected radio emission properties
M31. Relativistic leptons produced by DM annihilation due to WIMP annihilation in M31, in section III we de-
in its halo emit synchrotron radiation on radio frequen-
cies due to Andromeda’s magnetic field. Given the size
of the galaxy and its DM halo, the strength of its mag-
netic field and its close distance, large radio fluxes from
DM annihilation may be expected from M31. There-
fore M31’s radio properties may be exploited in order
to put strong upper limits on the annihilation cross sec-
tion. Surprisingly, almost no attempts in this direction
have been made in the past: in the literature only few
articles can be found dedicated to indirect searches in
M31. The most remarkable among them is [4]. These
authors obtained some constraints on WIMP parameter
space by Cherenkov ground based gamma observations.
However, those observations did not have enough sen-
sitivity to probe the relevant region of DM parameters,
and since then no significant progress has been made on
this object.
WecomputedtheexpectedradiofluxduetoDManni-
hilation(sectionII),andthencompareditwithavailable
data of radio observations of M31, which allowed us to
putupperlimitsontheannihilationcrosssection(section Figure 1. Optical image of M31 with the marked region of
III). For comparison with real observations we chose all interest (ROI),selected for ourpurposes. Circle radius is 5′.
appropriate radio surveys, which cover a wide range of
frequencies: VLSS(74MHz),WENSS(325MHz),NVSS
(1400 MHz) and GB6 (4850 MHz). The limits obtained rived the actual constraints by comparison of the pre-
can be considered as conservative because we did not dicted fluxes with observational data, and section IV
make any specific assumptions about the radio emission summarizes the results of our work.
other than the one from DM in the center of M31, and
allowedfor an unconstrained contribution from all other
unknown backgrounds. II. COMPUTING THE RADIO FLUX
We calculated the constraints for two annihilation
channels, specifically the ones annihilating into b¯b and A. General theory
τ+τ− pairs. We chose these channels among allpossibil-
itiesbecause,asexplainede.g. in[5],b¯bandτ+τ− nearly Inthissectionwepresentedtheprocedureforthecom-
present the channels with the softest and hardest lepton putationoftheradiofluxdensityfromthecenterofM31.
yields respectively. Any other case would therefore pro- Weneglectedherethepotentialabsorptionofradioemis-
duce radio fluxes at intermediate levels with respect to sion between the source and the observer, since our es-
these two. In this sense, these can be considered as the timates showed that it occurs at a negligible level (see
two limiting cases. appendixA).Inthecaseofanopticallytransparentemit-
In our analysis we considered only the central part of tingmedium,thetotalfluxdensityfromourROIinM31
M31, and specifically the bulge area of circular shape can be obtained just by integrating the local medium
withtheangularradiusα 5′aroundgalacticcenter(see emissivity j(ν,~r) over the volume of halo contained in
≈
fig. 1). We chose this specific region of interest (ROI) our ROI:
as a target of indirect searches because of the following
j(ν,~r)dV
considerations: i)theradioquietnessoftheM31nucleus, S = , (1)
which indicates low contamination in the radio band by R 4πd2
otherstandardastrophysicalprocesses;ii)theabsenceof
whered=785 25kpcisthedistancebetweenusandthe
any projected point source inside it; iii) sufficient halo ±
Andromedacenter([6]),~risthepositionvectorinsidethe
size to produce a relevantsignal. More details about the
M31 halo originating in the M31 center. We disregarded
ROI choice will be explained below.
here all (small) redshift effects. Then as a next step we
Through our paper we adopted the Hubble constant needed to compute the local emission coefficient at an
value H0 = 71 km/(sMpc), which was taken from the arbitrary position in M31 halo j(ν,~r). The synchrotron
·
WMAP7 data at http://lambda.gsfc.nasa.gov/. emissivityofleptonsproducedbyWIMPannihilationhas
This paperis organizedasfollows: sectionII describes the form:
3
mχc2 mχc2
j(ν,~r)= Z (cid:18)Pe+(Ee+,ν,~r)ddEnee++dEe+ +Pe−(Ee−,ν,~r)ddEnee−−dEe−(cid:19)=2 Z Pe(E,ν,~r)ddnEedE, (2)
mec2 mec2
where e+ and e− represent positrons and electrons re- injected leptons with the size of our ROI and concluded
spectively, P (E,ν,~r) is the synchrotron emission power that it could be done without significant effect on our
e
of one lepton with energy E ona frequency ν (measured results(seeappendix Bfordetails). Withthese twosim-
in erg/(sHz) in CGS), dne is the energy distribution of plifications we can easily solve the diffusion equation:
· dE
leptons - the number ofleptons per unitvolume per unit
energy range. We assumed electron and positron terms mχc2
dn 1
in eq. (2) to be equal to each other, which is reasonable e(E,~r)= q (E′,r)dE′. (6)
e
dE b(E,~r)
becausethesebothspeciesbehavesimilarlyinallrelevant Z
E
aspects. According to [7]:
At this point we should specify the functions b(E,~r)
π and q (E,r). The energy loss rate b(E,~r) is constituted
x e
P (E,ν,~r)= dθ′π√3sin2θ′r m cν F , (3) mainly by four differentcoolingprocesses: ICS emission,
e e e 0 sinθ′
Z0 (cid:16) (cid:17) synchrotronemission,bremsstrahlungandCoulombscat-
tering:
where
b(E,~r)=b +b +b +b . (7)
ICS sync brem Col
r = e2 , Using previousresults presentedinthe literature([7], [1]
e mec2 and [9]) we constructed each of these four terms respec-
 x= 2ν 1+ γνpl 2 3/2, tively as follows:
 νν0pl(~r=)3ν=0γ2π2eeπmB2(cid:16)mn(e~re)−c,p(cid:0)laνsm(cid:1)a(cid:17)frequency of ambient plasma, bICS =0.761 eVU·pchm−3 (cid:18)1 GEeV(cid:19)2, (8)
F(t) q1.25(648)1/12t1/3e−t synchrotron spectrum, here and in other loss terms numerical pre-factors like
 γ = m≈Eec2 −lepton Lorentz-fa−ctor. (4) 0co.7r6refsoplloonwdsinfrgomlostsheraetxeacfotranaallyeptitcoanl.comUpphutdaetnioonteosftthhee
This set of formulas describes the synchrotron power total energy density of radiation at relevant locations.
of a lepton in magnetic field B(~r) in the presence of an For the galactic center regionit is constituted mainly by
ambient plasma with concentration n. Integration over star light photons. We substituted for this quantity the
θ′ in eq. (3) represents the averaging over all possible fixed value U =8 eVcm−3, which is quoted in [1] as a
ph
·
random angles between lepton’s velocity and magnetic characteristic value for the MW center.
field. In order to proceed our computation of the radio
B(~r) 2 E 2
fluxexpectedweneededtoderivetheenergydistribution b =0.025 , (9)
sync
dne. For this purpose we applied the standard diffusion 1 µG 1 GeV
dE (cid:18) (cid:19) (cid:18) (cid:19)
equation:
n E
∂ dne = D dne + ∂ b(E,~r)dne +q (E,r), bbrem =4.701 cm−31 GeV, (10)
e
∂t dE ∇ ∇dE ∂E dE
(cid:18) (cid:19) (cid:18) (cid:19) where n denotes the concentration of ambient plasma
(5)
in cm−3.2 For this quantity we also used the constant
where D is the spatial diffusion coefficient, b(E,~r) is the
value n = 0.1 cm−3 relying on the results of study [10],
energylossrateforaleptonthroughvariousenergydissi-
wheresuchplasmaconcentrationwasderivedasatypical
pation mechanisms (measured in the units of energy per
averagevalue for the Galactic center region.
time), q (E,r) is the source function - how many elec-
e
trons (or positrons - but only one of these two species) n 1 E
are produced by WIMP annihilation per unit time per bCol =6.131 cm−3 1+ 75log m c2n . (11)
unit volume per unit energy range. We can simplify this (cid:18) (cid:18) e (cid:19)(cid:19)
equation assuming the stationary limit ∂ = 0, as com-
∂t
monly done for diffusion of DM annihilation products in
galaxiesandclusters(seee.g. [8]). Anotherusefulsimpli- 2 Thislosstermbbrem waspresentedinanother forminthework
[7]. However, wesuspect a mistake at this point in [7], because
ficationistheabsenceofthespatialdiffusionofannihila-
our verification of this loss rate did not confirm the expression
tionproducts,whichwillmakethefirstterminther.h.s. there. For this reason, we use here the expression precisely de-
vanish. We investigated the validity of such assumption rivedin[9](formula(6.74)there), whichappears tohavebetter
bycomparingthecharacteristicdiffusionlengthofnewly justification.
4
Alllossratesherearemeasuredintheunits10−16GeV/s.
For illustrative purposes, we showed in fig. 2 the energy
dependences of each term in eq. (8)-(11) over the rele-
Synchr.(Ρ’=0, z’=0)
vantrangeofleptonenergies. Wecanclearlyseethatthe
Synchr.(Ρ’=1 kpc, z’=0)
synchrotronlossestermdepends onthe locationthrough
B(~r). For this reason, we presented bsync on fig. 2 for ICS
three relevant trial locations, which are away from the Synchr.(Ρ’=0, z’=2 kpc)
M31 center by 0, 1 and 2 kpc. More details about B(~r)
distribution will be discussed in subsection IIC. In the
most general case, the quantities Uph and n depend on Bremsstr.
~r as well, indeed. However, at the current level of accu-
racyofourmodelwemade twosimplifying assumptions.
First, we assumed constant values for these quantities Coulomb
with conservative choices; and second, we adopted the
values for the MW relying on high similarity between
these two galaxies.
As for the source function q (E,r), it is computed es-
e
sentially as a product of number of leptons produced by
one WIMP annihilation per unit energy range dNe and
dE
number of annihilations per unit time per unit volume:
2
1 ρ (r) dN
q (E,r)= DM σv e, (12)
e
2 m h i dE Figure 2. The energy dependence of various types of lepton
(cid:18) χ (cid:19)
energy losses. Synchrotron losses are shown for three repre-
where ρ (r) is the local DM density, σv is thermally sentative locations around the M31 center. For more details
averagedDMannihilation cross section anhd mi is WIMP see subsections IIA and IIC. Color version of our article is
χ
mass. Annihilation yields of leptons dNe for two neces- available online.
dE
sary channels χχ b¯b, τ+τ− were taken from the web-
→
site http://www.marcocirelli.net/PPPC4DMID.html. Some other profiles potentially describe a DM distribu-
A detailed description of this resource is presented in tion in galactic halos even better than NFW - this was
[11]. We took the pre-computed annihilationyields from demonstrated e.g. in [14] for the Einasto density pro-
this resource in their last version, which includes elec- file. However, we did not use any other profiles except
troweakcorrections. Theimportanceofthesecorrections NFWfurtherdue tothe followingreasons. First,allcur-
for a yields computation was justified in [12]. Inclusion rently available models of the DM density distribution
of electroweak corrections leads to slightly different sec- in M31 were made in the frame of NFW profile. The
ondaryleptonyieldsfromannihilationintoprimaryprod- secondreasonis that amongallprofiles,the NFW yields
ucts with different polarizations: left- or right-polarized moderate emission fluxes due to DM annihilation (lower
tau-leptons, transversely- or longitudinally-polarizedW- than e.g. Moore profile - see [8]), and this makes our
bosons,etc. Butaspracticeshowed,thedifferenceinthe constraintsconservativewith respectto the assumptions
final radiation flux is very minor due to this splitting. on the adopted density profile. Also, the NFW profile is
For example, the final flux for the case of M31 DM halo specifiedbyonlytwofreeparameters,whilealargernum-
concentration c100 = 12, observational frequency ν = 74 berofparametersisrequiredbysomeotherprofiles. This
MHz, WIMP mass mχ =100 GeV and annihilation into choicethereforeminimizestheoverallnumberof“degrees
τ+τ− differs just by < 1% for different polarizations; of freedom” in our model.
whichismuchsmallerthanothertotalflux uncertainties Quantitatively the DM density distribution of our
presented. That’s why further we did not distinguish choice is the following:
different polarizations of primary annihilation products.
For the next step we needed to specify the DM den- ρ (r =50 pc), if 06r <50 pc;
DM
sity distribution ρ (r), which appears in formula (12). ρDM(r)= ρs , if r >50 pc. (13)
Next subsection isDMdedicated to the DM density profile (rrs(1+rrs)2
derivation. ItpresentstheexactNFWprofilewithaminormodifica-
tion,whichcanbecalledasaflatcore: wedonotextrap-
olate the density below the radialdistance of r =50 pc,
t
B. DM density distribution insidethiscentralregionweleaveDMdensityatthecon-
stant level of ρ (r = 50 pc). This choice is motivated
DM
As a model ofDM density distribution ρ (r) inM31 by the divergence of the NFW distribution towards the
DM
haloweusedthestandardNFWprofile,firstlyintroduced central point of the halo, which prevents to estimate the
in [13] and widely used since that due to its universality. DMdensityreliablyintheverycentralregion. TheNFW
5
profile providesreliable density estimate downto the ra- and numerous, that’s why they are able to substantially
dius 50 pc for a M31–size halo (see e.g. fig. 4 in [15]). increase the total flux due to DM annihilation from a
∼
That’s why in order to be conservative in our computa- whole halo - by 10-100 times (see e.g. [17]). But in our
tion of the medium emissivity, which depends quadrati- case we can neglect substructure contribution, because
callyonthelocalDMdensity,wetruncatedNFWprofile the main part of our expected flux comes from very cen-
at r = r = 50 pc and left the density at smaller radii tral region with a size about 1 kpc. As can be seen e.g.
t
on the level of ρ (r = 50 pc). We tested the depen- in [17], substructures do not survive so close to the cen-
DM
denceofourresultsonthespecificarbitrarychoiceofthe ter due to tidal disruption, and their contribution to the
truncation radius r . This study showed that a decrease total expected flux would be negligible at so small radii.
t
(increase) of r by a factor two (a factor four) produce Thus we did not need to include subhalos in our anal-
t
an increase (decrease) of 3% ( 10%) in final fluxes, ysis. Moreover, subhalos’ presence would only increase
≈ ∼
which is far below the overall level of accuracy of our theoverallflux,sothattheconstraintsobtainedhereare
model. to be considered conservative.
The DM density profile eq. (13) is completely defined Thus, at this point we have completely specified the
for a specific halo by two parameters - scaling radius r DM density distribution in M31 halo. Now we can move
s
and scaling density ρ . Thus, we needed to estimate to the next step - specifying the magnetic field distribu-
s
them for M31. These two parameters are unambigu- tion.
ously linked to another two halo parameters, which are
more meaningful - the halo mass M and its concentra-
∆
tion c . ∆ here denotes an overdensity: halo mass M C. Magnetic field
∆ ∆
by definition means the DM halo mass enclosed in the
sphere, average density inside which is equal to ∆ρ The emission due to DM annihilation, which we hope
crit
with ρ = 3H02 =9.5 10−30 gcm−3 being the current to detect, is generated through the synchrotron mecha-
criticalcrditensit8yπGof the U·niverse.· In our article we used nism. For synchrotron emission modeling it’s crucially
∆ = 100 due to some practical circumstances. How- important to know the magnetic field strength distri-
ever, the parameters of the actual DM density distribu- bution in our emitting volume. A global axis symme-
tion(13)donotdependindeedona∆choice. Halomass try of the large spiral galaxy M31 naturally suggests as-
and concentration are connected with the NFW profile sumption of two-parametric magnetic field distribution
parameters through the following relations: B(ρ′,z′), where B depends on the radial distance from
the centerρ′ inthe galacticplaneandthe verticalheight
1/3 r∆ z′ above the galactic plane.
r = 3M∆ , c r∆, M = ρ (r)4πr2dr, As for B dependence on the vertical coordinate z′,
∆ (cid:18)4π∆ρcrit(cid:19) ∆ ≡ rs ∆ Z0 DM eitx’sp(plazu′s/ibzl′e) twoitahsssuommee esxcaploenehnetigiahltdze′p.enSduecnhcedeBpen∼-
(14) −| | 0 0
dence is commonly used for large spiral galaxies - see
As for the determination of the relevant parameters for
e.g. [18]. In order to consider a reasonable radial de-
M31, we used the results of M and c from [6] and
100 100
pendence of B in the galactic plane we used the findings
[16]. While these are in a good agreement with each
other in determination of M ((1.2 0.3) 1012M vs. of [19],[20]. In [19] large-scale magnetic field in the disk
100 ⊙
(0.91 0.16) 1012M respectively),±they sh·ow a rather planewasmeasuredbetweenradialdistancesof6and14
⊙
± · kpc. These measurements yielded almost constant mag-
bigdiscrepancyinthe estimationofc (the mostprob-
100
netic field over the whole mentioned annular region of
ablevaluescitedbeing 12and 28respectively). The
≈ ≈ 7µG,whichistheexpectedvalueforgalacticdiscs(see
expected radiation fluxes from these two concentration
≈
e.g. [21]). As for the magnetic field strengthin the inner
values differ by about one order of magnitude. This is
region ρ′ < 6 kpc, we were unable to find reliable data.
whywedecidedtotreatthesetwocasesseparately. They
Somemeasurementsforthisregionwerereportedonlyin
canbeconsideredastwolimitingcasesyieldingthemost
[22] and only for two trial locations with poor justifica-
conservativeconstraintsfor c =12 andthe mostopti-
100
tion. Besides this information on the central field prop-
mistic constraints for c =28.
100
erties we can reliably assume probably only one thing -
AtthispointwewerereadytospecifytheNFWprofile
that the field grows towards the center. In such situa-
by estimation of its two parameters ρ ,r . For the first
s s
casec =12we tookM =1.2 1012M from[6]and tion we made a decision to introduce the following field
100 100 ⊙
· distribution (including the vertical dependence part):
obtained ρ ,r by formulas (14). For the second case
s s
c =28 we did not need to calculate ρ ,r - they were
100 s s ρ′ z′
taken directly from the table 2 in [16] as the best-fitting B(ρ′,z′)= B +B exp exp | | .
10 s −ρ′ − z′
values. (cid:18) (cid:18) 0(cid:19)(cid:19) (cid:18) 0 (cid:19)
(15)
Another potentially relevant question in this subsec-
Here B = 7 µG - the plateau value observed on
tion is substructures contribution in our flux. It is well 10
knownthatanyDMhalocontainsalotofsmallsubhalos ρ′ ≈ 6 − 14 kpc. The term Bsexp −ρρ′0′ presents
inside it (see e.g. [1]). These subhalos are very dense an exponentially decreasing central spik(cid:16)e with(cid:17)the char-
6
acteristic radial extent ρ′. Of course, such distribu-
0
tion has a non–physical plateau extending to infinity
B(ρ′ = ,z′ = 0) = 7 µG. But this plateau does not
∞
matter for our calculations: at distances of tens of kpc
DM density becomes very small, and a contribution of
these regions to the overall flux is tiny. The distribu-
tion in eq. (15) can be completely specified by two
free parameters: the vertical scale height z′ and the
0
central field value B(0,0) = B + B . After speci-
10 s
fying these two parameters, the radial scale length ρ′
0
becomes automatically set by smooth connection of the
central exponential spike B exp ρ′ to the plateau
s −ρ′0
B(ρ′ >6 kpc,z′ =0) 7 µG. (cid:16) (cid:17)
Now let’s discuss th≈e specification of the values of z′
0
and B(0,0). We estimated the scale height z′ using [23]
0
and[19]. The firstmentionedarticleoutlines suchagen-
eralpropertyofgalacticmagneticfields: the scaleheight
z′ istypicallyabout4timesgreaterthanthescaleheight
0
of a synchrotronemission,which arises from cosmic rays
in a whole galaxy volume. The last mentioned syn-
Figure3. Densityplotofthemagneticfielddistribution(15),
chrotronscaleheightwasreportedin[19]forM31witha
which was chosen for our computations. Parameter values
kvaplcu.eM≈o0r.e3ovkeprc,.sTuchhisawvoaulludeyfioerldz0′z0′is≈co0n.3firkmpced·4b=y t1h.2e hfieerldesctorrernegstphonBd(0t,o0)th=e5m0oµsGt.pFroorbambolerecdaesteaiwlsitsheetshuebsceecnttiroanl
authors of [19] on the basis of other independent con- IIC. Color version of our article is available online.
siderations: they conclude that the minimal expected
z′ value is not less than 1 kpc. We therefore chose
z0′ = 1.2 kpc for all our c≈alculations. This estimate re- model parameter. These alternative less realistic values
0
flectsaconservativechoice,sincelowervaluesseemtobe werechosentobeB(0,0)=15µGandB(0,0)=300µG.
unrealistic and possible greater values can only increase Theycanbeconsideredasapproximateboundariesofthe
the final expected flux and, hence, strengthen the final intervalwheretheactualfieldvaluelies. Thedetaileddis-
constraints. cussion of the final results variation due to uncertainties
As for the second necessary parameter - B(0,0), the in B(0,0) will be done in the subsection IIIB.
situation is much less certain. As we already men- At this point we have completely specified magnetic
tioned, we found in the literature only one attempt to field distribution. For illustration purposes we showed
measure the central field values in M31 - [22]. This the density plot of our distribution in eq. (15) on fig. 3
work reports the total field strengths for the two lo- for the case of B(0,0)=50 µG.
cations: B(ρ′ 0.3 kpc,z′ 0) = 15 3 µG and The magnetic field distribution in eq. (15) is written
B(ρ′ 0.9 kpc,≈z′ 0) = 19 ≈3 µG. These±values were inacylindricalcoordinatesystemattachedtoM31. Axes
obtain≈ed assuming≈energy equ±ipartition between cosmic x′y′ lie in the M31 disk plane. Another non trivial step
is to make a transformation between coordinate system
rayparticlesandmagneticfields. Suchassumption,how-
ever,hasbeenshowntoleadtoanunderestimationofthe attached to M31 and "laboratory" system with Oz axis
pointing from the M31 center to the observer. This is
magnetic field in similar setting. For example, the au-
necessaryforthe fluxcomputation,becausewearegoing
thors of [20] explored the central field properties of the
toperformanumericalintegrationovertheemittingvol-
MW and derived the stringent lower limit of 50 µG on
the field strengthbasedonmulti-wavelengthstudies of a ume along our line of sight, which is coincident with Oz
cosmic ray emission. At the same time, they mentioned axis and not coincident with Oz′. The two coordinate
that a simple equipartition assumption implies a central systems are shown on fig. 4. We performed the coor-
field strength of only 10 µG. As for the upper bound dinate transformation (ρ,ϕ,z) (ρ′,ϕ′,z′) using the
on a central field stren∼gth, [20] demonstrates that there formulas below, which were deri←ve→d in [24] (formula (6))
isnocertainlyknownvalueforit. Valuesoffewhundreds for cartesian coordinates (x,y,z) (x′,y′,z′):
←→
µG are mentioned to be absolutely possible.
x′ = xsinP +ycosP,
ConsideringthepossibleunderestimationofM31mag-
−
netic field due to the equipartition assumption ([22]), y′ = xcosPcosi ysinPcosi zsini, (16)
− − −
and the similarity between MW and M31, we chose z′ = xcosP sini ysinP sini+zcosi,
− −
B(0,0) = 50 µG as the most probable value. However,
weperformedallourcalculationsforalternativefieldval- where P and i represent the position angle and axis in-
ues as well (which are consideredas less probable)in or- clination respectively, and define the orientation of the
der to study the dependence of our final results on this M31 disc plane with respect to the sky plane. Follow-
7
ing [6], we assumed: P = 38◦,i = 78◦. From eq. (16), nates it’s easy to obtain the one-to-one correspondence
using relations between cylindrical and cartesian coordi- (ρ,ϕ,z) (ρ′,ϕ′,z′):
←→
1
ρ′ = ( ρcosϕsinP +ρsinϕcosP)2+(ρcosϕcosPcosi+ρsinϕsinPcosi+zsini)2 2 ,
− (17)
(z′ = ρcosϕcosPsini ρsinϕsinPsini+zcosi.
−(cid:0) − (cid:1)
Finally, we substituted eq. (17) into eq. (15) obtaining wegaininoverallreliabilityofcalculations. Also,wecan
therefore B as a function of (ρ,ϕ,z), which allowed to not increase α any higher than 5′, because it will lead
∼
compute the final flux. to non-desirable capturing of projected point sources.
As an intermediate step of calculations, we decided to
study also the contribution of leptons with different en-
D. Estimated fluxes
ergies into the total flux. For this purpose we computed
the part of the total flux produced by the leptons with
A final expected radio flux due to DM annihilation in initial energies at production below different thresholds
M31 halo is defined by eq. (1). In order to compute it E . We spanned the whole possible range of E val-
m m
we combined all relevant relations (2)-(4), (6)-(15), (17) ues from the lepton’s rest energy to the WIMP’s rest
and substituted them in eq. (1). The integration limits energy. The results are presented on fig. 7 for all four
forspatialcoordinatescanbe easilyfiguredoutfromthe frequencies used. These plots are generatedfor the most
fig. 4: relevant magnetic field model, WIMP mass m = 100
χ
GeV (relevance of this mass scale will be seen in the re-
06ρ6(d z)tanα, sults section) and the c =12 halo model. Essentially,
− 100
06ϕ62π, (18) we can see on fig. 7 that the main contribution to the

 r100 6z 6r100, total flux (& 50%) for the all frequencies and annihila-
− tion channels comes from a relatively narrow window of
whereαistheangularradiusofourROI.Fig. 5presents initialleptonenergiesbetween 0.01mχc2 and0.1mχc2,
∼
or1GeVand10GeV.Thisenergywindowofoneorderof
thedependenceofthetotalexpectedradiofluxatν =74
MHz on the angular radius of ROI for both selected an- magnitudewidthisconsideredtobenarrowwithrespect
to the whole energy range of produced leptons between
nihilation channels and three trial WIMP masses: m =
χ
10,100and1000GeV.Boundariesofeveryshadedregion 10−5mχc2 and1mχc2,whichspansmorethan5orders
∼
ofmagnitude. Thisimpliesthatthemostimportantfrac-
correspond to the two limiting cases of DM density dis-
tribution in the halo,which were discussed in subsection tion of the lepton population for our considerations has
anenergyabove 1 GeV.This fact, inturn, meansthat
IIB: c =12and c =28. Thus the actualflux is ex-
100 100 ∼
pectedtoliesomewhereinsideshadedregions. Thethick according to fig. 2 the main energy loss mechanisms for
suchpopulationwouldbesynchrotronandICSemission.
central lines present the algebraic averages between the
corresponding boundary curves (fluxes). The flux satu- Thus, these kinds of losses dominate in comparisonwith
rates around α 5′ and does not grow significantly at the other two - bremsstrahlungand Coloumblosses,and
larger radii. We≈also computed the expected signal-to- play a main role for final results. Meanwhile, it’s also
noise ratio (SNR) for the same set of parameters and useful to note that in the case of τ+τ− channel the to-
presenteditonfig. 6. Detailsaboutnoiselevelcomputa- talfluxisconstitutedbymoreenergeticleptonsthanthe
tionscanbefoundinsubsectionIIIA.Theseplotsjustify flux for b¯b channel, which is in complete agreement with
the difference between these two channels mentioned in
our choice of α. Basically,in this choice we are governed
by two counteracting effects. On one hand, fig. 6 clearly the Introduction- τ+τ− channelproduces harder lepton
energy spectrum.
suggeststousethesmallestαpossibleinordertoachieve
thebestSNR.Ontheotherhand,asdiscussedinsubsec- Another important result of our computations is the
tions 2.2-2.3, the accuracy of DM density and magnetic frequencydependenceofradiofluxfromM31halo. These
field strength estimates decreases when approaching the dependencies are shown on fig. 8 for three different
central point of the halo. As a matter of fact, we know WIMP masses, both annihilation channels, both limit-
almost nothing about real distributions of these quanti- ing halo models and the chosen ROI with α = 5′. As a
ties inside the region of size α 0.25′. Hence, our ROI generalfeature,wesee,forallcases,adecreaseofthera-
should be much larger than 0≈.25′ in radius in order to dio flux asthe frequency increases. The steepness ofthis
≈
encompass a volume where our overall knowledge of pa- decrease is, however, more pronounced for lower WIMP
rameters distribution is precise enough. Based on these masses. This behavior suggests to primarily leverage on
considerations, we chose α = 5′ as the optimal ROI ra- lowfrequencyobservationsforobtainingDMconstraints.
dius. With such choice we lose little in SNR, however, However, the sensitivity of radio surveys typically in-
8
Figure 4. Two coordinate systems (x,y,z) and (x′,y′,z′). The plane x′y′ is coincident with the M31 disc plane. Integration
overthe emitting volumeis performed in (x,y,z). More details are in subsection IIC.
mx = 10GeV mx = 10GeV
m = 100GeV
x m =100GeV
x
m =1000GeV
x
m =1000GeV
x
ΧΧ®bb ΧΧ®Τ+Τ-
Figure 5. The expected radio flux dependenceon theangular radius of the ROI.The expected flux lies inside shaded regions.
Upper and lower boundaries of these regions correspond to the different halo concentration values c = 28 and c = 12
100 100
respectively. Thick central lines present algebraic averages between corresponding limiting cases. Curves are shown for the
two annihilation channels χχ → b¯b, τ+τ−, three WIMP masses mχ = 10,100, 1000 GeV and the most probable magnetic
field distribution with the central value of 50 µG. The annihilation cross section is hσvi = 3·10−26 cm3/s, the observational
frequency is ν = 74 MHz. Formore details see subsection IID.
9
Figure 6. The expected SNR dependence on the angular radius of ROI. The expected SNR lies inside shaded regions. Upper
and lower boundaries of these regions correspond to different halo concentration values c = 28 and c = 12 respectively.
100 100
Curves are shown for the two annihilation channels χχ → b¯b, τ+τ−, three WIMP masses mχ = 10,100, 1000 GeV and the
most probable magnetic field distribution with the central value of 50 µG. The annihilation cross section is hσvi = 3·10−26
cm3/s, the observational frequency is ν = 74 MHz. Formore details see subsection IID.
ΧΧ®bb ΧΧ®Τ+Τ-
c =12 c =12
100 100
74MHz 74MHz
325MHz 325MHz
1400MHz 1400MHz
4850MHz 4850MHz
Figure 7. The part of total radio flux produced by leptons with energies at production below different thresholds. These
illustrationsaregeneratedforthefollowingsetofparameters: c =12,hσvi=3·10−26 cm3/s,B(0,0)=50µG.WIMPmass
100
is chosen to be mχ =100 GeV, which reflects themost relevant mass scale with respect to final results, as will be seen below.
We can note that the majority of the total flux are produced by leptons with initial energies above ∼0.01mχc2. More details
are in subsection IID.
10
9. Firstly, we notice that our ROI does not capture any
Table I.Main parameters of all radio surveysused.
contaminating point sources. On the last three images
(325, 1400 and 4850 MHz) some signal from M31 center
Beam
RMS noise is clearly visible. The 74 MHz image presents essentially
diameter
Survey Frequency, or FWHM level inside noise. For the precise flux upper limit derivations we
MHz the beam
2βi, σ(b), mJy neededto obtainthe measuredsignalvaluesci (i=1−4
arcseconds i denotes different observational frequencies) and the rms
(root of mean square or root of dispersion) noise values
σ forourROIforallimages(weassumeGaussiandistri-
i
VLSS 74 80 100 butions of noise levels inside a radio telescope beam and
inside the ROI). The measured signal values were easily
read by the Aladin. As for the rms noise values, we de-
WENSS 325 82 4.0
rived from basic principles that they can be estimated
as σ = σ(b)√N = σ(b) α, where N are the numbers of
NVSS 1400 46 0.5 i i i i βi i
beams contained in the ROI, β are the angular radii of
i
beams (or half of the FWHM for each survey) and σ(b)
i
GB6 4850 240 4.0 arethermsnoisesinsideabeamforeachsurvey(citedin
table I). Thus, we have obtained both necessary ingredi-
ents for a flux upper limit estimationfromthe ROI - the
measured signal values and the rms noise values.
creaseswith frequency, partiallycompensating the effect Having all the necessary information, we obtained the
of the lower signal. There is, therefore, no obvious op- limiting cross section values for cases of different param-
timal frequency to be used and that would lead to the eter sets using the following formula for the probability
strongest constraints possible. For this reason, we opted density of the noise values inside our ROI:
radioobservationsatverydifferentfrequenciesforactual
constraintsderivation. Thisprocedureanditsresultswill 1 (c w ( σv ))2
i i
p (n =c w ( σv ))= exp − h i ,
be described in the next section. i i i− i h i √2πσ − 2σ2
i (cid:18) i (cid:19)
(19)
where n - the noise values inside the ROI; c - the mea-
i i
III. OBTAINING CONSTRAINTS ON DM sured signal values from the ROI and w ( σv ) - the ex-
i
ANNIHILATION h i
pectedsignalsfromWIMPannihilation,whichislinearly
proportional to an unknown annihilation cross section
A. Comparison with radio surveys σv . Here we do not include any background radiation
h i
(this question will be further developed in the next sec-
In order to make a next step on the way to final con- tion). For our final constraints onto σv we decided to
h i
straints, we needed to obtain the upper limits on the use 99.73% confedence level, which corresponds to 3σ
actually observed radio fluxes from our ROI. Compari- Gaussian confedence level.
son of these upper limits with the expected fluxes com- In order to guide intuition on which experiment is the
putedaboveimmediatelyprovidestheupperlimitsonthe mostconstrainingfor a givenparticle mass,infig. 10 we
annihilation cross section for different values of WIMP showed our results on the cross section when each fre-
masses,since the expected fluxes is directly proportional quency is considered individually. There we plotted the
to σv . For this purpose we studied many major radio constraints obtained at 99.7% confidence level by four
h i
surveys conducted in the past and chose four of them different surveys used, for two halo models (described
(corresponding frequencies and detailed descriptions are in subsection IIB), two annihilation channels selected
provided in brackets): VLSS (74 MHz, [25]), WENSS χχ b¯b,τ+τ− andthemostprobablemagneticfielddis-
→
(325MHz, [26]), NVSS (1400MHz, [27]) andGB6 (4850 tributionwithB(0,0)=50µG.Wespannedtherangeof
MHz, [28]). All other surveys were rejected due to ei- WIMP masses between 6 GeV and 1000 GeV and com-
ther one of the following two reasons: too low resolution puted the limiting σv values for a subset of WIMP
h i
of a survey (i.e. the telescope beam size is larger than masses inside this range, interpolating results for mass
the ROI size) or not covering M31 area of the sky. In values in between. The general key properties of our
table I we outlined the main surveys’ parameters, which results can be summarized as follows: first of all, differ-
were collected from [25]-[28] and then used in calcula- ent halo models yield very different limiting annihilation
tions. The uncertaintiesofsurveyfrequencies,causedby cross sections (by about one order of magnitude) for all
finite bandwidth, do not exceed the level of 4%, and WIMP masses. This reflects the high sensitivity of our
∼
do not affect final results significantly. results to the halo model assumed. Another noticeable
Real sky images were viewed and analysed by the Al- factisthatτ+τ− annihilationchannelappearsto be less
adin software, which is described in [29]. The images of constrainedthanb¯bchannel,incompleteagreementwith
our ROI from all four radio surveys are shown on fig. the mentioned expectation that harder secondary spec-