Fermi-LAT Observations of γ-Ray Emission toward the Outer Halo of M31

Abstract

The Andromeda galaxy is the closest spiral galaxy to us and has been the subject of numerous studies. It harbors a massive dark matter halo, which may span up to ~600 kpc across and comprises ~90% of the galaxy's total mass. This halo size translates into a large diameter of 42° on the sky, for an M31-Milky Way (MW) distance of 785 kpc, but its presumably low surface brightness makes it challenging to detect with γ-ray telescopes. Using 7.6 yr of Fermi Large Area Telescope (Fermi-LAT) observations, we make a detailed study of the γ-ray emission between 1-100 GeV toward M31's outer halo, with a total field radius of 60° centered at M31, and perform an in-depth analysis of the systematic uncertainties related to the observations. We use the cosmic-ray propagation code GALPROP to construct specialized interstellar emission models to characterize the foreground γ-ray emission from the MW, including a self-consistent determination of the isotropic component. We find evidence for an extended excess that appears to be distinct from the conventional MW foreground, having a total radial extension upward of ~120-200 kpc from the center of M31. We discuss plausible interpretations of the excess emission, but emphasize that uncertainties in the MW foreground-and in particular, modeling of the H i-related components-have not been fully explored and may impact the results.

Keywords: Galaxy: halo; astroparticle physics; cosmic rays; dark matter; galaxies: individual (M31); gamma rays: diffuse background.

Figure 37.

Fractional energy residuals resulting from varying the index of the IC components using a PL scaling. Otherwise, the fit is performed in the standard way. The left shows the results for the TR and the right is for FM31.

Figure 38.

Left: the index of the H i-related emission as a function of Galactocentric radius. The black circles give the baseline index for the M31 IEM corresponding to the GALPROP prediction (~2.75). The cyan circle is the best-fit index for the local annulus obtained in the TR (using the M31 IEM), which is consistent with the GALPROP prediction. The red squares show the results for scaling the index of the H i-related components in FM31. The middle ring, A6, has the smallest radial extension—and likewise, it has the largest error bars. We also repeat the fit using the IG IEM, which only has one outer ring. The results for the IG IEM are shown with blue diamonds, and they are qualitatively consistent with the M31 IEM. For comparison, we also show other measurements. The purple upward-pointing triangles are from Acero et al. (2016). For the local ring, the fit includes all longitudes and 10° < |b| < 70°; for the outer Galaxy (last two rings), the fit includes all latitudes and 90° < l < 270°. The gray rightward-pointing triangles are from Yang et al. (2016). The fit is performed in the latitude range |b| < 5°. The green dashed band is also from Yang et al. (2016), and it shows the 1σ average photon index (above 2 GeV) in the region 10° < |b| < 15° and 90° < l < 150°, which corresponds to the M31 direction. Last, the brown dashed curve is a model fit from Recchia et al. (2016), which is based on nonlinear CR propagation in which transport is due to scattering and advection off self-generated turbulence. These other studies find evidence for a gradual softening toward the outer Galaxy. There is clearly a significant anomaly in FM31. Right: fractional energy residuals resulting from scaling the index of the H i-related components, for both the M31 IEM and the IG IEM.

Figure 39.

The isotropic component includes unresolved diffuse extragalactic emission, residual instrumental background, and possibly contributions from other Galactic components that have a roughly isotropic distribution. The spectrum has a dependence on the IEM and the ROI used for the calculation, as well as the data set. For the IG IEM (which uses the isotropic IC sky maps), we calculate the All-Sky (solid black line) isotropic component in the following region: |b| ⩾ 30°, 45° ⩽ l ⩾ 315°. We also calculate the isotropic component in the different sky regions as follows. North: b ⩾ 30°, 45° ⩽ l ⩽ 315° (orange dashed line). South: b ⩽ −30°, 45° ⩽ l ⩽ 315° (green dashed line). East: |b| ⩾ 30°, 180° ⩽ l ⩽ 315° (blue dashed line). West: |b| ⩾ 30°, 45° ⩽ l ⩽ 180° (purple dashed line). The calculations are performed using a log parabola (LP) scaling for the diffuse components. In addition, we calculate the isotropic spectrum at high latitudes (|b| ⩾ 50°), scaling just the normalizations of the diffuse components. The brown squares show the official FSSC isotropic spectrum (iso_P8R2_CLEAN_V6_v06). The gray band is our calculated isotropic component systematic uncertainty for the IG IEM, as shown in Figure 8.

Figure 40.

Flux (upper panel) and fractional count residuals (lower panel) for the fit in FM31 with the IG IEM. The H ii and Bremsstrahlung components are fixed to their GALPROP predictions. The normalizations of the IC, H i-related, and H₂-related components are fit to the γ-ray data in FM31, as well as 3FGL sources within 20° of M31, along with additional point sources that we find using our procedure. The fit is performed with the high-latitude isotropic component fixed to its nominal value (1.0). The bottom panel shows the fractional residuals, and the blue band shows the corresponding fractional residuals for the baseline fit (with IC scaled) with the M31 IEM. For reference, the residuals (data–model) are also plotted in the upper panel (faint gray band).

Figure 41.

Spatial count residuals (data–model) resulting from the fit in FM31 with the IG IEM for three different energy bands, as indicated above each plot. The energy bins are chosen to coincide with the excess observed in the fractional residuals. The color scale corresponds to counts/pixel, and the pixel size is $0_{.}^{°}2\times 0_{.}^{°}2$. The images are smoothed using a 1° Gaussian kernel. This value corresponds to the PSF (68% containment angle) of Fermi-LAT, which is ~1° at 1 GeV. For reference, the position of M33, $(l,b)=({133}_{.}^{°}61,-{31}_{.}^{°}33)$, is shown with a yellow triangle.

Figure 42.

The maps show the difference between the spatial residuals resulting from the baseline fit and the spatial residuals resulting from the 3FGL optimized fit. For the optimized fit, the PL spectral models are replaced with LogParabola spectral models. Three energy bins are shown, just as in Figure 15. Green crosses show 3FGL sources with TS ⩾ 25, and slanted green crosses show 3FGL sources with 9 ⩽ TS < 25. For the baseline fit, numerous 3FGL sources with PL spectral models were overmodeling in bins 1 and 3, and undermodeling in bin 2, as seen in Figure 16. As seen here, in bins 1 and 3, the 3FGL overmodeling is deeper for the baseline fit, resulting in the surrounding blue regions, and in bin 2, the 3FGL undermodeling is more severe for the baseline fit, resulting in the surrounding red regions; i.e., numerous 3FGL sources show improvement in the spatial residuals with the optimized fit.

Figure 43.

All 3FGL sources in FM31 with a PL spectral model are fit with a LogParabola spectral model. The spectral parameters for each source (norm, α, β, E_b) are initially set to the corresponding values for the respective PL spectra, with β initially set to zero. Optimization of the 3FGL sources leads to marginal improvement in the fractional energy residuals, and most notably for the high energy deficit in the last few energy bins. The corresponding differences for each energy bin are reported in Table 19. For the baseline fit, the likelihood value is −log L = 143349; for the optimized 3FGL fit, it is −log L = 143308.

Figure 44.

The top panel shows the best-fit spectra resulting from the FSSC IEM using the Clean data class (with the Galactic diffuse index fixed). The bottom panel shows the resulting fractional count residuals. Black squares are for the Clean class (with the Galactic diffuse index fixed), green circles show the same fit but with the index of the Galactic diffuse component freed, and blue triangles are for the UltraCleanVeto (UCV) class (with Galactic diffuse index fixed). All components are fit in FM31, including the isotropic.

Figure 45.

Spatial count residuals (data-model) resulting from the fit in FM31 with the FSSC IEM for three different energy bands, as indicated above each plot. The energy bins are chosen to coincide with the excess observed in the fractional residuals. The color scale corresponds to counts/pixel, and the pixel size is $0_{.}^{°}2\times 0_{.}^{°}2$. The images are smoothed using a 1° Gaussian kernel. This value corresponds to the PSF (68% containment angle) of Fermi-LAT, which at 1 GeV is ~1°. For reference, the position of M33, $(l,b)=({133}_{.}^{°}61,-{31}_{.}^{°}33)$, is shown with a yellow triangle.

Figure 46.

The M31-related components (not including the arc template) are added to the model and fit with the FSSC IEM. For this fit, the normalization of the isotropic component is held fixed to its best-fit value obtained in the baseline fit (1.04). All other components are fit simultaneously in the standard way. Two variations of the fit are performed. In one variation, the M31-related components are given PL spectral models (dashed blue curves). In the second variation, the M31-related components are fit with a power law per every other energy band over the range 0.3–300 GeV (dashed gray curves). The free parameters include an overall normalization, as well as the index in each respective energy bin. Corresponding results for the M31 IEM are shown with solid purple curves.

Figure 47.

TS map before and after the additional point sources are included in the model. Note that M31 is modeled with an elliptical template based on the IRIS 100 μm map of the galaxy The region shown is a 14° × 14° square, centered at M31 (white circle). The color scale corresponds to the TS value (2Δ log L), as calculated by gttsmap. Overlaid on the initial TS map are the positions of the additional point sources that we find with procedure. Point sources with TS ⩾ 25 are shown as red crosses, and sources with 9 ⩽ TS < 25 are shown as angled red crosses.

Figure 48.

The blue solid line shows a smooth NFW halo appropriate for warm DM models that do not produce significant structure below the dwarf galaxy scale. The parameters of the NFW profile are as follows: mass = 10¹² M_☉, concentration = 11.2, R_virial = 210 kpc, R_scale = 18.9 kpc, γ = 1.0. The green (lower) dashed line, labeled NFW + Substructure (Low), shows the expected DM for a typical ΛCDM cosmology with thermal WIMP DM. For the NFW + Substructure (Low): overall boost factor = 4.2, substructure fraction = 13%, minimum halo mass = 10⁻⁶ M_☉. The corresponding MW contribution along the line of sight is shown with a thin dashed purple line. The black (upper) dashed line, labeled NFW + Substructure (high), shows a scenario in which DM is produced very cold, such that the minimum mass structures form with very high concentrations. These smallest structures would dominate the annihilation signal. For the NFW + Substructure (High), we use results from Gao et al. (2012).

Figure 1.

Observed counts (left) and saturated counts (right) for a 60° radius centered at M31, and an energy range of 1–100 GeV. The green dashed circle (21° in radius) corresponds to a 300 kpc projected radius centered at M31, for an M31–MW distance of 785 kpc, i.e., the canonical virial radius of M31. Also shown is M31 ’s population of dwarf galaxies. M31 and M33 are shown with cyan triangles, and the other dwarfs are shown with 1 ° green circles, each centered at the optical center of the respective galaxy. The sizes of the circles are a bit arbitrary, although they roughly correspond to the point-spread function (PSF; 68% containment angle) of Fermi-LAT, which is ~1° at 1 GeV. Most of the MW dwarfs are not detected by Fermi-LAT, and so we do not necessarily expect the individual M31 dwarfs to be detected. The primary purpose of the overlay is to provide a qualitative representation of the extent of M31 ’ s outer halo, and to show its relationship to the MW disk. Note that ~3 dwarfs (which are thought to be gravitationally bound to M31) reach as for as ~300 kpc, with one dwarf (And XXVIII) reaching as far as ~360 kpc, as seen in the figure.

Figure 2.

The LIS for CR protons (top), He (middle), and all electrons (e⁻ + e⁺; bottom). The latest AMS-02 measurements from Aguilar et al. (2014, 2015a, 2015b) are shown with red squares. The green dashed line shows the results from Boschini et al. (2017, 2018a), which we employ GALPROP and HelMod together in an iterative manner to derive the LIS. We adopt their derived GALPROP CR parameters, and the LIS for our IEM (M31 IEM: solid black line) are roughly the same. The thin dotted black line shows the LIS modulated with HelMod (Boschini et al. 2017, 2018a). Yellow triangles show the Voyager 1 p and He data in the local interstellar medium (Cummings et al. 2016). Voyager 1 electron data are below 100 MeV and therefore are not shown. In addition, we show the LIS for the (“Yusifov”) IEM in Ajello et al. (2016), which we use as a reference model in our study of the systematics for the M31 field (see Appendix B.2).

Figure 3.

The total interstellar emission model (IEM) for the MW integrated in the energy range 1–100 GeV. The color corresponds to the intensity, and is shown in logarithmic scale. The intensity level is for the initial GALPROP output, before tuning to the γ-ray data. The map is shown in a Plate Carrée projection, and the pixel size is 0.25 deg/pix. The model has contributions from π°-decay, (anisotropic) IC emission, and Bremsstrahlung. Overlaid is the ROI used in this analysis. From the observed counts (Figure 1) we cut an 84° × 84° ROI, which is centered at M31. The green dashed circle is the 300 kpc boundary corresponding to M31’s canonical virial radius (of ~21°), as also shown in Figure 1. We label the field within the virial radius as field M31 (FM31), and we label the region outside (and south of latitudes of $-{21}_{.}^{°}57$) as the tuning region (TR). Longitude cuts are made on the ROI at l = 168° and l = 72°, as discussed in the text. For reference, we also show the region of the Galactic Center (GC), which corresponds to a 15° × 15° square centered at the GC.

Figure 4.

Schematic of the eight concentric circles that define the annuli (A1–A8) in the IEM, as described in the text. The ranges in Galactocentric radii are reported in the legend. Note that the full extension of A8 is not shown. Only A5–A8 contribute to the Galactic foreground emission for the field used in this analysis.

Figure 5.

Gas-related components of the IEM (π°-decay related to H i, H ii, and H₂, and Bremsstrahlung emission) integrated in the energy range 1–100 GeV. The components correspond to different annuli, as indicated above each plot. The color corresponds to the intensity, and is shown in logarithmic scale. The intensity level is for the initial GALPROP outputs, before tuning to the γ-ray data. The maps are shown in a Plate Carrée projection, and the pixel size is 0.25 deg/pix. Overlaid is the ROI used in this analysis, as well as the GC region (see Figure 3).

Figure 6.

Anisotropic Inverse Compton (AIC) components of the interstellar emission model for the MW in the energy range 1–100 GeV. The color corresponds to the intensity, and is shown in logarithmic scale. The intensity level is for the initial GALPROP outputs, before tuning to the γ-ray data. The map is shown in a Plate Carrée projection, and the pixel size is 0.25 deg/pix. The IC A6 and A7 components are highly degenerate, and so we combine them into a single map A6+A7. Overlaid is the ROI used in this analysis, as well as the GC region (see Figure 3). Note that we use the anisotropic IC maps as our default component. Unless otherwise stated, all reference to the IC component implies the anisotropic formalism.

Figure 7.

The IEM employs the anisotropic IC sky maps, as discussed in the text. For comparison, we show the differential flux ratio (AIC/IC) between the anisotropic (AIC) and isotropic (IC) inverse Compton components (all-sky). The top figure shows the spatial variation of the ratio at 1 GeV. The bottom figure shows the energy dependence of the ratio for four different spatial points, including M31. The ratio is close to unity toward the GC, increases with Galactic longitude and latitude, and reaches maximum at midlatitudes toward the outer Galaxy. Note that we use the anisotropic IC maps as our default component. Unless otherwise stated, all reference to the IC component implies the anisotropic formalism.

Figure 8.

The spectrum of the isotropic component has a dependence on the IEM and the ROI used for the calculation, as well as the data set. For the M31 IEM (which uses the AIC sky maps), we calculate the all-sky (solid black line) isotropic component in the following region: |b| ⩾ 30°, 45° ⩽ l ⩽ 315°. We also calculate the isotropic component in the different sky regions, as follows. North: b ⩾ 30°, 45° ⩽ l ⩽ 315° (orange dashed line). South: b ⩽ −30°, 45° ⩽ l ⩽ 315° (green dashed line). East: |b| ⩾ 30°, 180° ⩽ l ⩽ 315° (blue dashed line). West: |b| ⩾ 30°, 45° ⩽ l ⩽ 180° (purple dashed line). See Table 2 for the corresponding best-fit normalizations. Magenta triangles show the all-sky isotropic component for the M31 IEM derived using the isotropic IC formalism. The brown squares show the official FSSC isotropic spectrum (iso_P8R2_CLEAN_V6_v06). The gray band is our calculated isotropic systematic uncertainty for the IG IEM, which uses the isotropic IC formalism (see Appendix B.2).

Figure 9.

Total model counts for the full ROI. For the tuning region (TR), we mask within the 300 kpc circle and latitudes above $-{21}_{.}^{°}57$, as discussed in the text.

Figure 10.

Flux (upper panel) and fractional count residuals (lower panel) for the fit in the TR. The H ii component is fixed to its GALPROP prediction. The normalizations of all other diffuse components are freely scaled, as are all 3FGL sources in the region. The residuals show fairly good agreement over the entire energy range.

Figure 11.

Correlation matrix for the fit in the TR. For brevity, IC A6–A7 is labeled as ICA67 and the isotropic component is labeled as Iso.

Figure 12.

Spatial count residuals (data–model) resulting from the fit in the TR for three different energy bands, as indicated above each plot. The energy bins are chosen to coincide with an excess that is later observed in the fractional energy residuals for the fit in FM31, as discussed in the text. The color scale corresponds to counts/pixel, and the pixel size is $0_{.}^{°}2\times 0_{.}^{°}2$. The images are smoothed using a 1° Gaussian kernel. This value corresponds to the PSF (68% containment angle) of Fermi-LAT, which is ~1° at 1 GeV.

Figure 13.

The TS map is calculated after the baseline fit in FM31 (tuned). Overlaid are the additional point sources found by our procedure. Red crosses represent new sources with TS ⩾ 25, and red slanted crosses represent new sources with 9 ⩽ TS < 25.

Figure 14.

Flux (upper panel) and fractional count residuals (lower panel) for the fit in FM31 (tuned). The H ii and Bremsstrahlung components are fixed to their GALPROP predictions. The normalizations of the IC (A5 and A6–A7) and isotropic components are held fixed to the values obtained in the tuning region. The normalizations of the H i- and H₂-related components are fit to the γ-ray data in FM31, as well as 3FGL sources within 20° of M31, in addition to point sources that we find using our procedure. Note that the top of FM31 has contribution from IC A8, and its normalization is also freed in the fit. The fractional residuals show an excess between ~3–20 GeV, reaching a level of ~4% (error bars show 1σ statistical error). Above and below this range, the data are overmodeled as the fit tries to balance the excess with the negative residuals. This is in contrast to the fit in the TR, which shows fairly good agreement over the entire energy range. For reference, the residuals (data–model) are also plotted in the upper panel (faint gray band).

Figure 15.

Spatial count residuals (data–model) resulting from the fit in FM31 (tuned) for three different energy bands, as indicated above each plot. The energy bins are chosen to coincide with the excess observed in the fractional residuals. The color scale corresponds to counts/pixel, and the pixel size is $0_{.}^{°}2\times 0_{.}^{°}2$. The images are smoothed using a 1° Gaussian kernel. This value corresponds to the PSF (68% containment angle) of Fermi-LAT, which is ~1° at 1GeV. For reference, the position of M33, $(l,b)=({133}_{.}^{°}61,-{31}_{.}^{°}33)$, is shown with a yellow triangle.

Figure 16.

Same residual maps as shown in Figure 15. Here, we show the maps in gray scale, and intentionally saturate the images to bring out weaker features. Overlaid are the point sources in the region. Crosses show sources with TS ⩾ 25, and slanted crosses show sources with 9 ⩽ TS < 25. Fermi 3FGL sources are shown in green, and new sources found in this analysis are shown in red.

Figure 17.

Top row: H i column density contours for A5, A6, and A7, as indicated above each plot. For reference, a yellow circle ($0_{.}^{°}4$) centered at M31 is overlaid, and a yellow triangle is overlaid at the position of M33. The units are 10²⁰ cm⁻², and the levels are indicated on the maps. Middle row: the same H i column density contours are overlaid on the residual maps for FM31. The maps are integrated over the entire energy range 1–100 GeV. The residual emission is observed to be correlated with the column densities. In addition, the column densities of A6 and A7 are observed to be correlated with the major axis of M31 (the position angle of M31 is 38°). Bottom row: the same maps as for the middle row, but for a 5° radius centered at M31. Contours for the IRIS 100 μm map of M31 are overlaid. The levels shown range from 6 to 22 MJy sr⁻¹. Also overlaid are the regions corresponding to the two main cuts (space and velocity) that are made on the underlying gas maps when constructing the MW IEM, as detailed in the text. Last, we overlay the 3FGL sources (magenta crosses) in the region with TS ⩾ 25. In particular, we consider the two point sources located closest to the M31 disk, because we are interested in the true morphology of the M31 emission. The source located to the right of the disk (3FGL J0040.3+4049) is a blazar candidate and has an association. The source located to the left of the disk (3FGL J0049.0+4224) is unassociated.

Figure 18.

Additional freedom is given to the baseline fit. The IC components are fit simultaneously with the other contributing diffuse components and point sources. The isotropic component remains fixed to its value obtained in the TR (1.06).

Figure 19.

Correlation matrix for the FM31 baseline fit with the IC components scaled.

Figure 20.

Fractional residuals calculated in different spatial regions. The field is evenly divided into top, middle, and bottom regions. Each slice is then further divided into right and left halves. The regions are indicated above each plot. Black data points show the residuals resulting from the baseline fit (which is over the entire field, with IC scaled in addition to the other contributing components). We then rescale the diffuse components in the different subregions, masking the rest of the region and keeping the point sources fixed to their baseline values (green data points). This is done to allow for a spatially varying spin temperature and/or CR and ISRF densities, which would in turn change the normalizations of the γ-ray components. Even in these smaller regions, the diffuse components are unable to flatten the residuals, with the exception of the bottom right, which is fairly flat.

Figure 21.

The first two panels show the spatial count residuals integrated between 1–100 GeV, resulting from the baseline fit (see Figure 18). In order to construct a template for the large arc extending from the top left corner to the projected position of M33 (arc template), we divide the total residual map into positive residuals (left) and negative residuals (middle). The maps show the geometry used to help facilitate the template construction (the green axes, circle, and ellipse), as detailed in the text. The corresponding geometrical parameters are given in Table 7. The resulting arc template is shown in the far right panel. In addition to fitting the full arc template, we also perform a variation of the fit in which the arc template is divided into a north component (arc north: $b>-{16}_{.}^{°}5$) and a south component (arc south: $b⩽-{16}_{.}^{°}5$), where the spectral parameters of each component are allowed to vary independently. The cut is made right below the bright emission in the upper-left corner, and it allows the north component to be at a different distance along the line of sight than the south component, as discussed in the text. The cyan triangle shows the projected position of M33.

Figure 22.

Spectra and fractional energy residuals resulting from the arc fit. Left: The full arc component is given a PL spectrum, and the normalization and index are fit simultaneously with the other components in the region, just as for the baseline fit. Black dashed lines show the H i A5 (top), A6 (bottom), and A7 (middle) components from the baseline fit (not the arc fit). Note that A7 has a greater radial extension than that of A6, and it likewise has a greater overall flux. Correspondingly, the gray markers (squares, circles, and triangles) show the H i A5–A7 spectra resulting from the arc fit. The blue solid line is the best-fit spectrum for the arc template. The bottom panel shows the remaining fractional residuals. For reference, the residuals (data-model) are also plotted in the upper panel (faint gray band). Right: The arc template is given additional freedom by dividing it into north and south components. The arc components are given PLEXP spectral models, and the spectral parameters (normalization, index, and cutoff) are freely scaled with the other components. Downward-pointing blue and green triangles give upper limits. Bands give the 1σ error. The arc template is unable to flatten the excess between ~3–20 GeV.

Figure 23.

The correlation matrix for the arc north (AN) and south (AS) fit.

Figure 24.

Spatial count residuals resulting from the arc fit. To give a sense of the deviations, here we show the fractional residuals, where we divide by the model counts for each pixel. The residuals are integrated in three energy bins, just as for the residuals in Figure 15. We show residuals from the arc north and south fit with the PLEXP spectral model. Residuals for the full arc fit with the PL spectral model are very similar. As expected, the arc structure no longer dominates the residuals. The position of M33 is indicated with a yellow triangle, and the center of M31 is indicated with a $0_{.}^{°}4$ open circle.

Figure 25.

The average local (A5) emissivity per H atom. The solid gray curve comes from the baseline fit with IC scaled, and it gives the proper estimate of the emissivity in FM31. The dashed gray curve comes from the arc fit with PL spectral model: it only includes the contribution from the H i A5 component, but not the emission associated with the arc. The blue data points (squares) are from Casandjian (2015), and the corresponding error bars are systematic +statistical. The fit includes absolute latitudes between 10–70°. The data points for the different regions (red circles, green upward-pointing triangles, and yellow rightward-pointing triangles) are from Ackermann et al. (2012c), and the corresponding error bars are statistical only (1σ). The teal band shows the total uncertainty (statistical+systematic) from the same analysis (from the erratum). The different regions are among the nearest molecular cloud complexes, within ~300 pc from the solar system. We also plot (as black leftward-pointing triangles) the measurements from Abdo et al. (2009a), as determined from a midlatitude region in the third Galactic quadrant.

Figure 26.

Top panel shows the dust temperature map for FM31, and the bottom panel shows the dust reddening map from Schlegel et al. (1998), as discussed in the text. Overlaid are contours for the arc template. Contours for the IRIS 100 μm map of M31 are also overlaid in the top panel. The cyan triangle shows the (projected) position of M33.

Figure 27.

Left: FM31 residuals from the baseline fit (with IC scaled) with the Loop III shell plotted over it. The two lines correspond to two somewhat different positions and radii obtained from continuum and polarization observations (Vidal et al. 2015). The shell radius is approximate and the shell itself can be several degrees thick. The shaded area gives an idea of the error associated with the parameters of the shell. Right: M31’s virial radius (300 kpc) is shown with a cyan dashed circle, and cyan triangles show the positions of M31 and M33. The gray circles show Loop III at the top and Loop II at the bottom. Loop IIIs (which is only visible in polarization) is shown with a dashed–dotted magenta circle.

Figure 28.

M31-related components are added to the model, along with the arc template and standard baseline components. The left panel is for the full arc template with PL spectral model, and the right panel is for the north and south arc templates with PLEXP spectral model, just as in Figure 22. Black dashed lines show the best-fit spectra for the H i A5 (top), A6 (bottom), and A7 (middle) components. The black dashed-dotted line shows the isotropic component, which remains fixed to its best-fit value obtained in the tuning region, just as for all other fits. The best-fit spectra of the remaining components are similar to that shown in Figure 18, and are omitted here for visual clarity. Downward-pointing triangles give upper limits. Bands give the 1σ error. The bottom panel shows the remaining fractional residuals, which are fairly flat over the entire energy range, and likewise show a normal distribution with a mean of zero.

Figure 29.

A systematic excess can be observed between ~3–20 GeV at the level of ~3–5%. Systematic overmodeling is also present above and below this range. We note that there is one model for which the signal can be flattened (shown with green circles). This results from using the FSSC IEM (intended for point source analysis) and fitting both the isotropic and Galactic diffuse (including the index) in the signal region. The FSSC IEM is not intended for extended source analysis, and this result illustrates how the application of an improper IEM for analysis of largely extended emission can alter the physical results. The M31 IEM is our benchmark model. The different models are as follows. Black squares: FSSC IEM, fitting the isotropic and Galactic diffuse (with index fixed) in the signal region, using Clean data, corresponding to the fit in Figure 44. Blue upward-pointing triangles: same as for the black squares, but using UltraCleanVeto (UCV) data; see Appendix B.3 for details. Green circles: same as for the black squares, but also freeing the index of the Galactic diffuse. Orange diamonds: M31 IEM baseline fit, varying the index of the IC components A5–A8 using a power law scaling, corresponding to the fit in Figure 37. Purple rightward-pointing triangles: M31 IEM baseline fit, varying the index of the H i-related components A5–A8 using a power-law scaling, corresponding to the fit in Figure 38. Note that, in this case, FM31 shows a significant anomaly in the index of the gas-related emission toward the outer Galaxy, as is clearly shown in Figure 38. Blue band: M31 IEM baseline fit, corresponding to the fit in Figure 18. Green band: M31 IEM tuned fit, corresponding to the fit in Figure 14. Pink band: M31 IEM arc fit, corresponding to the fit in Figure 22 (this is our primary model). Black band: inner Galaxy (IG) IEM, corresponding to the fit in Figure 40.

Figure 30.

The fractional count residuals calculated over the different spatial regions corresponding to the spherical halo and far outer halo components, as indicated above each plot. Note that these are the residuals before adding the M31-related components, and they correspond to the spatial residuals shown in Figure 24, resulting from the baseline fit with the arc north and south templates. The goal here is to further examine the symmetry of the residual emission associated with the M31-related components. We consider the northern and southern regions of the templates, where the cut is made at the midpoint of FM31 along the horizontal direction (parallel to the Galactic plane), corresponding to a latitude of $-{21}_{.}^{°}5$. The first column shows the residuals calculated over the entire region, for the spherical halo and far outer halo, respectively. The second column shows the residuals in the north, and the third column shows the residuals in the south.

Figure 31.

The best-fit spectra resulting from the symmetry test fit, where the spherical halo and far outer halo templates are divided into north and south components, and the spectral parameters for each component are allowed to vary independently. The cut is made at the midpoint of FM31 along the horizontal direction (parallel to the Galactic plane), corresponding to a latitude of $-{21}_{.}^{°}5$. The northern components are shown with square markers, and the southern components are shown with circle markers. Downward-pointing triangles give upper limits. Also overlaid are the spectra for the full component fit (with arc north and south), as shown in Figure 28.

Figure 32.

Correlation matrix for the symmetry test fit. In addition to the standard components, the fit includes components for the arc north and south (AN and AS), inner galaxy (not shown here), spherical halo north and south (SHN and SHS), and the far outer halo north and south (FHN and FHS).

Figure 33.

Left: radial intensity profile for the M31-related components. Red square markers show the results from the north and south arc template with PLEXP. The profiles for the PL arc fit are basically the same. Purple circle markers show the results from the fit with the M31-related templates divided into north and south components (from Figure 31). For reference, we compare the radial profile to expectations for DM annihilation in the line of sight. Note that this also includes the contribution from the MW’s DM halo in the line of sight, which has not been accounted for in our analysis and may be at least partially embedded in the isotropic component and Galactic diffuse components. Likewise, the M31-related components may contain a significant contribution from the MW’s extended halo. Details regarding the DM profiles are given in Appendix C. Right: spectral shape comparison to the Galactic center excess (for an arbitrary normalization), as observed in Ajello et al. (2016). Also shown is a prediction for CRs interacting with the ionized gas of the circumgalactic medium from Feldmann et al. (2013). Note that the prediction is for an MW component, but we are primarily interested in a spectral shape comparison.

Figure 34.

Residual maps showing the structured emission integrated in the energy range 1–100 GeV. The color scale corresponds to counts/pixel, and the pixel size is $0_{.}^{°}2\times 0_{.}^{°}2$. The images are smoothed using a 1° Gaussian kernel. This value corresponds to the PSF (68% containment angle) of Fermi-LAT, which is ~1° at 1 GeV. Maps are shown in the cubehelix color scheme (Green 2011). In the top row, contours for the IRIS 100 μm map of M31 are overlaid, and three zoom levels (2°, 7°, full field) centered at M31 are shown. The white circle (1°) shows the position of M33. The bottom row shows two zoom levels (1°, 3°) centered at M33, and the H i integrated intensity map (units of K) of M33 is overlaid. In the third panel, we show the M31 zoom 0 map, rescaled in order to provide a sense of the relative intensity toward the MW disk. We stress that these maps have not subtracted any Galactic H i-related emission.

Figure 35.

Pixel distribution of the smoothed residual map (1 GeV–100 GeV) after removing the H i-related components, as shown in Figure 34. The yellow dashed lines are at 0 and 4 counts.

Figure 36.

The structured γ-ray emission in FM31 is overlaid with some M31-related objects observed at other wavelengths. We stress that this is only done as a qualitative gauge of M31’ s outer halo. In the figure we have not subtracted any Galactic H i-related emission, and we do not expect the M31-related observations to outshine the MW emission, as discussed in the text. Contours for the IRIS 100 μm map of M31 are overlaid. The solid cyan circle ($0_{.}^{°}4$) shows the boundary of the FM31 inner galaxy component, and the black dashed circle ($8_{.}^{°}5$) shows the outer boundary of the FM31 spherical halo component, as detailed in Section 3.4. Overlaid are H i emission contours from the HI4PI all-sky survey based on EBHIS and GASS (Bekhti et al. 2016), integrated over the velocity range −600 km s⁻¹ ⩽ V_LSR ⩽ −95 km s⁻¹. M31’s confirmed globular clusters are shown with black stars. M31’s population of dwarf galaxies is shown with open black triangles. The M31 cloud can be seen (albeit obscured by globular clusters). We note the serendipitous enclosure by the spherical halo of the M31 cloud, as well as a majority of M31’s globular cluster population and dwarf galaxies. H i contours corresponding to M33 can be seen in the lower-left corner. The hook-shaped gas cloud to the right of M33 is Wright’s cloud. The red gas contours toward the top of the map are clouds of Complex H. The black H i contours toward the top of the field correspond to the plane of the MW, and likewise for the bright (white) γ-ray emission. To the far right of the field, a bright arm of emission extends to higher latitudes. Although not considered when making the overlay, the M31-related observations can be seen to trace the left boundary of the arm. This may be an observational bias, due to foreground gas and dust. We stress that these maps have not subtracted any Galactic H i-related emission.