Two Emission Mechanisms in the Fermi Bubbles: A Possible Signal of Annihilating Dark Matter

We study the variation of the spectrum of the Fermi Bubbles with Galactic latitude. Far from the Galactic plane (|b|>30 degrees), the observed gamma-ray emission is nearly invariant with latitude, and is consistent with arising from inverse Compton scattering of the interstellar radiation field by cosmic-ray electrons with an approximately power-law spectrum. The same electrons in the presence of microgauss-scale magnetic fields can also generate the the observed microwave"haze". At lower latitudes (b<20 degrees), in contrast, the spectrum of the emission correlated with the Bubbles possesses a pronounced spectral feature peaking at 1-4 GeV (in E^2 dN/dE) which cannot be generated by any realistic spectrum of electrons. Instead, we conclude that a second (non-inverse-Compton) emission mechanism must be responsible for the bulk of the low-energy, low-latitude emission. This second component is spectrally similar to the excess GeV emission previously reported from the Galactic Center (GC), and also appears spatially consistent with a luminosity per volume falling approximately as r^-2.4, where r is the distance from the GC. We argue that the spectral feature visible in the low-latitude Bubbles is the extended counterpart of the GC excess, now detected out to at least 2-3 kpc from the GC. The spectrum and angular distribution of the signal is consistent with that predicted from ~10 GeV dark matter particles annihilating to leptons, or from ~50 GeV dark matter particles annihilating to quarks, following a distribution similar to the canonical Navarro-Frenk-White (NFW) profile. We also consider millisecond pulsars as a possible astrophysical explanation for the signal, as observed millisecond pulsars possess a spectral cutoff at approximately the required energy. Any such scenario would require a large population of unresolved millisecond pulsars extending at least 2-3 kpc from the GC.


I. INTRODUCTION
Data from the Fermi Gamma-Ray Space Telescope have revealed a pair of large gamma-ray lobes extending approximately 50 • north and south of the Galactic Center [1]. These lobes, known as the Fermi Bubbles, are visible in gamma-rays between ∼1-100 GeV and have a markedly harder spectrum (dN/dE ∝ E −2 ) than the gamma-ray emission associated with the Galactic Disk.
The Bubbles were originally studied as a possible gamma-ray counterpart to the WMAP haze [2], a spectrally hard microwave excess in the inner Galaxy most clearly visible in WMAP 's 23 and 33 GHz frequency bands. The haze was first discovered in 2003 [3], and has been studied over the past decade as a possible signature of a new hard electron population in the inner Galaxy [4][5][6], producing microwave synchrotron radiation in the Galactic magnetic field. Recently, the existence of the microwave haze has been confirmed by the Planck experiment [7], whose data indicate a strong degree of spatial coincidence between the microwave haze and the gammaray Bubbles, further supporting the hypothesis that these signals have a common origin. Perhaps the simplest possibility is that the gamma rays arise from inverse Compton scattering (ICS) by the same hard electron population that produces the haze via synchrotron.
The question of the origin and nature of the Fermi Bubbles has been the subject of much debate. One key question is whether these gamma-rays are produced by a hadronic [8][9][10] or leptonic [11][12][13][14][15] mechanism, i.e. whether they arise from the scattering of energetic protons on the gas of the interstellar medium, or from the ICS of photons from the interstellar radiation field by energetic electrons. An example of a hadronic scenario was proposed by Aharonian and Crocker [8], in which the Bubbles are billion-year-old reservoirs of energetic protons, which were injected as a result of star formation in the Galactic Center and are confined within the Bubbles by magnetic fields. Leptonic scenarios have garnered somewhat more attention in the literature, and provide a straightforward link with the spatially correlated emission observed in the microwave and radio [16] (hadronic scenarios will also generate synchrotron emission through the electrons produced in charged pion decays, but since these electrons will diffuse after being produced, the connection between their spectrum and spatial distribution and that of the gamma-rays is not as straightforward as in the leptonic scenario). In such leptonic models, electrons are accelerated by a shock or a series of shocks and/or by Fermi acceleration in turbulent magnetic fields behind the shock front [17]; the shock(s) may be fueled by accretion onto the supermassive black hole at the Galactic arXiv:1302.6589v1 [astro-ph.HE] 26 Feb 2013 Center, by starburst activity, or by some other mechanism.
In this article, we reexamine the Fermi Bubbles and the variation of their spectrum with Galactic latitude. Far from the Galactic plane (|b| > ∼ 30 • ), the observed gamma-ray spectrum is nearly invariant with latitude and fairly flat over the energy range observed by Fermi. This spectrum can be well explained by inverse Compton scattering of cosmic microwave background (CMB), infrared, and starlight photons by a population of GeV-TeV electrons with an approximately power-law spectrum (dN e /dE e ∼ E −3 e ). Furthermore, we find that this same population of cosmic ray electrons leads to synchrotron emission of the same amplitude as the observed microwave haze, if microgauss-scale magnetic fields are present in the high-latitude regions of the Bubbles. The success of this simple and self-consistent picture provides strong support for a leptonic origin of the high-latitude emission from the Fermi Bubbles.
At latitudes closer to the disk, however, a leptonic origin of all the emission associated with the Bubbles does not appear possible. The gamma-ray spectrum of the Fermi Bubbles at latitudes within approximately 20 • of the Galactic plane possesses a peak at energies of a few GeV, and cannot be generated by inverse Compton scattering of starlight, infrared, or cosmic microwave background radiation by any realistic steady-state electron population. Furthermore, no realistic spectrum of cosmic ray protons is capable of accounting for the gamma-ray spectrum observed at these low latitudes.
Gamma-ray emission with a similar spectrum has been previously identified from the region surrounding the Galactic Center (GC) [18][19][20][21][22]. Proposed origins for this excess include annihilating dark matter [18][19][20]22], a population of millisecond pulsars [18,19,[22][23][24], or cosmic ray interactions with gas [18,19,22,25,26]. In this paper, we show that the non-inverse Compton component of the emission from the Fermi Bubbles identified in this study is spectrally and morphologically consistent with being the extended counterpart of this GC excess, revealing that this emission is not confined to the GC, but extends out to at least ∼2-3 kiloparsecs from the Galactic plane. The morphology of this signal is consistent with originating from annihilating dark matter distributed according to a generalized Navarro-Frenk-White (NFW) profile with an inner slope of ρ ∝ r −1.2 , where r is the distance to the GC (the GC excess, in isolation, favors a power-law slope in the range of 1.2-1.4 [18,19,22]). The spectral shape of this signal can be accommodated by dark matter particles of mass ∼10 GeV annihilating to leptons, or by ∼50 GeV particles annihilating to quarks. In either case, the normalization of the observed signal requires an annihilation cross section on the order of σv ∼ (6 − 8) × 10 −27 cm 3 /s, similar to that expected of a thermal relic of the Big Bang.
The remainder of this paper is structured as follows. In Sec. II we describe the analysis used to extract the spectra from various regions of the Fermi Bubbles. In Sec. III we demonstrate that the high-latitude regions of the Bubbles can be accounted for with an approximately power-law spectrum of GeV-TeV electrons, which can simultaneously produce the observed microwave haze as synchrotron. In Sec. IV we turn our attention to the low-latitude regions of the Fermi Bubbles, and show that their spectra cannot be accounted for by inverse Compton scattering. Instead, an additional mechanism is required, capable of producing a spectrum which peaks strongly at energies of a few GeV. In Sec. V we compare this signal with that previously observed from the Galactic Center. In Sec. VI, we study the residuals remaining when the (best-fit) known backgrounds are subtracted, and demonstrate how the few-GeV peak appears in these residuals. In Sec. VII, we show that the emission observed at low-latitudes is better described by a spherically symmetric, NFW-like morphology than by a flat-brightness distribution confined to the regions of the Fermi Bubbles. In Sec. VIII we discuss possible interpretations of this signal, including annihilating dark matter and a population of gamma-ray pulsars. We summarize our results and draw conclusions in Sec. IX. This paper also includes five appendices, which describe various cross checks of our results and other supplementary material.

II. EXTRACTING THE LATITUDE-DEPENDENT SPECTRUM OF THE GAMMA-RAY BUBBLES
In our analysis, we employ the publicly available Pass 7 Version 6 data release from Fermi, with 4.5 years of photon data. 1 We apply a standard zenith angle cut to exclude the Earth limb, rejecting events with zenith angles greater than 100 • . We also employ the recommended diffuse analysis cuts on the data quality, nominal science configuration and rocking angle: DATA QUAL==1, LAT CONFIG==1, ABS(ROCK ANGLE)<52. Throughout, we use the class of events designated ULTRACLEAN, but have confirmed that employing the SOURCE or CLEAN classes does not alter our conclusions.
We generate skymaps for 30 log-spaced energy bins spanning the range from 0.3 GeV to 300 GeV, binning the photons on an NSIDE=256 HEALPix grid, and smooth all maps to 2 • FWHM (full width at half maximum). As in Ref. [1], we use front-converting events only (which have better inherent angular resolution) at energies below 1 GeV. We follow the prescription from Ref. [1] for point source subtraction, using the Fermi 1-year source catalog and masking the very brightest sources.
For each energy bin, we fit the skymap as a linear combination of templates, maximizing the Poisson likelihood. The Gaussian errors on the fit coefficients are computed from the likelihood by ∆ ln L = 1/2, and do not take into account the systematic error in the event that the linear combination of templates is a poor description of the data. Further details of the fitting procedure may be found in Ref. [1] and in Appendix B of Ref. [2]. We employ several different template combinations to test the robustness of our results to the foreground model.
In the Galactic disk, there is a substantial population of unsubtracted point sources, as well as bright diffuse emission; consequently, we mask the inner disk. Throughout our study, we will show results for masks with |b| < 1 • , |b| < 2 • and |b| < 5 • , to test the dependence of our results on this parameter.
To determine the spectrum of the Fermi Bubbles as a function of latitude, we divide the standard spatial template for the Bubbles (as defined in Ref. [1]) into subregions by (absolute) latitude: |b| < 10 • , 10 • < |b| < 20 • , 20 • < |b| < 30 • , 30 • < |b| < 40 • , and 40 • < |b| < 50 • (see Fig. 1). We smooth all templates to the scale of the maps. We fit separately for the spectrum in each of these latitude bands, varying the degree of masking of the Galactic Disk, and with a range of template-based models for the known backgrounds. All fits include an isotropic offset, to absorb residual cosmic-ray contamination and isotropic diffuse emission. The two primary possibilities we consider for the additional templates are: • Diffuse model : We take the Fermi diffuse model from version P6V11 of the Fermi Science Tools, smooth it to match the maps we are using, interpolate to the appropriate energies, and perform the fit using only this template (in addition to the latitude-sliced Bubbles templates and the isotropic template, which are universal to all subtraction methods). This version of the diffuse model has been adjusted to fit the data assuming no contribution from the Bubbles; consequently it may absorb some of the Bubbles-correlated emission at the cost of oversubtraction in neighboring regions. It was also designed primarily to model the emission at energies 50 GeV, and is not recommended for use at very low latitudes, |b| < 3 • . However, since our signal extends to quite high latitudes and the energies of greatest interest are at 50 GeV, these latter caveats do not pose severe problems for our study.
• Low-energy template: We employ the Schlegel-Finkbeiner-Davis (SFD) dust map [27] as a template for emission from cosmic-ray protons scattering on the gas (see Refs. [1,2] for a discussion). We take the Fermi data at 0.5-1.0 GeV (where the Bubbles are less pronounced [1]) and subtract the SFD dust map to obtain an approximate template for emission from inverse Compton scattering by cosmic-ray electrons; this is the dominant contribution to the diffuse background after the dust/gas-correlated emission has been removed. We then fit the higher-energy data using this template, the SFD dust map, and a flat template for the large soft-spectrum structure known as Loop I. This method avoids the use of complicated models and minimizes the use of external maps, but by construction cannot probe the Bubbles spectrum at energies around or below 1 GeV, and does not take into account spectral variation in the various emission components with position in the Galaxy.
Each of these templates has been discussed in greater depth in Ref. [1]; we refer the reader to that work for further details. We also employ the same normalization convention as Ref. [1]; the coefficients of the SFD dust map, 0.5-1.0 GeV map, and Fermi diffuse model are multiplied by the average value of these templates within the entire region defined by the Bubbles (and outside of the mask). The other templates are flat (in projected intensity) within the regions that they are non-zero (over a given latitude range within the Bubbles, for example).
We employ the "diffuse model" fit for our primary results, and use the "low-energy template" fit as a check for possible systematic errors introduced by our use of the Fermi diffuse model. We have also tested the "simple disk" model as employed in Ref. [1], where the ICS emission is described by a simple geometric diffused disk template. As this template proved inadequate for modelling the data close to the Galactic plane, and the results obtained using it were found to depend strongly on the degree of masking of the disk, we have relegated discussion of this model to Appendix A.
In Fig. 2 we show the extracted spectrum for each of the fitted templates, masking only the region within one  Fig. 1), for the two foreground models we employ (see text). The Galactic Disk is masked for |b| < 1 • in each case. The left and center panels employ the "diffuse model" fit, for the entire sky in the left panel and the southern hemisphere in the center panel. The right panel employs the "low-energy template" fit over the entire sky (see text for the details of the fitting procedures).
degree of the Galactic plane, |b| < 1 • . We show results found using the "diffuse model" and the "low-energy template". In the center frame, we show the fit restricting to the southern sky (b < 0), which we might expect to be less contaminated by bright features such as Loop I. As expected, the error bars are larger in this case due to lower statistics, but the results are not otherwise significantly altered. In Fig. 3, we show the spectra extracted for the gamma-ray Bubbles, and the dependence on the degree of masking of the disk, in each range of Galactic latitude. For our two preferred template models, the results are largely stable to changes in the mask.
While the gamma-ray spectra extracted using the lowenergy template appear somewhat different from those derived using the diffuse model, this is natural and expected, particularly at low energies, since part of the emission associated with the Bubbles is included in the low-energy template itself. Additionally, at high energies and low latitudes the low-energy template fit yields a significant amount of emission roughly flat in E 2 dN/dE, which is nearly absent in the diffuse model fit; these issues are discussed in further detail in Appendix B.
We note that the spectrum is almost invariant from |b| = 20 • − 50 • . This suggests that the electrons responsible for the observed emission in any leptonic scenario must either be accelerated in situ or instead travel from the inner Galaxy very rapidly, avoiding significant energy losses (the distance over which TeV electrons propagate via standard diffusion without significant energy losses is considerably less than the 5 or more kpc to which this angular range corresponds). In contrast, a pronounced change in the Bubbles' spectrum is observed at lower latitudes. In an attempt to quantify the significance of this transition, we have compared the quality of the fit found using five separate latitude-sliced Bubbles templates to that found using only a single Bubbles template. Even conservatively limiting our analysis to the cleaner south-ern bubble, and masking within 5 • of the disk, we find that the five-Bubbles-templates model is favored over the single Bubbles template at the level of approximately 16σ. However, it is important to note that this is a formal significance, accounting only for statistical error; there is a degree of unavoidable and unaccounted-for systematic error in that neither model is a "good fit", in the sense of describing the sky to the level of Poisson noise.

III. COSMIC RAY ELECTRONS AS THE SOURCE OF THE HIGH-LATITUDE GAMMA-RAY BUBBLES AND SYNCHROTRON HAZE
Following Blumenthal and Gould [28], we employ the full Klein-Nishina formula to compute the spectrum of inverse Compton emission from an arbitrary electron population. For the problem at hand, we need to consider scattering with the CMB as well as with starlight and infrared radiation. In our calculations, we adopt the interstellar radiation model of Ref. [29]. At energies below ∼3 × 10 −3 eV, the CMB dominates the energy density, while starlight is important at higher energies.
The gamma-ray spectra observed from various regions of the Fermi Bubbles are shown again in Fig. 4 (as found using the diffuse model template fit, and masking within 1 degree of the disk). To determine whether these gamma-rays could be the product of inverse Compton scattering, we take an arbitrary (binned) spectrum of electrons and compare the resulting inverse Compton emission to that shown in Fig. 4. In Fig. 5, we plot the electron spectrum which provides the best possible fit to the gamma-ray spectrum for each latitude range (and error bars around the best fit). The solid line in each frame of Fig. 4   The spectrum extracted for the gamma-ray Bubbles in ten-degree latitude bands: in order from the top row, The left and center panels use the "diffuse model" template fit (see text); in the center panels, the fit is restricted to b < 0 in addition to the masking. The right panels use the "low-energy template" approach (see text). The different colors show different choices for the latitude cut to remove the Galactic Disk: Where the 1σ error bars overlap with zero, we instead plot downward-pointing arrows corresponding to the 3σ upper limits on the emission.

FIG. 4:
The gamma-ray spectrum of the Fermi Bubbles, broken into different regions by Galactic latitude (see Fig. 1). The solid lines denote the best-fit spectrum of inverse Compton emission, as calculated from the central values of the electron spectra shown in Fig. 5. At high latitudes, the spectrum of the Fermi Bubbles is consistent with originating entirely from the inverse Compton scattering of GeV-TeV electrons, while at lower latitudes inverse Compton scattering alone cannot account for the observed emission. The dashed line in each frame denotes the spectrum of inverse Compton scattering that would be predicted from a spectrum of electrons the same as that required to generate the inverse Compton scattering spectrum observed in the highest latitude region (|b| = 40 • − 50 • ), as discussed in Sec. IV.
FIG. 5: The spectrum of electrons within the volume of the Fermi Bubbles that is best able to fit the observed gamma-ray spectrum (see Fig. 4) through inverse Compton scattering alone. At high latitudes, the electron spectrum yields a good fit to the observed gamma-ray spectrum, while at low latitudes no spectrum of electrons is able to produce a gamma-ray spectrum consistent with observations. power-law spectrum of electrons (dN e /dE e ∝ E −3 e ) between GeV and TeV energies can produce a spectrum of gamma-rays consistent with that observed from the Fermi Bubbles. At |b| = 40 • − 50 • , for example, the best-fit electron spectrum shown in Fig. 5 provides an excellent fit (χ 2 = 21.8 over 29-7 degrees-of-freedom) to the observed gamma-ray spectrum.
As noted in the Introduction, the inverse Compton interpretation of the Fermi Bubbles is supported by the observation of spatially correlated emission at both microwave [7] and radio [16] wavelengths. For a given magnetic field strength (or distribution of magnetic field   . Red data points indicate the spectrum of the haze in the southern sky within each given latitude range, as derived from Ref. [6]; the average over l is performed for −5 • < l < 15 • , as in Fig. 5 of Ref. [6]. At somewhat lower frequencies, the green triangles and dashed lines connecting them denote the spectral index (not amplitude) for polarized emission between 2.3 GHz and 23 GHz, as extracted from S-PASS [30] data by Ref. [16]; the averaging is performed over −5 • < l < 15 • and restricted to |b| > 15 • to avoid depolarization associated with the Galactic plane. This comparison demonstrates that in the presence of microgauss-scale magnetic fields, the electron population responsible for the observed ICS emission can also account for the observed synchrotron haze. Also shown for comparison as blue data points is the spectrum of the microwave haze as a whole, as given by Ref. [7]. As these results are not binned by latitude, however, they should not be directly compared to the predicted synchrotron spectra. strengths), we can calculate the spectrum of synchrotron emission that is predicted to result from the spectrum of electrons shown in Fig. 5. For simplicity, we take the central values for the extracted spectrum and assume an isotropic electron population; we do not extrapolate to energies above the highest bin. In Fig. 6, we show the resulting synchrotron spectra for a range of magnetic fields from 0.1-100 µG, assumed to be uniform throughout the region of the Bubbles template in question. For comparison, the magnetic field strength in the local region of the Milky Way is thought to be of order a few µG [31]. One should keep in mind, however, that magnetic fields in localized regions and filaments may be much higher than the spatially averaged value.
To compare the predicted synchrotron spectrum to observations, we employ the WMAP 7 results of Ref. [6], and the spectral index of the polarized emission between 2.3 GHz and 23 GHz, as obtained from S-PASS [30] data and described in Ref. [16]. 2 We also show the recent results combining data from WMAP and Planck [7], but these results were not binned by latitude so should not be directly compared to the predicted synchrotron curves.
The amplitudes of the gamma-ray Bubbles and microwave haze appear consistent with arising from a common population of cosmic ray electrons for magnetic fields on the order of a few µG. The spectrum deduced from WMAP is somewhat harder than expected from the corresponding gamma-ray data, but this may be the result of contamination by the CMB (which would be expected to produce the appearance of spectral hardening). In contrast, the spectral index for the polarized emission found by the authors of Ref. [16] appears consistent with expectations from the gamma-rays. If the spectral hardness observed in WMAP data is confirmed by Planck, however, it could suggest a scenario in which most of the synchrotron emission is produced in regions with higher than average magnetic fields. For example, if magnetic fields as high as ∼ 30µG are present in even ∼1% of the volume of the Bubbles, then this discrepancy could be ameliorated. We caution that the success of this picture does not rule out hadronic scenarios; we focus on leptonic scenarios in this work because their consistency can be checked in a straightforward and model-independent way, since the microwaves and gamma-rays originate from the same steady-state electron population. In hadronic scenarios, in contrast, the gamma-rays provide a probe of the proton CR spectrum, but the microwaves probe the electron spectrum after diffusion and cooling, and so the consistency may depend on diffusion parameters in addition to the magnetic field. Consequently, we leave a careful study of consistency in the hadronic case for future work.
In the previous section, we demonstrated that the spectrum observed from the high latitude (|b| > ∼ 30 • ) regions of the Fermi Bubbles can be accounted for by inverse Compton scattering of an approximately power-law spectrum of cosmic ray electrons. The same electrons, in the presence of µG-scale magnetic fields, can also account for the observed amplitude of the WMAP haze. At lower Galactic latitudes, however, we find that inverse Compton scattering alone cannot account for the observed gamma-ray emission. In particular, we find that the gamma-ray spectrum climbs rapidly between 0.3 and 1 GeV at low-latitudes (see Fig. 4), and this rise cannot be accounted for by the inverse Compton scattering of any physically realistic spectrum of electrons. Quantitatively, for the choice of binning used in Fig. 5, we find that an entirely inverse-Compton origin for the gammarays observed from the |b| = 1 • −10 • or |b| = 10 • −20 • regions of the Fermi Bubbles yields best-fits of χ 2 = 100.4 and 95.7, respectively, each over 29-7 degrees-of-freedom. From this, we are forced to conclude that a non-negligible fraction of the gamma-ray emission observed from the low-latitude regions of the Fermi Bubbles is not the result of inverse Compton scattering. 3 If we assume that the spectral shape of electrons present throughout the volume of the Fermi Bubbles does not change significantly with latitude, and that the gamma-ray emission from the highest latitude region (|b| = 40 • − 50 • ) is dominated by the products of inverse Compton scattering, we can calculate the inverse Compton contribution from each of the lower latitude regions. In each frame of Fig. 4, the dashed line denotes this predicted contribution from inverse Compton scattering. Not surprisingly, this component makes up most or all of the total observed emission at high latitudes, but only a small fraction at low latitudes. Note that the FIG. 7: The gamma-ray spectrum of the Fermi Bubbles after subtracting a contribution from inverse Compton emission, derived using the electron spectrum (up to normalization) found in our best-fit to the |b| = 40 • − 50 • region. This illustrates the characteristics of the additional (non-inverse Compton) component of the gamma-ray emission from the Fermi Bubbles, which is quite bright at low Galactic latitudes. We caution that these extracted spectra are subject to a number of systematic uncertainties, such as those associated with the interstellar radiation field model, and due to uncertainties and variations in the electron spectra throughout the volume of the Bubbles. These extracted spectra can, however, be taken as indicative of the broad spectral features of the non-inverse Compton component of the Bubbles emission. Shown as dashed lines is the predicted contribution of gamma-rays from the annihilations of 10 GeV dark matter particles (to τ + τ − ) distributed according to a generalized NFW profile with an inner slope of γ = 1.2, as described in Sec. V. We remind the reader that the backgrounds are largest near the disk and thus there are significant systematic uncertainties in the spectrum from the low latitude (|b| = 1 • −10 • ) region, especially at low energies.
shapes of these dashed lines differ with latitude only as a result of variations in the radiation field model.
To better isolate the spectrum of the non-inverse Compton emission, we subtract this estimate of the inverse Compton contribution from the total emission of the Fermi Bubbles in each latitude range. This residual spectrum is plotted in Fig. 7. While one should keep in mind that the error bars shown in this figure do not take into account any variations in the spectral shape of cosmic ray electrons throughout the volume of the Bubbles, this provides us with what is likely to be a reasonable estimate of the spectrum and intensity of the non-inverse Compton emission exhibited in the low-latitude regions of the Fermi Bubbles.
The spectra shown in Fig. 7 exhibit some rather distinctive features. In particular, this spectral component peaks strongly at energies of ∼1-5 GeV, and has no statistically significant presence above ∼10 GeV. Furthermore, the intensity of this component is a very strong function of Galactic latitude, being more than an order of magnitude brighter at low latitudes than at intermediate latitudes. These spectral and morphological characteristics are quite similar to those exhibited by the gamma-ray emission previously observed from the inner few degrees surrounding the Galactic Center [18][19][20][21][22]. In the following section, we will discuss this comparison in more detail.

V. COMPARISON WITH GAMMA-RAY EMISSION FROM THE GALACTIC CENTER
In previous studies of Fermi data, multiple sets of authors have identified the presence of a bright and spatially extended gamma-ray source around the Galactic Center, peaking at energies of a few GeV [18][19][20][21][22]. In particular, Ref. [19] reports that the morphology of this source implies a luminosity per volume that scales as r −2.4 to r −2.8 , where r is the distance to the Galactic Center. Similar profiles were found to provide good fits in Refs. [18,21,22]. Each of these studies also found that the spectrum of this spatially extended emission peaks strongly at energies of a few GeV, very similar to the peak found in this study in the low-latitude emission from the Fermi Bubbles.
This comparison strongly suggests that the non-inverse Compton emission we have identified in the low-latitude regions of the Fermi Bubbles is in fact the more spatially extended counterpart of the gamma-ray signal previously reported from the innermost few degrees around the Galactic Center. To further explore this comparison, we plot as a dashed line in each frame of Fig. 7 the contribution predicted by an annihilating dark matter model found in Ref. [32] to provide a good fit to the Galactic Center signal. In this model, 10 GeV dark matter particles annihilate to tau lepton pairs (with σv = 2 × 10 −27 cm 3 /s to τ + τ − , or σv = 6 × 10 −27 cm 3 /s if annihilations proceed equally to τ + τ − , µ + µ − and e + e − ), and are distributed following a generalized Navarro-Frenk-White profile, with an inner slope of γ =1.2 (the classic NFW profile corresponds to an inner slope of γ = 1), a scale radius of 20 kiloparsecs, and normalized such that the local dark matter density is 0.4 GeV/cm 3 . Specifically, the functional form is: From this comparison of both the spectrum and morphology of these signals, we conclude that the gammaray component we have identified within the low-latitude regions of the Fermi Bubbles is the more spatially extended continuation of the gamma-ray signal previously observed in the Galactic Center. Furthermore, we confirm that the dark matter models previously shown to be capable of accounting for the gamma-ray emission in the inner Galaxy are also capable of providing an explanation for the non-inverse Compton emission observed in the low-latitude spectrum of the Fermi Bubbles.

VI. MORPHOLOGY OF THE LOW-LATITUDE SIGNAL
One approach to extracting the morphology of the few-GeV spectral feature discussed in the previous sections is to examine the residual sky maps produced when the background templates (for example, the diffuse model template and isotropic background), multiplied by their best-fit coefficients, are subtracted from the data. Equivalently, these maps are produced by taking the residuals of the template fit, and re-adding the latitude-sliced Bubbles templates, weighted by their best-fit coefficients. This largely cancels out structure due to mismatches between the shape of the Bubbles templates and the actual excess, although such mismatches may still bias the fit. (This procedure has been recommended for extracting the spectrum of an excess in Ref. [7].) Fig. 8 shows these "residual" maps (in average E 2 dN/dE) for the "diffuse model" fit, masked within 5 degrees of the Galactic plane, in four energy bands spanning the range from 1 − 20 GeV. The high-latitude Fermi Bubbles are visible in all bins with comparable brightness, as expected since the high-latitude spectrum is approximately flat in E 2 dN/dE. At low latitudes, there is a pronounced excess around the Galactic Center, not disk-like in shape, which is clearly visible in the 1 − 2 and 2 − 5 GeV maps, but rapidly fades at higher energies.
Having generated these residual maps, we can study the emission binned in b, within a narrow longitude range (|l| < 5). As an illustrative exercise, we subtract the 10 − 50 GeV (residual) map from the 1 − 10 GeV (residual) map, where the excess is most pronounced (we ignore energies higher than 50 GeV because the P6V11 Fermi diffuse model was not designed for studies at those energies and, due to low photon statistics, little information would be gained from higher-energy photons in any  case). For more details, see Appendix E. Since the maps are given in E 2 dN/dE, a zero result indicates an average spectrum of dN/dE ∝ E −2 between these two bins. The results for the southern hemisphere, where there are fewer bright sources and local features, are shown in Fig. 9. While we find that at high latitudes the spectrum is consistent with dN/dE ∝ E −2 , the lower latitude emission (|b| < 10 • − 15 • ) reveals significant additional emission at low energies. The error bars shown in this figure, computed from the standard deviation of the pixel values in each ∆b = 2 • bin, are quite large and non-negligibly correlated, but provide a sense of the uncertainty in the rate at which the signal falls off away from the Galactic plane.

VII. FITTING THE SPECTRAL BUMP WITH AN NFW PROFILE
So far, we have focused our attention on the region of the sky occupied by the Bubbles. If annihilating dark matter is responsible for the signal observed at low latitudes, however, this signal should be distributed around the Galactic Center with approximate spherically symmetry, and will not be restricted to the Bubbles. In this section, we include an additional template in our analysis, corresponding to the signal expected from annihilating dark matter distributed according to a somewhat steepened NFW profile (with an inner slope of γ = 1.2).
In Figs. 10 and 11, we show the spectra associated with the NFW template itself, and with the latitude-sliced Bubbles in the presence of the NFW template, respectively. The normalization of the NFW template corresponds to its brightness at 5 • from the Galactic Center.
With the inclusion of this additional NFW template, we find that the low-latitude (non-inverse Compton) emission prefers the morphology of the NFW template over being absorbed into the low-latitude portion of the Bubbles. The significant improvement in the fit is driven both by the spherical symmetry of the NFW profile, and the fact that the imposed Bubbles templates are flat in brightness rather than falling off with increasing Galactocentric radius (see also Appendix D). With the inclusion of this NFW template in the fit, the spectrum of the emission associated with the Bubbles is not only largely constant at high latitudes (as expected), but also has a similar (near-flat) spectral shape at low latitudes, at least where it can be detected (in several cases there is no significant detection of Bubbles-correlated emission in the lowest-latitude band, after the addition of the NFW template). The significant GeV-scale feature identified in Sec. IV is no longer absorbed by these Bubbles templates, but instead is present in the spectrum of the spherically symmetric NFW template. The spectrum as-sociated with the NFW template peaks at ∼1-2 GeV, and falls off above ∼10 GeV (there is little or no evidence for emission above 10 GeV in the spectrum of the NFW template, at least in the diffuse model fit). The total flux (integrated over solid angle) associated with this template is ∼ 1 keV/cm 2 /s, corresponding to a total luminosity within the solar circle of ∼ 10 37 erg/s. The regions of the sky that are most important for discriminating between the NFW template and the Bubbles templates are those near the inner Galaxy but outside the Bubbles themselves. Unfortunately, such regions are also near the Galactic plane. As a result, the extracted spectrum associated with the NFW template depends somewhat on the details of the fit, including the mask of the Galactic plane, especially below ∼1 GeV where the angular resolution of Fermi is somewhat poor. In the simplest case, in which we employ the "diffuse model" fit and add the NFW template, the spectrum we extract below 1 GeV is considerably softer than that extracted from the Bubbles templates (see the left frame of Fig. 10). Masking |b| < 5 • hardens the low-energy spectrum, as does adding a separate template to absorb soft-spectrum emission from Loop I, or restricting the fit to the southern hemisphere where the structures associated with Loop I do not contaminate the signal (see left center and right center frames of the figure). Given these results, we find it plausible that the low-energy emission (below 1 GeV) associated with the NFW template can be attributed in large part to contamination from the disk and Loop I.
In Fig. 12, we show how the spectral feature at a few GeV is reallocated from the lowest-latitude slices of the Bubbles to the NFW template, once the NFW template is added. In order to minimize background contamination, we examine the spectra masking |b| < 5 • and restricted to the southern hemisphere. As discussed above, we find that the "bump" is much better fitted by the NFW template and there is no remaining significant evidence for emission correlated with the low-latitude Bub-  The spectrum extracted for the gamma-ray Bubbles in ten-degree latitude bands, now including an NFW template in the fit: in order from the top row, 40 • < |b| < 50 • , 30 • < |b| < 40 • , 20 • < |b| < 30 • , 10 • < |b| < 20 • , and |b| < 10 • . The left and center panels use the "diffuse model" template fit (see text); in the center panels, the fit is restricted to b < 0 in addition to the masking. The right panels use the "low-energy template" fit (see text). In all cases an additional (squared, projected) NFW template is added to the fit, as in Fig. 10 bles. Quantitatively, we find that the fit including the NFW template improves over the five-Bubbles-templates model at the level of approximately 12σ. Again, though, we reiterate that this is a formal significance derived including statistical errors only, whereas there are significant systematics originating from the fact that the models we are using do not fully describe the data. We have performed fits using generalized NFW profiles with a variety of inner slopes, γ. Fitting over both hemispheres and masking within 1 degree of the Galactic plane, we find that γ 1.2 provides the best fit (for this reason, we have adopted γ =1.2 as our default value). If we restrict our fit to the southern hemisphere and mask within five degrees of the Galactic plane, however, we find that a somewhat steeper distribution is preferred, γ ∼ 1.5 − 2. This may reflect a slope which varies with distance from the GC, although systematic uncertainties prevent us from making any strong statement to that effect. We show the ∆ ln L values for various choices of γ, in both these cases, in Fig. 13.
We note that it is not surprising that a spherically symmetric signal centered around the Galactic Center would first become apparent in an analysis of the region of the Fermi Bubbles. By their shape, the Bubbles exclude both most of the background from the disk, and the arc structures associated with Loop I, making them reasonably well designed for extracting even a spherically symmetric excess (note the similarities between the Bubbles regions and the dark matter regions-of-interest as described in Refs. [33][34][35]).

VIII. INTERPRETATION
In the existing literature, a number of possibilities have been discussed for the origin of the extended gamma-ray emission observed from the region of the Galactic Center, including annihilating dark matter [18][19][20]22], a population of millisecond pulsars [18,19,[22][23][24], and cosmic ray interactions with gas [18,19,22,25,26]. In this section, we discuss these possibilities in light of the new information provided in this study, including the evidence that this emission is not confined to the inner Galaxy, but instead extends out to at least 2-3 kpc from the Galactic plane.

A. Diffuse Emission Mechanisms
As we have previously argued in Sec. IV, no realistic spectrum of cosmic ray electrons can produce the observed spectral features of the GeV emission identified in this study. Cosmic ray protons scattering with gas also fail to provide a reasonable explanation for the observed gamma-rays. Firstly, the morphology of the observed emission is not highly correlated with the distribution of gas in the Milky Way, as would be predicted in such a scenario. Secondly, the spectral shape of the observed gamma-ray emission requires a peculiar spectrum of cosmic ray protons, sharply peaked at approximately 25 GeV. In the upper left frame of Fig. 14, we compare the observed spectrum of this emission to the spectral shape predicted from proton collisions with gas for a proton spectrum given by a delta function at 25 GeV (solid ), or by a broken power-law, following E −2 p below 25 GeV and E −3 p at higher energies (dashed ). While neither provides a particularly good fit, it is clear that from this comparison that a very sharply peaked spectrum of protons would be required to potentially account for this emission. And while we cannot absolutely rule out the existence of such a strongly peaked cosmic ray proton spectrum, this rather extreme requirement further disfavors cosmic ray origins for the observed low-latitude emission.
Note that in each frame of Fig. 14, we make comparisons to the observed gamma-ray spectrum (after subtracting the contribution from inverse Compton), as extracted using the diffuse model template and from the |b| = 10 • − 20 • regions of the Bubbles. We have chosen to compare to this region over that extracted from lower latitudes, or from the NFW-template, as this region provides a spectrum that is more robust to contamination from the Galactic Disk. In particular, a comparison of the low-energy spectrum extracted from the |b| < 10 • Bubbles template or (especially) the NFW templates varies considerably depending on how the disk is masked, and on whether we consider both hemispheres, or only the southern sky. On the other hand, the spectral shape extracted at higher latitudes depends more strongly on the subtraction of the inverse Compton component, which is negligible compared to the bump for |b| < 10 • , but increasingly substantial as one moves farther from the Galactic plane. However, the consistency in the spectrum of the bump between the 10−20 • band and even higher latitudes supports the simple picture where the electron spectrum generating the ICS component (or the proton spectrum, in a hadronic scenario) is essentially invariant.

B. Millisecond Pulsars
Instead of relying on diffuse emission mechanisms, a large population of unresolved gamma-ray point sources could, in principle, be responsible for the observed emission. In particular, the spectrum of gamma-ray pulsars has been observed to peak at GeV energies, and it has been suggested that a collection of ∼10 3 such objects may be able to account for the extended gamma-ray emission observed from the Galactic Center [22][23][24]. Of particular interest are millisecond pulsars, which are thought to be formed as parts of low-mass X-ray binary systems. This could help to accommodate the very concentrated distribution of the gamma-ray emission observed around the Galactic Center (the stellar distribution in the inner Galaxy scales as n ∝ r −1.2 , but objects formed through FIG. 14: Comparisons of the observed gamma-ray spectrum of the low-latitude emission, after subtracting the contribution from inverse Compton scattering (see Fig. 7) to that predicted from the scattering of cosmic ray protons with gas (upper left), millisecond pulsars (upper right), and dark matter annihilations (lower ). For proton-gas collisions, the solid and dashed lines denote cosmic ray proton spectra which take the form of a delta function at 25 GeV or a broken power-law following E −2 p below 25 GeV and E −3 p at higher energies, respectively. To accommodate the spectral shape of the observed gamma-ray emission, the cosmic ray proton spectrum throughout the inner several kiloparces of the Fermi Bubbles must peak very strongly at approximately 25 GeV. The spectrum shown for pulsars is that corresponding to the average millisecond pulsar observed by the Fermi collaboration [36,37]. For annihilating dark matter, we show results for two models: 10 GeV particles annihilating to tau leptons (dashed) and 50 GeV particles annihilating to bb. In each case, we have adopted a generalized NFW profile with an inner slope of γ = 1.2, and normalized the signal to a local density of 0.4 GeV/cm 3 and an annihilation cross section of σv = 2 × 10 −27 cm 3 /s (τ + τ − ) or σv = 8 × 10 −27 cm 3 /s (bb).
interactions between stars could plausibly be distributed as steeply as the square of this distribution) [22,38]. The binary companions of such pulsars could also act as a tether, explaining why they do not free-stream out of the Galactic Center as a result of velocity kicks. Furthermore, millisecond pulsars are spun up through accretion, and can thus produce high luminosities of gamma-ray emission over much longer timescales than other types of pulsars.
Unfortunately, relatively little is known about the gamma-ray spectrum of millisecond pulsars. Fermi has reported spectra from only eight millisecond pulsars, which together yield an average spectrum that is well fit by dN γ /dE γ ∝ E −1.5 γ exp(−E γ /2.8 GeV) [36,37]. In the right frame of Fig. 14, we show this spectral shape (with arbitrary normalization) compared to the observed emission from the low-latitude Bubbles. This does not provide a good fit, especially at low energies (although this is also where the potential systematic errors are most significant). One could imagine, however, that the eight millisecond pulsars reported by Fermi may not be representative, perhaps being biased toward the brightest or most locally common examples of such objects. To address this issue, it has been suggested that the gamma-ray spectra of globular clusters (which are thought to contain large numbers of millisecond pulsars) could provide a more reliable determination of the average spectrum from millisecond pulsars [23]. At present, however, the error bars on the gamma-ray spectra of globular clusters are quite large, and (on average) do not appear to favor a much harder spectral index than is observed from individual resolved pulsars [39] (see also Fig. 9 of Ref. [19]). Future data from Fermi could potentially be useful in further testing the possibility that millisecond pulsars produce a gamma-ray spectrum compatible with the low-latitude, GeV emission under consideration here.
In addition to their gamma-ray spectra, the distribution of millisecond pulsars within the Milky Way is not well constrained empirically. That being said, one expects the formation of millisecond pulsars to roughly follow that of stars, possibly with an additional preference for regions of very high stellar density, such as in globular clusters. To account for the observed morphology of this signal, however, there would have to be a significant number of millisecond pulsars in the halo, at distances of at least a few kpc outside the Galactic plane. If this emission is generated by such objects, they would require a distribution that is very different from that observed among other stellar populations. And while such distributions have been proposed [40], such distributions are constrained by the small number of millisecond pulsars resolved by Fermi and by the high degree of isotropy observed in the gamma-ray background [41].

C. Annihilating Dark Matter
In the previous two sub-sections, we have described some of the possibilities and challenges involved in explaining this gamma-ray signal with astrophysical sources or mechanisms. Annihilating dark matter can provide a simple explanation for the sharply peaked spectrum and distinctive morphology of this emission. In the lower frame of Fig. 14, we compare the gamma-ray spectrum predicted in two dark matter scenarios to the observed spectrum. First, we consider a 10 GeV dark matter particle species annihilating to tau leptons. We also show the spectrum resulting from a 50 GeV particle annihilating to bb which, given the systematic uncertainties in the extraction of the spectrum, we also consider to be a viable possibility. For a generalized NFW profile (γ = 1.2) normalized to a local density of 0.4 GeV/cm 3 , we require an annihilation cross section of σv = 2 × 10 −27 cm 3 /s in the 10 GeV τ + τ − case, or σv = 8 × 10 −27 cm 3 /s in the case of a 50 GeV particle annihilating to bb. As annihilations to electrons, muons, or neutrinos do not significantly contribute to the gamma-ray spectrum, the total annihilation cross section in the leptonic case may be larger by a factor of a few or several.
We note that the dark matter distribution required to fit the observed signal is well supported by the results of numerical simulations. As has been known for some time, N-body simulations which model the evolution of cold dark matter without baryons typically find halos of the NFW form, with inner slopes of ρ ∝ r −1 (γ = 1) [42,43]. Hydrodynamical simulations, which include the effects of baryonic processes involved in galaxy formation and evolution, in many cases predict the steepening of this inner slope, typically from γ = 1 to γ = 1.2 to 1.5 (see Ref. [44] and references therein). Supporting this picture further are observations of the Milky Way's rotation curve and microlensing optical depth, which together find that the dark matter distribution is best fit by γ = 1.3 (although with large uncertainties) [45]. 4 The annihilation cross section required to normalize the signal in question is also an attractive value from a theoretical perspective. To be produced thermally in the early universe in an abundance equal to the measured dark matter density, a (∼10-50 GeV) dark matter particle must have an annihilation cross section of σv ≈ 2 × 10 −26 cm 3 /s, when thermally averaged over the process of freeze-out [48]. In many typical dark matter models, however, the annihilation cross section in the low velocity limit (as is relevant for annihilation in the halo) is smaller than this value by a factor on the order of a few. Neutralinos, for example, when annihilating to leptons through the t-channel exchange of sleptons typically exhibit σv FO /σv v=0 ∼ 5. While neutralino annihilations through p-wave amplitudes (which are suppressed at low velocities) significantly contribute to the relic abundance calculation, they do not factor into the current annihilation rate. We also note that annihilation cross sections in the range required here (σv ∼ (2 − 8) × 10 −27 cm 3 /s) are currently consistent with all observational and experimental constraints, including those derived from dwarf galaxies [49,50], the cosmic microwave background [51][52][53], the spectrum of cosmic-ray antiprotons [54,55], and from searches for monophoton-plus-missing energy events at LEP [56].
From a dark matter model building perspective, there are a number of possibilities one could consider. In the case of light dark matter particles which annihilate primarily to leptons, the annihilations could proceed through the t-channel exchange of new charged particles which carry lepton number (such as sleptons), or thorough the exchange of a leptophillic Z [57] or a leptophillic higgs [58]. Alternatively, one could imagine a scenario in which the dark matter itself carries lepton number. The dark matter could also annihilate to a pair of light force carriers which interact with the Standard Model only through kinetic mixing the with photon. The decays of the light force carrier then yield combinations of mesons and charged leptons, leading to a gamma-ray spectrum that is very similar to that predicted from tau decays [59].
One additional consideration is that if the ∼10 GeV dark matter particles also annihilate to e + e − at a rate similar to τ + τ − , then they would also be capable of explaining the anomalous synchrotron emission observed from the Milky Way's radio filaments [60], and possi-bly a significant fraction of the isotropic radio background [61,62]. Such annihilations could also lead to a feature in the cosmic ray positron spectrum, potentially identifiable by AMS [63].

IX. DISCUSSION AND CONCLUSIONS
In this study, we have identified a new component of gamma-ray emission from the low-latitude regions of the Fermi Bubbles. This emission does not appear to be compatible with originating from either inverse Compton scattering or proton collisions with gas (for any physically realistic spectrum of cosmic rays). The spectrum of this emission peaks at energies of a few GeV and is distributed with an approximately spherically symmetric morphology about the Galactic Center, with a luminosity per volume that is consistent with a squared (generalized) NFW profile with an inner slope of ρ ∝ r −1.2 (although the best-fit slope depends somewhat on the region chosen for the fit). The broad features of this signal are robust to the various details of our analysis, such as the choice of templates and the degree of masking around the disk.
The diffuse gamma-ray signal we have identified here is consistent both in spectrum and morphology with being the more spatially extended counterpart of the emission previously observed from the Galactic Center [18][19][20][21][22]. In addition to further confirming the existence of this gamma-ray source, the results presented here are important in that they confirm that the morphology of this signal extends to at least ∼ 2 − 3 kpc from the Galactic Center (whereas previous studies of the Galactic Center were sensitive only to the inner few hundreds of parsecs). This provides us with valuable information with which to test various interpretations of the observed emission. Millisecond pulsars, in particular, have been proposed as a source of the emission observed from the Galactic Center. In such a scenario, the significant high-latitude component of the signal reported in this paper would seem to require a large high-latitude population of faint unresolved millisecond pulsars that does not trace known stellar populations. In contrast, the morphology of the signal is consistent with that predicted a priori from annihilating dark matter.
If annihilating dark matter is responsible for the emission reported here, we find that a halo profile of the generalized NFW form, with an inner slope of ρ ∝ r −1.2 , is preferred. The spectral shape of the emission is well described by an approximately 10 GeV particle annihilating to tau leptons, with a normalization corresponding to an annihilation cross section of σv ∼ 2 × 10 −27 cm 3 /s (in addition to annihilations proceeding to electrons, muons and/or neutrinos, each of which provide comparatively few gamma-rays). We note that dark matter in this mass range may also be capable of accommodating a number of other anomalous signals reported from direct and indirect detection experiments (for a summary, see Ref. [32]). Alternatively, the annihilations of a ∼40-50 GeV dark matter particle to quarks can also provide a reasonable fit, requiring a normalization of σv ∼ 8 × 10 −27 cm 3 /s. At first appearance, the gamma-ray spectra as extracted using the diffuse model and low-energy template methods are quite different (see the left and right frames of Fig. 3). This comparison can be misleading, however, as the spectra extracted using the low-energy template do not describe the emission present in the Bubbles region after subtraction of some physical model -rather, they describe the degree to which the Bubbles-correlated emission is harder than the ordinary ICS emission associated with the disk (which dominates the 0.5-1 GeV skymap after removal of the dust-correlated emission), with a pivot point at 0.5-1 GeV. In contrast, the diffuse model template fit extracts the excess over the diffuse model spatially correlated with the Fermi Bubbles. Of course, the diffuse model is fitted to the data, and so may also soak up emission physically associated with the Bubbles.
In this Appendix, we provide a direct comparison of these techniques. Comparing the gamma-ray spectra extracted using the two methods, we see that there are discrepancies at both low and high energy (with the latter being particularly pronounced at low latitudes). These discrepancies have different origins; we will discuss them in turn.
At high energies and low latitudes, the low-energy template fit returns a roughly flat (in E 2 dN/dE) Bubblescorrelated spectrum, similar to the spectrum obtained at high latitudes. In contrast, the diffuse model fit detects virtually no emission. This indicates that the diffuse model is "soaking up" this hard-spectrum emission in some way. The diffuse model has, of course, been adjusted to fit the data, and has the potential to oversubtract Bubbles-correlated emission in trying to fit the data without the freedom to include the Bubbles explicitly. Given the information to hand, it is not possible to distinguish between this scenario and one where the Fermi diffuse model is adequately capturing the physics of this spectral component, which should be assigned to the ordinary Galactic emission even though it is harder than the norm. The amplitude of this hard spectral component, at low latitudes, is not of great interest for the purposes of this study; to facilitate checking agreement of the few-GeV spectral feature between the two methods, we subtract a component with dN γ /dE γ ∝ E −2 γ from the results for the low-energy template fit, designed to remove this high-energy emission.
The low-energy discrepancy is simpler to deal with, as it is an entirely natural and expected consequence of the different methodologies being employed. Since the low-energy template already includes some contribution from the Bubbles, to obtain the true Bubbles spectrum we should re-add a component with the same spectrum as the low-energy template (dN γ /dE γ ∝ E −2.65 γ , to a good approximation). The normalization of this component is a priori unknown, but to check consistency between the two methods, we can normalize it to match the diffuse model result at 0.5-1 GeV.
Having chosen these two free parameters (the amplitude of the high-energy spectrum and the amplitude of the correction due to the low-energy template), we compare in Fig. 16 the spectra extracted using the diffuse model template (black) to that extracted using the lowenergy template, after performing the adjustments described in this Appendix (blue). After these corrections, the low-energy template results are almost identical to those derived using the diffuse model template.
FIG. 16: A comparison of the gamma-ray spectra extracted from various latitude regions of the Fermi Bubbles using the diffuse model template (black) to that extracted using the low-energy template, after performing the corrections described in the text (blue). After accounting for these corrections, the low-energy template results are almost identical to those derived using the diffuse model template.
unexpected, and motivates the addition of a template to absorb this additional emission.  ), for the diffuse model foreground template, restricting the fit to b < 0 (top left), b > 0 (top right), l < 0 (bottom left), l > 0 (bottom right). The Galactic Disk is masked for |b| < 1 • in all cases.
FIG. 18: The spectrum of the "bubble complement" region, defined as √ b 2 + l 2 < 20 • and that does not lie within the Bubbles templates, as extracted from a fit including the complement template, the five latitudinally sliced component Bubbles templates, and background templates as described in Fig. 3. Here, we have used the Fermi diffuse model. The different frames correspond to different choices for the latitude cut to remove the Galactic disk: |b| < 1 • , |b| < 2 • , |b| < 5 • . (Note that this substantially changes the region over which the complement spectrum is averaged, and thus may truly modify the result.) The dashed line in each frame denotes the average spectrum expected in this region from the dark matter model and distribution as shown in Fig. 7  The residual emission after re-adding the latitude-sliced Bubbles templates with their best-fit coefficients, as a function of energy. Emission is E 2 dN/dE averaged over |l| < 5 • and in bins of width ∆b = 2 • ; the error bars describe the pixelto-pixel scatter within each bin. The upper panel shows the fit using the diffuse model, the lower the fit using the low-energy template, both masked at 5 • from the plane. The absolute value of the residual emission after subtraction of all templates (dotted line) and the emission after re-adding the latitude-sliced Bubbles templates with their best-fit coefficients (red stars), for 1-10 GeV. The total emission is shown in black diamonds. Emission is E 2 dN/dE averaged over |l| < 5 • and in bins of width ∆b = 2 • ; the error bars describe the pixel-to-pixel scatter within each bin. The second and fourth panels show the residual emission after subtraction of all templates (not its absolute value). The upper two panels show the fit using the diffuse model, the lower the fit using the low-energy template; in both cases the fit was performed for |b| > 5 • .