Constraints on the Intergalactic Magnetic Field from Bow Ties in the Gamma-ray Sky

TL;DR

This study uses gamma-ray observations to place new constraints on the intergalactic magnetic field, excluding certain strengths and coherence lengths based on the absence of expected gamma-ray halos around blazars.

Contribution

It introduces a novel method using angular power spectra to constrain the intergalactic magnetic field by analyzing gamma-ray halos, improving previous limits.

Findings

Excludes IGMF with coherence length >100 Mpc at >3.9σ

Limits IGMF strength to between 10^{-16} and 10^{-15} G at current epoch

Provides a direct non-detection constraint on gamma-ray halos

Abstract

Pair creation on the cosmic infrared background and subsequent inverse-Compton scattering on the CMB potentially reprocesses the TeV emission of blazars into faint GeV halos with structures sensitive to intergalactic magnetic fields (IGMF). Previous work has shown that these halos are then highly-anisotropic and extended. If the coherence length of the IGMF is greater than the inverse-Compton cooling length of the pairs, then the orientation of the gamma-ray halo will be correlated with the direction of the magnetic field which is unknown and expected to change for each source. In order to constructively add each source we then use angular power spectra which are insensitive to the jet orientation. By looking at known GeV blazars detected by Fermi, we exclude the existence of an IGMF with coherence lengths $>100$ Mpc at greater than $3.9\sigma$ with current-day strengths in the range $10^{-16}$ to $10^{-15}$ G, and at 2$\sigma$ from $10^{-17}$ to $10^{-14}$ G. This provides a direct measurement of the non-existence of gamma-ray halos, providing an important check on previous results.

Equations

$${\mathcal{P}}_m=\frac{1}{N_{\mathrm{src}}}\sum _{\mathrm{src}}\frac{1}{N^2}{\left|\sum _j{e}^{im{\theta }_j}\right|}^2,$$

$$n_j^s={\bar{n}}_j+\sum _i{f}_{ji}^s{N}_i$$

$$\sum _i\left(f_{ji}^{1.2}-f_{ji}^{4.7}\right)N_i=\left(n_j^{1.2}-n_j^{4.7}\right).$$

$$N^{-1}\frac{dN}{d{\log }_{10}{F}_{35}}=0.45{\left(\frac{F_{35}}{10^{-9}\mathrm{ph}{\mathrm{s}}^{-1}\mathrm{c}{\mathrm{m}}^{-2}}\right)}^{-1.1}.$$

$$w(F_{35})=\{\begin{array}{ll}\exp [-12.88({\log }_{10}{F}_{35}/F_{\mathrm{th}}{)}^2]& F_{35}<F_{\mathrm{th}}\\ 1& \text{otherwise.}\end{array}$$

$$F_{\mathrm{th}}=4.27\times 10^{-10}|b/90^{\circ }{|}^{-0.33}\mathrm{ph}\mathrm{c}{\mathrm{m}}^{-2}{\mathrm{s}}^{-1}.$$

$$p_1(F_{35})=1-p_0(F_{35})=1-e^{-F_{35}ϵ}\approx F_{35}ϵ,$$

$${\mathcal{P}}_m^h={\mathcal{P}}_m^l+{\mathcal{P}}_m^{ss}$$

$${\delta }_m=\sqrt{\frac{{\Delta }_m^4}{{\Delta }_m^2+{\Sigma }_m^2}},{\sigma }_m^2=\ln \left(\frac{{\Sigma }_m^2+{\Delta }_m^2}{{\Delta }_m^2}\right),$$

Full text

BL Lacertae objects: general — gamma rays: diffuse background — gamma rays: general — infrared: diffuse background — plasmas — radiation mechanisms: non-thermal

††software: astropy (Price-Whelan et al., 2018), matplotlib (Hunter, 2007), scipy (Jones et al., 2001–)

1 Introduction

Astronomical magnetic fields appear on scales ranging from the terrestrial to those of galaxy clusters, the latter spanning many Mpc. Within all of these systems magnetic fields are believed to play a variety of important roles. These include providing additional pressure components (Beck, 2016), mediating the magneto-hydrodynamic turbulence (Goldreich & Sridhar, 1995), providing and limiting the effective viscosity that transports angular momentum in accretion disks (Balbus & Hawley, 1991) and mediates collapse (Hennebelle & Falgarone, 2012; Crutcher, 2012), contributing to the acceleration of high energy particles, e.g., cosmic rays (Blandford & Eichler, 1987), and dictating their subsequent propagation (Hanasz & Lesch, 2003; Pakmor et al., 2016) as well as diffusion of thermal energy and momentum in weakly collisional and collisionless systems (Schekochihin et al., 2009).

Observed astrophysical magnetic fields are likely amplified from weak primordial seed fields. These seed fields may originate from Biermann battery processes (Kulsrud & Zweibel, 2008), electroweak or quantum chromo-dynamics phase transitions in the early universe, or during cosmic inflation. Biermann battery processes can naturally generate magnetic fields when temperature and density gradients are not parallel, resulting in electric fields that have a curl and providing a source for the magnetic field. Typically, these processes have a correlation length that reflects the length scales associated with the temperature and density gradients (Kulsrud & Zweibel, 2008; Widrow et al., 2012), and generate small scale ($\sim 10$ kpc) fields ($B\sim 10^{-20}$ G), which can be amplified by the action of a dynamo to $\mu$G fields observed in galaxies today (Schober et al., 2013; Beck, 2016; Pakmor et al., 2017).

By their nature, fields generated in primordial processes may be present in the intergalactic medium as the IGMF. While the general nature of these fields in still unknown, causal constraints limit the coherent length scales of primordial fields generated from post-inflation mechanisms, e.g., phase transitions, to small comoving scales, typically $<10$ kpc (Widrow et al., 2012). Small scale fields may seed large scale fields through the inverse cascade up to $\sim 0.1$ Mpc, i.e., the sonic scale of $10^{4}$ K gas over the age of the universe (Ryu et al., 2012).

To generate large-scale ($\gtrsim 1$ Mpc) volume-filling fields, seed fields must have been generated prior to or during inflation (Biermann, 1950; Turner & Widrow, 1988) or from second-order interactions between the photon and electron fluid prior to recombination (Ichiki et al., 2006). In the former case, additional physics of super-adiabatic magnetic amplification (Tsagas & Kandus, 2005) is needed to keep these fields relevant today. This leaves the strength of these inflationary fields unconstrained (Widrow et al., 2012). In the latter case, the expected comoving field strength is $B\sim 10^{-24}$ G.

Probing these primordial fields is difficult, but some constraints can be set on a large-scale IGMF. The best direct upper limits on the strength of a putative large-scale IGMF to date come from the impact of such fields on the Cosmic Microwave Background (CMB), roughly $10^{-9}$ G, followed closely by limits arising from Faraday rotation measurements, roughly $10^{-9}$ to $10^{-8}$ G depending weakly on the coherence length of the field $\lambda_{B}$ (Neronov & Vovk, 2010). Indirect limits on large-scale IGMFs can be obtained from the rotation measures of galaxy clusters, constraining the IGMF on scales exceeding a Mpc to be less than $10^{-12}$ G (Dolag et al., 2005; Neronov & Vovk, 2010). In contrast, indirect lower-limits can be placed by the observed absence of inverse-Compton gamma-ray cascades.

Very-high energy gamma rays (VHEGRs), with energies exceeding $E>100$ GeV, lie above the pair-creation threshold with the infrared background (Gould & Schréder, 1967; Salamon & Stecker, 1998). As a result, the mean free path of VHEGRs in the nearby universe is $10^{2}(1+z)^{-4.5}(E/6~{}{\rm T}{\rm eV})^{-1}$ Mpc (Broderick et al., 2012). That is, while VHEGRs can traverse intergalactic distances, they are absorbed on cosmological distances, producing highly relativistic beams of electron-positron pairs (i.e., Lorentz factors of order $10^{6}$). This redshift-dependent absorption of high energy gamma rays has been observed (Ackermann et al., 2012a) and subsequently used to probe the infrared background (Domínguez et al., 2011; Gilmore et al., 2012; Domínguez et al., 2013).

In the absence of any other process, e.g., collective beam-plasma instabilities (Broderick et al., 2012), these beams will inverse-Compton cool by upscattering cosmic microwave background (CMB) photons to GeV energies. That is, this chain of events would effectively reprocess the TeV emission of AGN into GeV emission. Thus, bright, nearby, extragalactic TeV sources should exhibit a GeV bump that is correlated with their TeV spectral energy distribution (SED). Stringent limits have been placed by Fermi on the existence of such a GeV bump already, effectively excluding such a feature in several TeV sources (Neronov & Semikoz, 2009; Neronov & Vovk, 2010).

Presuming that inverse-Compton cooling dominates the evolution of the beams, these non-detections provide indirect, circumstantial evidence for the existence of a strong IGMF, with $B>10^{-15}$ G to $10^{-17}$ G, depending on the particular extragalactic background light model employed (Neronov & Vovk, 2010; Taylor et al., 2011; Vovk et al., 2012). This assumes that the duty cycle of the gamma-ray blazars is comparable to that observed for the radio mode in AGN, roughly $10^{6}$ yr; shorter duty cycles could result in GeV echos that lag behind the initial VHEGR emission by long times, resulting in a weaker limit: $B>10^{-17}$ G (Dermer et al., 2011; Taylor et al., 2011). Such an IGMF would deflect the resulting pairs, whose inverse-Compton emission would then be highly beamed away from us, naturally removing the GeV bump. However, at the same time, this would open the door to the direct detection of the IGMF through the observation of GeV halos surrounding intrinsically gamma-ray bright objects, for which the inverse-Compton emission is beamed into our line of sight.

A number of searches for inverse-Compton halos using the publicly available Fermi observations have been performed, yielding contradictory claims but no convincing evidence for their presence (Ando & Kusenko, 2010; Neronov et al., 2011; Ackermann et al., 2013; Chen et al., 2015). These have focused on the large radial extent of the gamma-ray halos, typically subtending many degrees. However, this is generally degenerate with the point spread function (PSF) of Fermi’s Large Area Telescope (LAT), which until recently was highly uncertain (Neronov et al., 2011; Ackermann et al., 2012b). As result, even when detections are claimed they remain plagued by large systematic uncertainties.

Previous work in Broderick et al. (2018) stacked, rotated Fermi images of Fanaroff-Riley class I and II objects from the VLA FIRST survey (Helfand et al., 2015), to limit the strength of the IGMF to $<10^{-15}G$ for large correlation length magnetic fields. A key difference in this work is that we use Fermi blazars and their full spectral information in the GeV–TeV range (see Tiede et al., 2016) to attempt to measure gamma-ray halos. This is different to the work by Broderick et al. (2018) who statistically sampled the blazar SEDs using the distribution of Fermi and VHEGR spectra. Due to this, the work in this paper provides an additional and independent constraint on the IGMF strength.

The layout of this paper is as follows: In Section 2, we review our method first described in two companion publications (Broderick et al., 2016; Tiede et al., 2016), where we used angular power spectra to coherently stack angular structure from gamma-ray halos. In Section 3 we describe a model for the gamma-ray sky background, which adds power at low multipoles. Section 4, details the main results followed by a discussion of the impact of various potential systematics. Finally, Section 5 summarizes the conclusions of this work.

2 Method

In two companion publications (Broderick et al., 2016; Tiede et al., 2016) we have demonstrated that these inverse-Compton halos are expected to be highly anisotropic (see also Neronov et al., 2010; Long & Vachaspati, 2015). This originates from different physical origins for a small-scale, tangled IGMF ($\lambda_{B}\ll 3$ Mpc) and a large-scale, uniform IGMF ($\lambda_{B}>100$ Mpc). For the latter, of interest here, this is due to extreme beaming of the inverse-Compton emission coupled with the geometric conditions upon the gyration of the particles. This typically produces a bilateral, narrowly beamed halo, the two lobes being due to the opposite gyration directions of the electrons and positrons (Broderick et al., 2016).

In Tiede et al. (2016) we demonstrated that existing Fermi observations should permit the detection of a large-scale IGMF with strengths comparable to those implied by the absence of GeV bumps in nearby TeV sources, i.e., $10^{-16}$ to $10^{-15}$ G. This makes use of stacked one-dimensional angular power spectra that amplify the impact of anisotropy111Note that Broderick et al. (2018) used the VLA survey so that the jet direction is known, so stacking angular histograms coherently was possible.:

[TABLE]

where $\theta_{j}$ is the polar angle of the $j$th photon defined relative to the position of the central source and an arbitrary direction, and $N$ the number of photons for the source. For this we make use of the Pass 8R2_V6 ULTRACLEANVETO event class and limit ourselves to gamma-ray energies between 1 GeV and 100 GeV to avoid large variations in the instrument response. We further apply an energy-dependent and conversion-location-dependent mask (i.e., Front and Back), equal to the 68% inclusion region of the Pass 8R2_V6 PSF, thus excising the region most dominated by the direct emission from the source.

The existence of a bilateral, jet-like feature appears primarily as an excess of quadrupolar power; with a sufficient number of photons it would produce a sawtooth-like structure in the power spectrum in which even multipoles exhibit an excess. Based upon the source SEDs and Monte Carlo modeling of the Fermi source population we have identified sets of Fermi objects that are optimized for the detection of the anisotropic halo feature. These necessarily depend upon the assumed strength of the IGMF, with larger IGMFs producing observable halos in more sources and therefore permitting a larger sample that may be profitably stacked. These optimized source lists may be found in Table 1 of Tiede et al. (2016), with the number of sources ranging from 4 to 18 for current-day strengths of $B_{0}=10^{-17}$ G to $10^{-14}$ G (note that the source-frame field strength is given by $B=B_{0}(1+z)^{2}$, though due to the small redshifts of the gamma-ray sources employed the two differ by at most a factor of 2). For each, an uncontaminated region about each source has been identified via visual inspection. Note that as a result these fields do typically contain many dim sources identified in the Fermi Large Area Telescope Third Source Catalog (3FGL Acero et al., 2015).

Critically, the anisotropy signal is unlikely to be confused with a number of anticipated systematic contaminants. In principle, these may arise from unknown physics at the TeV source, the LAT response, or unmodeled structures in the gamma-ray sky. In practice, chief among the potential contaminants is the degree-scale substructures in Galactic contribution to the diffuse gamma-ray background. To assess the degree to which these contribute to the stacked angular power spectra we now construct a model of the gamma-ray background.

3 Modeling the background contributions to the power spectra

In this section we compute the contributions of the photons from background and contaminating 3FGL sources to the stacked angular power spectra.To assess, quantify, and guide the modeling of said contribution, we first identified 507 fields with Galactic latitudes more than $20^{\circ}$ from the Galactic plane. The locations of these are shown in Fig. 1. Also shown for comparison are the source fields (red). Note that while the background fields appear to nearly cover the sky, the bright point sources in the Fermi sky have effectively been excised according to the criteria used in Tiede et al. (2016) to determine the size of the source fields. That is, all of the bright sources fall into the gaps between the source and background fields. Note that these fields are not devoid of 3FGL sources – sources with fluences less than roughly 40 photons do appear in these background fields and in the source fields. The stacked power spectrum of these 507 fields comprise the background power spectrum that we sought to successfully model.

The model we constructed consists of four components:

3FGL sources 2. 2.

Completion of the 3FGL below the survey sensitivity limit 3. 3.

Large-scale gradients ($>9^{\circ}$) 4. 4.

Small-scale gradients ($\approx 2^{\circ}$)

We will discuss the construction of each individually.

3.1 Inspecting the Fermi sky

The construction of all of these is facilitated by the direct inspection of the Fermi sky, obtained from 239557417 s (the official start of the Fermi data release) to 487012885 s in mission elapsed time (MET). A density map of the 1-100 GeV photons in the Front- and Back-converted Pass 8R2_V6 ULTRACLEANVETO events class is smoothed with a Gaussian kernel. The size of this kernel is informed by typical structure of the point sources, set by the Pass 8R2_V6 PSF, below which no discernible structure appears. Thus, as our finest resolution map we choose a kernel full-width, half-max (FWHM) of $1.2^{\circ}$, which is shown in the top row of Fig. 2. The multitude of point sources that comprise the 3FGL are immediately evident.

3.2 Component 1: 3FGL Sources

Because the period from which the 3FGL was constructed is roughly half that over which these maps are generated, both exposure differences and source variability will produce significant differences in the anticipated number of photons from each 3FGL object. Therefore, we begin by refitting the individual 3FGL source photon fluences.

We do this by comparing the observed gamma-ray maps smoothed on $1.2^{\circ}$ and $4.7^{\circ}$ scales to the expected response from the individual point sources. That is, we begin by computing the observed photon density at each of the 3FGL source locations, given by $n^{1.2}_{j}$ and $n^{4.7}_{j}$, respectively, for the $j$th source. We then simulate the corresponding anticipated fractional response of each individual source at each of the 3FGL positions, i.e., the number density of photons located at the $j$th source location associated with $N_{i}$ photons from the $i$th source are $f_{ji}^{1.2}N_{i}$ and $f_{ji}^{4.7}N_{i}$ for the two smoothing scales, respectively. In doing this we employ the Pass 8R2_V6 ULTRACLEANVETO PSFs and assume the spectral index listed in the 3FGL.

In addition to the point sources in the 3FGL there exists a diffuse component associated with gamma-ray emission from Galactic leptons and hadrons. Unlike point sources, this contribution should be similar at both smoothing scales. Thus, we have

[TABLE]

for some fixed set of diffuse background photon densities, $\bar{n}_{j}$, at both smoothing scales. This may be concurrently solved for the $\bar{n}_{j}$ and $N_{i}$, yielding the set of linear equations for the latter,

[TABLE]

We solve this for the $N_{i}$, and subsequently verify that the implied photon density associated with the diffuse background is positive. While this method for estimating the 3FGL source fluences does account for the differing exposures and source variability, it does not account for the possibility of new sources. It is also poorly suited to correcting the fluences of extended sources. We do not anticipate either of these limits to present significant obstacles.

With the individual source fluences we then construct a realization of the 3FGL, smooth it over the appropriate scale, and subtract it from the observed photon density. This is shown for a smooth scale of $1.2^{\circ}$ for Front- and Back-converted Pass 8R2_V6 ULTRACLEANVETO events in the second row of Fig. 2. Apart from a handful of extended Galactic sources, and up to Poisson fluctuations in the number of sources for the Back-converted map, point sources have been reasonably well removed. In practice, because the background fields specifically exclude bright sources, the details of the 3FGL removal are relatively unimportant. Simply subtracting the 3FGL sources assuming the 3FGL fluxes and fixed exposure produces qualitatively identical results.

3.3 Component 2: Completion of the 3FGL

The second component, the unresolved completion of the 3FGL, can only be inserted in a statistical sense. Ultimately, a complement of sources with the appropriate fluxes are selected and randomly placed within the background/source fields in a manner consistent with the $9.4^{\circ}$-smoothed background. The source fluxes are chosen from the power-law distribution shown in Fig. 3:

[TABLE]

This is supplemented with a detection efficiency, which we assume has a log-normal cutoff, i.e.,

[TABLE]

While for the $|b|>20^{\circ}$ population the threshold flux, $F_{th}$, is well described by a single number, in detail it does depend on Galactic latitude. This is clear in Fig. 4, in which there is a trend in the lower envelope of 3FGL fluxes with Galactic latitude, scaling roughly as $|b|^{-0.33}$. Therefore, we set

[TABLE]

Note that, despite the fact that the assumed flux distribution produces a formally infinite flux (and is thus unphysical), the condition that photons are observed in integer numbers, i.e., at least one, limits the number of photons observed in practice. That is, for a flux $F_{35}$ the probability that at least a single photon will be observed is

[TABLE]

where $\epsilon$ is the exposure. Thus, at small fluxes, the probability of observing at least one photon is linearly suppressed. Therefore, any flux distribution flatter than $F_{35}^{-2}$ produces a convergent number of photons.

In practice, the contribution of the unresolved completion of the 3FGL is negligible.

3.4 Component 3 and 4: Gradient fields

The third component consists of randomly drawing the remaining background photons with a spatial distribution that follows the $9.4^{\circ}$-smoothed map. This modifies the dipole component significantly, though again makes a negligible difference at $m>1$ as shown in Fig. 5.

Finally, the fourth component consists of an attempt to model the degree-scale fluctuations in the diffuse Galactic background. Stacked angular power spectra generated from realizations of the background fields using smoothed maps of the gamma-ray sky are shown in Fig. 5. Unsurprisingly, larger smoothing scales correspond to less power at large $m$. Specifically, at smoothing scales of $1.2^{\circ}$ the power above $m=5$ is noticeably suppressed. Beyond smoothing scales of $2.4^{\circ}$ this suppression extends to $m\geq 3$; by $9.4^{\circ}$ only the dipole power remains. This suggests that any missing power at $m=2$ is most sensitive to a few degree-scale structures in the diffuse gamma-ray background. Therefore, to clearly identify these we generate difference maps between high ($1.2^{\circ}$) and low ($9.4^{\circ}$) smoothing scales. These are shown in the third row of Fig. 2.

Generally, these exhibit a latitude dependent amplitude, falling rapidly with distance from the Galactic plane. This supports their association with Galactic structures. A similar conclusion holds regarding their impact on the stacked angular power spectra. Figure 6 shows a direct comparison between the observed and simulated stacked power spectra, including only the first three components, from the background fields assuming different latitude cuts. When restricted to fields with $|b|>60^{\circ}$ the two are very similar, as anticipated by the large reduction of the Galactic diffuse emission at high latitudes.

In contrast to their absolute amplitudes, the fractional amplitudes in the photon densities are roughly fixed, however (shown in the bottom row of Fig. 2). The PDF of the fractional density fluctuations is well described by the normal distributions, shown in Fig. 7. The mean is offset significantly from zero. This is a natural result from imperfect point source subtraction, likely arising from the failure to subtract unidentified sources. This is supported by the PDF of fractional density fluctuations constructed without any point source subtraction, shown in the dotted blue lines in Fig. 7, which also shows a noted shift to negative fractional flux fluctuations.

The impact of the degree-scale structures in the Galactic diffuse emission is to increase the power at low $m$. We quantify this by generating two stacked angular power spectra for the background fields from over-resolved realizations of the gamma-ray sky that include only the initial three components: 3FGL, its completion, and large-scale gradients (Figure 8). These differ in the maps used to generate the last component, the large-scale gradients; one using the map smoothed on $9.4^{\circ}$ and the other smoothed on $1.2^{\circ}$. The latter map is statistically indistinguishable from the stacked angular power spectrum obtained from the data directly. The difference between the two is entirely due to the otherwise unresolved degree-scale structures of interest.

3.5 Constructing the Background Power Spectra

We model the angular power due to these small-scale background fluctuations statistically as a multipole-dependent, log-normal additional component with means and standard deviations chosen to reproduce the difference in power and difference in the variances of the low- and high-resolution angular power spectra. That is, we assume that the power in the high- and low-resolution angular power spectra, $\mathcal{P}_{m}^{h}$ and $P_{m}^{l}$, respectively, are related via

[TABLE]

where the small-scale structure contribution $\mathcal{P}_{m}^{ss}$ is a log-normal random variable, i.e., $\mathcal{P}_{m}^{ss}=\delta_{m}e^{\sigma_{m}X}$ for a zero-mean, unit-variance Gaussian random variable $X$. The values of $\delta_{m}$ and $\sigma_{m}$ are fully defined by the means and variances of the $\mathcal{P}_{m}^{h}$ and $\mathcal{P}_{m}^{l}$; it is straightforward to show that choosing at each $m$

[TABLE]

where $\Delta_{m}$ is the difference in the mean powers and $\Sigma_{m}^{2}$ is the difference in the variances at high and low resolution at a given $m$, uniquely reproduces the means and variance of $\mathcal{P}_{m}^{h}$ and $\mathcal{P}_{m}^{l}$. This is done independently for angular power spectra constructed using Front-converted events, using Back-converted events, and using both Front- and Back-converted events, resulting in slightly different choices for $\delta_{m}$ and $\sigma_{m}$. For the last case, the resulting $\delta_{m}$ are shown in Figure 8 for the red and blue points. This is compared to the power spectra for the Galactic background of the 3FGL source subtracted Fermi sky (black points).

In summary, we have modeled the background power with a combination of the contaminating 3FGL sources, a random realization of the unresolved completion of the 3FGL source population, the large scale ($>9^{\circ}$) gradients, and random realizations of the small-scale gradients. At very high Galactic latitudes ($|b|>50^{\circ}$) this last component becomes insignificant, strongly implicating the Galactic contribution as its source. The resulting background power spectra is shown in Fig. 9, subdivided into their origins.

4 Results

4.1 Comparing the power spectra to the background

The stacked angular power spectrum for the optimized source set for $B_{0}=10^{-15}$ G is shown in Fig. 10. As anticipated by the power spectrum from the background fields, there is excess power at low $m$. However, this is entirely consistent with that expected for 15 background fields, adjusted via a constant shift for the different Poisson noise level (red points in Fig. 10). Thus, there is no evidence for any excess power at low multipoles.

This is further evident in Fig. 11, which shows the probability distribution of the simulated and observed power about the simulated median. Not only is each value consistent with being drawn from the background, the distributions about the mean is itself similar. Conversely, the ability to reproduce the statistics of the observed stacked angular power provides additional confidence in our ability to adequately model the background model. This remains the case regardless of which optimized source sample is used.

Additionally, we expect different low $m$ to be correlated when pair halos exist (Tiede et al., 2016). While we have made no attempt to quantify the discrepancies, we found that the cross-correlations are consistent with the expected null distribution, and qualitatively look discrepant in the model that accounts for inverse-Compton halos (see Appendix A).

4.2 Quantifying the lack of gamma-ray halos

To quantitatively assess the limits upon the IGMF strength that may be placed by the non-detection of excess structure, we employ our background model to generate multiple realizations of optimized Fermi source lists, with and without inverse-Compton halos present. From these we identify the one-sided confidence level at which a given value of the power may be excluded, i.e., the confidence level (CL), in each direction. Because the stacked powers at different multipoles are strongly correlated in presence of inverse-Compton halos (see Appendix A), we focus on the quadrupole exclusively when reporting constraints. In principle, we can improve upon these by including information at multiple $m$; in practice, this serves to obscure the nature of the test for marginal benefit. Note that because the optimized source lists vary with $B_{0}$, this procedure must be repeated for each IGMF strength.

In the top panel of Fig. 13 we show the 2$\sigma$ (95.5%) and 3$\sigma$ (99.8%) one-sided confidence levels for each multipole with (red) and without (purple) an inverse-Compton halo with $B_{0}=10^{-15}$ G. Even without modeling and including the background in the power spectra, we can exclude a large-scale IGMF with strength $B_{0}=10^{-15}$ G by 2$\sigma$. When the background is included we can exclude our IGMF model at 3.9$\sigma$. Fig. 12 shows this result does not depend on which detector, front or back, is used, suggesting this is not from telescope systematics.

Moving to $B_{0}=10^{-18}$ G, $10^{-17}$ G, $10^{-16}$ G, $10^{-14}$ G, and $10^{-13}$ G, we generally find similar results. For $B_{0}=10^{-18}$ G and $10^{-13}$ G we used the optimized source list described in Broderick et al. (2016) for $B_{0}=10^{-17}$ G and $10^{-14}$ G, respectively (for their power spectra see Fig. 13). These may be excluded at 1.9$\sigma$, 2.4$\sigma$, 4.0$\sigma$, 2.0$\sigma$, and 1.6$\sigma$, respectively. That is, for present-day field strengths between $10^{-16}$ G and $10^{-15}$ G a large-scale IGMF can be excluded at more than 3.9$\sigma$. At greater than 2$\sigma$ present-day field strengths slightly more than $10^{-18}$ G to $10^{-14}$ G a large-scale IGMF can be excluded. These regions are shown in Fig. 14 in comparison to the limits in Ackermann et al. (2018) and Broderick et al. (2018). Beyond those limits the inverse-Compton halos either become sufficiently small that they are confused with the source (small $B_{0}$) or extended that they lie beyond the regions about the sources considered (large $B_{0}$).

5 Conclusions

Combined with the SED lower limits, and making modest assumptions regarding the source duty cycles, we can exclude at more than 3.9$\sigma$ the existence of an IGMF of $10^{-16}<B_{0}<10^{-15}$ G, and at 2$\sigma$ field strengths of $10^{-17}<B_{0}<10^{-14}$, with $\lambda_{B}>100$ Mpc. The constraints presented in this work complement the results of Broderick et al. (2018) and improve the limit of large-scale magnetic fields to be $<10^{-17}{\rm G}$ at 2$\sigma$. However, in this paper we used known gamma-ray bright sources, i.e. FSRQ and BL-Lac objects, from the Fermi 3FGL source list. As a result, while the constraints on the IGMF in this paper are statistically weaker than the Broderick et al. (2018), they require less assumptions about the nature of AGN and their high-energy emission. Namely, the high-energy spectrum for each source is measured and we don’t require AGN unification. Therefore, these results form an independent test of the IGMF when compared to Broderick et al. (2018) and is more direct. Furthermore, we have demonstrated that even if the orientation of the gamma-ray halo is unknown i.e. the magnetic field orientation for each source is unknown, we can still constructively stack sources. This work strongly constrains models of magnetogenesis that invoke inflationary processes in the early universe, which generate fields on large scales.

It is nevertheless worth pointing out that all of the above limits can be avoided if inverse-Compton cascades are preempted by other physical processes. Chief among these are large-scale beam plasma instabilities in the intergalactic medium. The linear growth rates for these far exceed the inverse-Compton cooling times, a conclusion that appears to be robust even when higher order perturbative effects are included (Broderick et al., 2012; Schlickeiser et al., 2012; Chang et al., 2014, 2016). This may have important implications for the gamma-ray sky (Broderick et al., 2014a, b), the thermal history of the intergalactic medium, and cosmic structure formation in low-density environments (Chang et al., 2012; Pfrommer et al., 2012; Puchwein et al., 2012; Lamberts et al., 2015). However, it does remain debated whether these instabilities saturate at sufficiently large amplitudes to efficiently extract the energy of the pairs in the nonlinear regime (Miniati & Elyiv, 2013; Sironi & Giannios, 2014; Chang et al., 2014; Kempf et al., 2016; Chang et al., 2016; Shalaby et al., 2017; Vafin et al., 2018, 2019).

A.E.B., P.T. and M.S. receive financial support from the Perimeter Institute for Theoretical Physics and the Natural Sciences and Engineering Research Council of Canada through a Discovery Grant and Alexander Graham Bell scholarship (PT). Research at Perimeter Institute is supported by the Government of Canada through Industry Canada and by the Province of Ontario through the Ministry of Research and Innovation. C.P., M.S., and E.P. acknowledge support by the European Research Council under ERC-CoG grant CRAGSMAN-646955. PC is supported by the NASA ATP program through NASA grant NNH17ZDA001N-ATP. AL acknowledges support by the Programme National des Hautes Energies (France). E.P. acknowledges support by the Kavli Foundation.

Appendix A Power correlations

Generally the powers at different multipoles are highly correlated when inverse-Compton halos are present (Tiede et al., 2016). As a result, we have not made an attempt to leverage the observed broad inconsistency with the halo model expectation at low $m$. Nevertheless, in Fig. 15 we show the joint probability distributions of the angular power in various low $m$ multipoles in comparison to the observed values. From the comparison when inverse-Compton halos are present (Fig. 15, left) it is clear that the observed powers lie well outside the joint distribution for all instances shown. That is, the observations are even more highly inconsistent than implied by the quadrupolar comparison. In contrast, when inverse-Compton halos are absent, i.e., our null case (Fig. 15, right), the discrepancy disappears.