Observing small-scale \gamma-ray anisotropies with the Cherenkov Telescope Array

Observing small-scale -ray anisotropies with the Cherenkov Telescope Array

[    [

Disentangling the composition of the diffuse -ray background (DGRB) is a major challenge in -ray astronomy. It is presumed that at the highest energies, the DGRB is dominated by relatively few, still unresolved point sources. This conjecture has recently been supported by the measurement of small-scale anisotropies in the DGRB by the Fermi Large Area Telescope (LAT) up to energies of . We show how such anisotropies can be searched for with the forthcoming Earth-bound Cherenkov Telescope Array (CTA) up to the TeV range. We investigate different observation modes to analyse CTA data for small-scale anisotropies and propose the projected extragalactic large-area sky survey as the most promising data set. Relying on an up-to-date model of the performance of the southern CTA, we find that CTA will be able to probe anisotropies in the DGRB from unresolved point sources at a relative amplitude of at energies above and angular scales . Such DGRB anisotropies have not yet been ruled out by the Fermi-LAT. The proposed analysis would primarily clarify the contribution from blazars and misaligned active galactic nuclei to the very-high-energy regime of the DGRB, as well as provide insight on dark matter annihilation in Galactic and extragalactic density structures. Finally, it constitutes a measurement with complementary systematic uncertainties compared to the Fermi-LAT.

a,b,c]M. Hütten c]and G. Maier \affiliation[a]Max-Planck-Institut für Physik, Föhringer Ring 6, D-80805 München, Germany \affiliation[b]Humboldt-Universität zu Berlin, Newtonstraße 15, D-12489 Berlin, Germany \affiliation[c]DESY, Platanenallee 6, D-15738 Zeuthen, Germany \emailAddmhuetten@mpp.mpg.de \emailAddgernot.maier@desy.de \keywordsgamma ray detectors, gamma ray experiments

1 Introduction

The last decade has marked a major breakthrough in -ray astronomy. The Fermi Large Area Telescope (LAT, [1]), operating since 2008, has resolved thousands of astrophysical -ray sources above energies of 100  [2]. Ground-based instruments like H.E.S.S., MAGIC, and VERITAS revealed over a hundred sources emitting -rays in the very-high-energy (VHE) regime beyond 100 .111http://tevcat2.uchicago.edu/ [3] The presence of diffuse Galactic emission along the Milky-Way plane and its origin in cosmic-ray interaction with the interstellar gas and radiation fields in the Galaxy are well established [4, 5, 6]. The radiation processes giving rise to the diffuse Fermi-bubbles [7] and Loop-I structure [8], extending far above the Galactic disk, are still under debate. Nonetheless, the Fermi-LAT has detected millions of GeV photons, which cannot be attributed neither to resolved localised sources or to known diffuse structures [9, 10]: These -rays compose the so-called diffuse -ray background (DGRB)222The DGRB is synonymously denoted as isotropic -ray background (IGRB), extragalactic diffuse -ray background (EDGB), or unresolved -ray background (UGRB) throughout the literature. and, by definition, comprise all -rays outside the Galactic plane with almost isotropic distribution and not assigned to known -ray sources on the sky [11, 12]. Unresolved, far-away extragalactic -ray sources, distributed isotropically on the sky, definitely constitute a major fraction of the DGRB [13, 14, 15, 16, 17, 18]. However, also unresolved sources of Galactic origin at high latitudes [19, 20] and unaccounted large-scale Galactic diffuse emission [12] might contribute to the DGRB. Annihilation of dark matter (DM) particles in the Milky Way halo and throughout the cosmic web might provide an additional, still unidentified -ray source [21, 22, 23]. The relative contributions of numerous Galactic and extragalactic source classes to the DGRB are currently under debate [12].

The prevalence of extragalactic sources in the DGRB is further supported by the recent observation of its spectral softening below [10], in agreement with an attenuation due to pair-production losses from interaction with the extragalactic background light (EBL, [24]). If the DGRB is of mostly extragalactic origin, it must necessarily become more anisotropic at higher energies. Any known physical process generating -rays can only occur in the surroundings of galaxies or galaxy clusters, which appear virtually point-like on the sky for current -ray telescopes. The -rays arriving from these far-away sources are redshifted by the cosmic expansion and attenuated by the EBL, with the latter effect being stronger the higher the -ray energy [25]. For these reasons, only relatively close or spectrally hard objects are expected to emit -rays that are observable at VHE energies, independent of the specific production mechanism at source. This conjecture is in agreement with the identified extragalactic sources in the third Fermi-LAT full-sky catalogue of hard-spectra -ray sources (3FHL) [26] and the angular directions and redshifts of known VHE -rays emitters.

Focusing on such anisotropic signatures of the DGRB to decipher its origin, the approach of an angular power spectrum (APS) analysis was first applied to data from the Fermi-LAT several years ago. An anisotropic component to the DGRB on angular scales (“small-scale anisotropies”) was in fact observed [27]. The amplitude of these anisotropies was found compatible with being constant on all investigated angular scales, suggesting that the APS signal is caused by unclustered point sources.

A recent analysis up to based on more than six years of Fermi-LAT data confirmed this finding and indicates a spectral hardening of the underlying sources with increasing energy [28]. The detected APS is in agreement with the prediction for unresolved blazars [29, 30] and misaligned active galactic nuclei (AGN) [31]. However, it has been suggested that also Galactic sources, in particular high-latitude millisecond pulsars [19], or relic annihilation of Galactic and extragalactic DM [32, 33, 34, 35, 36, 37, 38, 39] might significantly contribute to small-scale DGRB anisotropies.

Further knowledge about the contributors to the DGRB can be obtained by probing its APS in the VHE regime. Due to the large cosmic-ray background, it is highly challenging to measure the absolute intensity of diffuse -rays with Earth-bound -ray detectors. However, the measurement of -ray anisotropies from ground may be much less constrained by the background. It was found by [40] that current instruments and their available data sets provide only limited capabilities for probing VHE anisotropies in the DGRB. However, with the future Cherenkov Telescope Array (CTA, [41]), a relatively large field of view (FOV) of a ground-based instrument will be combined with the best angular resolution ever achieved in -ray astronomy. With CTA, the authors of [40] have shown that a promising sensitivity to small-scale VHE -ray anisotropies can be reached. In the last years, optimization studies for the array layout and detailed studies on the expected performance of CTA have been performed, and dedicated observing plans for the first decade of operation have been drawn [42]. Therefore, this paper presents a refined assessment of CTA’s capability to resolve small-scale anisotropies in the DGRB.

This article is organised as follows: In section 2, we introduce the concept of angular power spectra and describe our likelihood-based analysis method of event data APS. In section 3, we present the instrumental model for the CTA event sampling (§ 3.1) and outline our setup of a CTA extragalactic survey (§ 3.2). We investigate the APS characteristic of the CTA cosmic-ray background and argue to prefer data from a shallow large-area survey for a study of -ray anisotropies (§ 3.3). Section 4 then presents our analysis and results on the CTA sensitivity to small-scale anisotropies. Our findings are finally discussed in section 5 with respect to existing data and expected DGRB anisotropies caused by different source classes. We conclude in section 6.

2 Likelihood-based analysis of event data angular power spectra (APS)

We quantify the anisotropies in the DGRB by decomposing the spatial distribution of -ray-like events into its auto-correlation angular power spectrum (APS), . A square-integrable function on the sphere can be expressed as a linear combination of spherical harmonics ,


where in the following, describes -ray intensities or the dimensionless maps of binned -ray event numbers. The are defined as the covariance of the uncorrelated coefficients ,


For a statistically isotropic field with , the ensemble average of the ,


provides an unbiased estimator for . With the common normalisation it is , and we will denote this dimensional intensity APS in the following with a superscript, . Throughout this paper, we assume the residual cosmic-ray background for Earth-bound -ray observation to be isotropic with no intrinsic power (). However, the shot noise from disjoint total events, binned in pixels, adds a noise power (indicated by the subscript N for noise) to the measurement,


where is the unmasked part of the sky. With this, we reconstruct the full-sky equivalent signal APS, , from a measured APS, , to333Note that the wording “signal” and “noise” are throughout this paper not used in the usual meaning of -ray signal and cosmic-ray background, but to distinguish an intrinsic physical anisotropy (in both - and cosmic rays) against the shot noise from events binned in discrete pixels.


where the power is unfolded by the beam suppression of a radially symmetric point spread function (PSF) [43],


and are the Legendre polynomials of the -th order, and . Finally, we make use of the approximation


to account for the sampling of the APS on a limited sky patch. More precisely, a masking of with a window in angular space results in a convolution in -space [44, 45, 46],


with the convolution kernel ,


At sufficiently large , it is and Eq. 7 holds [45]. However, like in Euclidean Fourier transformation, sharp window edges in angular space cause global spectral leakage artefacts in -space. In § 3.3, we will discuss these effects and the amount of spectral leakage for various CTA field of view shapes.

Due to its large residual background, CTA is unable to determine the absolute level of the diffuse -ray intensity, , which we are ultimately interested in probing for anisotropies. Therefore, we will work with the dimensionless fluctuation APS, . For event maps, the fluctuation APS is connected to the intensity APS via


Under the assumption that consists of two components, , where the background events contain zero intrinsic small-scale anisotropy (despite the shot noise power), the fluctuation APS of the DGRB component can be estimated as


Note that the estimated number of -ray events in the data set must be calculated with some a priori knowledge, e.g, using a spectrum measured by the Fermi-LAT and the expected CTA instrumental response. This is what we will later do to interpret the sensitivity towards anisotropies in CTA data in the context of the DGRB.

We finally use a maximum likelihood (ML) approach to estimate the significance of an APS detection and to fit the signal parameters. We tailor our analysis to the assumption that the DGRB is dominated by unresolved point sources, which show an APS constant in multipole [47].444A power-law scaling of the APS within the accessible multipole range, with can be assumed to probe extended source classes; in particular for the expected APS from Galactic DM subhalos [36, 37, 48, 39, 49]. It is straightforward to correspondingly generalise the here described likelihood analysis, respecting a couple of caveats (e.g., adopting the correct test statistic discussed in the later § 3.3). This constant power is indicated by the subscript P for Poisson power:


We pursue an analysis binned in multipoles . A binned analysis turns out be necessary to suppress remaining masking artefacts in the range of analysis, i.e., a correlation between neighbouring multipoles. We follow the findings by [28] and use the unweighted arithmetic mean in each bin, with . For an isotropic field, Eq. 2, the probability density of follows a distribution and the error on the binned is [43]


where we assume . We find that a logarithmic binning with bins per decade is able to eliminate the multipole correlation due to the masking (which breaks the assumption of Eq. 2) and to assure the correct error estimation according to Eq. 13.

With these definitions, we construct the likelihood function


where denotes the ensemble of reconstructed measured multipoles .

We then maximise the logarithm of the likelihood , Eq. 14, under the constraint and with allowed to vary, and use the ratio of these maximised log-likelihoods as test statistic,


to quantify the significance of some signal being present in the data. Here, and are the ML estimators for the signal APS amplitude and its variance according to Eq. 13; denotes the estimated variance under the constraint of no signal present in the data.

3 CTA instrumental characteristic and large-area survey model

3.1 Instrumental performance of the southern CTA

The Cherenkov Telescope Array (CTA) will be the next-generation Imaging Atmospheric Cherenkov Telescope (IACT) array, consisting of two separate arrays to be erected at the Paranal site (Chile) on the southern hemisphere and on La Palma (Spain) on the northern hemisphere. In this paper, we base the performance of CTA on the recently published prod3b-v1 instrumental performance simulation555http://www.cta-observatory.org/science/cta-performance and solely on the southern array, comprising 4 Large-Size Telescopes (LSTs), 25 Mid-Size Telescopes (MSTs), and 70 Small-Size Telescopes (SSTs). We rely on the southern CTA for the sake of simplicity and to consider the presence of SSTs (which are planned to be part only of the southern array), yielding a higher sensitivity in the TeV regime. In contrast to the earlier prod2 characteristics, the prod3b performance estimation provides reliable results for Cherenkov light hitting the camera planes under an incidence angle off-axis of the telescopes’ pointing direction, crucial for the study of the CTA survey performance. Note that little has changed with respect to the projected on-axis CTA performance () between prod2 (as used, e.g., in our previous study [49]) and prod3b.

The instrumental performance characteristics are the result of extensive Monte-Carlo (MC) simulations. These simulations comprise the modelling of all steps in the Cherenkov light emission and detection, from the shower development in the atmosphere to the response of the photosensors in the cameras. Because of the high computational costs, calculations are done in relatively coarse bins of the parameters (elevation and azimuth angle of the pointing, -ray energy, off-axis angle,…) such that we apply suitable analytical descriptions of the various instrument response quantities fitted to their discrete tabulated values. Such a fitting approach instead of interpolation is necessary for this study, where we intend to avoid artefacts in the APS which result from the numerics of the simulation (and would not be present in a CTA data set). We assume that all observations are taken at a zenith angle of , with equal proportions of southern and northern pointings. Besides the physics in the atmosphere and the instrumental layout, the performance characteristics depend on the chosen analysis cuts. If not explicitly stated differently, all event simulations are based on the CTA performance according to an event selection optimised for of observation and without spatial cuts, other than removing events at from the camera centre.

Figure 1: Effective areas for the prod3b southern CTA off-axis performance at different reconstructed -ray energies . The dashed and dotted curves additionally show the earlier CTA effective area models from [40].

We describe the dependence of the effective -ray collection area, , and the cosmic-ray residual background rate on the angle , using a higher-order Gaussian function:


While previous studies, e.g., [50, 40], have modelled the off-axis behaviour of the CTA effective area and background rate with a Gaussian function (), the modelling by Eq. 16 takes into account the plateau shape of the responses around the pointing position. Figure 1 shows the off-axis dependence of for the southern CTA, comparing the MC results for six bins between and with the fit by Eq. 16. This clearly reveals a deviation from a Gaussian shape, in particular at the lowest energies. Figure 1 also displays the Gaussian CTA effective area model in the earlier study on APS measurements with IACTs [40]. Figure 2 shows our off-axis fitting approach for some exemplary prod3b residual background rates (with the abscissa in log-scale, such that a Gaussian scaling would be represented by a parabola). When speaking about the CTA field of view (FOV) in the following, we refer to the energy-dependent off-axis shape of the residual background rate. We define the FOV radius, , via

Figure 2: Background rates for the prod3b southern CTA off-axis performance at different reconstructed -ray energies .

In table 1, we show the on- and off-axis event rates above various energy thresholds for the CTA instrumental response used throughout this paper.666Note that in [49], we have presented a similar table (Table 2) for the prod2 on-axis rates. It can be seen that the FOV radius according to Eq. 17 increases with energy, which causes the background over the full FOV to decrease more slowly with energy than the on-axis background rate alone. We also display the expected event rate from a DGRB spectrum after [10], using their foreground model B, which leaves the largest fraction of unassociated -rays to the DGRB:


In table 2, we show the CTA instrumental performance according to the previous work from [40], to be compared with our model in table 1. As shown in figure 1, the authors of [40] assumed a much wider FOV (both for the effective areas and background rates), however, with smaller on-axis effective areas. Combined, they obtain rather similar -ray event rates compared to our prod3b model. Also, their assumption of background rates between above has remained valid. However, we will argue in § 3.3 that the non-Gaussian and substantially smaller compared to [40] FOV in the prod3b model provides a major obstacle for analysing deep-field observations for -ray anisotropies.

DGRB event rate Background rate -rays/background
Energy FOV radius on-axis total FOV on-axis total FOV ratio
threshold on-axis total FOV
3.6 66.2
1.1 37.0
0.24 12.4
0.056 5.0
Table 1: DGRB event rates and background rates (without dead time correction), integrated over energies above different lower thresholds up to . The DGRB rates correspond to the DGRB intensity spectrum according to Eq. 18. Events in the “total FOV” comprise all events within the angular direction from the camera pointing position.
DGRB event rate Background rate -rays/background
Energy FOV radius on-axis total FOV on-axis total FOV ratio
threshold on-axis total FOV
0.99 100
0.099 10
Table 2: Model of the CTA performance from [40], according to their tables 1 & 2, , and the assumption of 100/10 Hz background rates above 100/300 GeV. Again, the DGRB rates are computed for the model B from [10]. Note that for a Gaussian off-axis acceptance, .

The finite CTA angular resolution causes a suppression of the APS at angular scales below the resolution, expressed by Eq. 6. The CTA PSF worsens at lower -ray energies and for large angles of the incident -ray in the camera field, i.e. it is . We model as a two parameter King function [51]


owing to the fact that CTA will exhibit a PSF with non-Gaussian tails (see the later figure 4). The two parameters and are again fitted to the MC simulations of the instrumental performance.

3.2 A model of the CTA extragalactic survey

With its large FOV, CTA will be the first IACT array to perform a large-area sky survey [42]. Within the CTA extragalactic survey key science project [42], it is planned to observe of the sky outside the Galactic plane within the first decade of observation. For a total observation time of about with both arrays, such a survey is projected to reach an average sensitivity to fluxes greater than ( the flux of the Crab Nebula) above for point sources with a spectrum similar to the one of the Crab Nebula [42]. In this paper, we adopt the same survey field as done in [49] for the extragalactic survey, namely, a circular region around the Galactic south pole, i.e., . As in [49], we assume that the whole area is covered by the southern array in , while [42] propose to raster of the survey field in with the southern array, and the remaining area in with the northern CTA. In this work, we consider a system dead time of 5%, reducing our effective total observation time to . Note that this choice is fairly conservative, and CTA is targeted to reach system dead times smaller than . On the other hand, we do not consider a loss of observation time due to the telescope slewing between each pointing.

The FOV size of the instrument is critical for the average survey exposure on each spot on the sky and the homogeneity of the exposure. The exposure homogeneity is in particular dependent on the pattern of the pointing distances of the individual survey observations. For the pointing pattern, we adopt a grid relying on the HEALPix pixelisation scheme [52]. The HEALPix tessellation facilitates equally spaced pointing positions on the sphere, where the sphere curvature becomes significant for an area as large as . This large-area survey telescope spacing strategy is illustrated in figure 3. In the remainder of this work, we consider two different grid spacings: A HEALPix grid with results in an average distance between the telescope pointings of .777 is defined as the square root of the “pointing pixel” area with size . Hereby, is covered by 3136 individual pointings with an observation time of each. This results in an average on-axis equivalent exposure of above with a homogeneity of . In contrast, a HEALPix grid with is considered, yielding with 12,416 single observations over each, and the same average exposure above with .

Because of the lower sensitivity of the northern CTA above 100  compared to the southern array, our simplified large-scale survey setup finally serves as a fairly realistic description of what can be achieved with CTA according to [42]: Using the cssens ctool [51], we find for a dead time corrected on-axis observation with and our selected CTA instrument response a survey sensitivity of of the Crab nebula flux above , not more than better than envisaged by [42].888Requiring a test statistic of and without applying trials corrections.

(a) ()
(b) ()
Figure 3: Large-area survey tiling strategy on in the HEALPix scheme as adopted in this work. For display purpose only, we show a spacing of () in the upper panel. In our default survey setup, we use (lower panel) and compare to . The green quadrilaterals show the HEALPix pixels in which the individual telescope pointings are centred. The red circles have a diameter of and show the CTA field of view diameter at about . Note that all HEALPix pixels cover the same solid angle, however, they are not congruent.

As the angular resolution degrades with offset from the camera centre, we also have to consider an average angular resolution whose homogeneity is, like the exposure, affected by the observation pattern. We approximate an average PSF under the assumption that each spot on the sky is equally observed under all incidence angles from the camera centre up to . In a finite energy interval , the effective PSF additionally scales with the energy distribution of events from the -ray intensity spectrum . Altogether, we model the effective CTA survey PSF in a finite energy interval as:


The -ray spectrum is chosen according to the source class hypothesised in the analysis. Figure 4 shows the averaged PSFs, , and their corresponding APS attenuation factors at different energies.

The PSF modelling determines the highest multipole , up to which we perform our analysis to avoid the impact of an increasing systematic uncertainty about the PSF at larger . This choice of together with our chosen map resolution also allows us to ignore the power suppression of the data due to the finite bin size (“pixel window”, see right panel of figure 4). Finally, we do not consider the finite energy resolution and energy bias of the CTA instrument in this study. As we investigate the CTA performance in only very coarse energy intervals, the energy resolution of the instrument is negligible.

Figure 4: Left panel: Southern CTA angular resolution in terms of the PSF averaged over , in angular space at different energies. Right panel: Corresponding window functions in multipole space. A Gaussian curve and its multipole transformation are given for comparison (dashed lines). Also, the pixel window function of a finite HEALPix grid with resolution is shown (dotted line on the right).

3.3 APS characteristic of the cosmic-ray background in CTA observations

As shown in table 1, CTA suffers from a large irreducible cosmic-ray background. This background constitutes the dominant challenge of power-spectral methods for Earth-bound -ray detectors, in contrast to the similar analyses already performed with Fermi-LAT data.

Figure 5: APS of the events from an intrinsically isotropic background in our CTA survey field on and for an observation spacing of . The APS is shown for four different energy thresholds. The dashed vertical line indicates one example of our lower cuts in -space. The figure shows the mean and the 68% credible intervals (shaded bands) based on 50 simulations each. The grey-dotted lines mark the estimator for each energy threshold.

To study the APS characteristic of this irreducible background, we use the gammalib and ctools frameworks [51]999http://cta.irap.omp.eu/gammalib/ and http://cta.irap.omp.eu/ctools/ to simulate hundreds of samples of the background events in our CTA extragalactic survey setup with , based on the background rates as presented in § 3.1. Subsequently, we collect the sampled events in maps of spatial bins in the HEALPix scheme, and compute the APS of the binned event maps. Figure 5 presents the resulting fluctuation APS above different energy thresholds. The shot noise estimator, Eq. 4, is indicated by the dashed horizontal lines in figure 5. It can be clearly seen that the window of the quarter-sky survey field dominates the multipole range at . According to section 2, we assume (i) the approximation to be valid in a regime unaffected by the window above some , (ii) the same for the error described by Eq. 13, and (iii) a signal hypothesis with only one additional free parameter in the maximisation of the likelihood Eq. 14, namely, . If all these conditions are satisfied, the distribution of TS values from many samples of our background maps must follow a distribution with one degree of freedom [53],


Requiring the reproduction of this statistic, we determine


as lower limits in -space for our analysis.101010For the integrated sensitivity, we use the most conservative . As mentioned in section 2, we find that a binning of the signal and its error in multipoles is crucial to reproduce Eq. 21. This procedure worked well for all energy up to , whereas at higher energies, window artefacts can no longer be removed by binning for . For these highest energies, we use empirically obtained from our MC simulations, but we emphasise that the results for this bin should be treated with some caution. The cuts finally define the angular scales below which CTA will be able to probe anisotropies.

Possible deviations of the available data set from the standard survey setup considered in the previous paragraph alter the background APS. Therefore, we have studied (i) the presence of 19 circular exclusion regions with a Gaussian mask and (corresponding to the 19 up-to-date known VHE -ray sources in the chosen survey field), (ii) a coarser survey grid with , and (iii) a 10% variation of the background rate between the individual observation pointings in the survey due to different weather conditions and calibration uncertainty. The latter is a rather conservative estimate, as CTA is required to provide a much more stable data rate over time on all timescales. Figure 6 compares the absolute (left) and relative (right) APS of these varied setups to our standard case from figure 5. The dashed vertical lines in these panels show the satisfying the background hypothesis in our standard survey setup. We conclude that neither the exclusion of some dozens of regions in the survey nor a varying event rate significantly pollute the APS above the threshold according to Eq. 22. However, a coarser survey pattern with in fact adds significant contamination to multipoles above .

Figure 6: The impact of various deviations from the standard survey simulation onto the residual background APS of events with . The relative change of power at each multipole compared to the spectrum presented in figure 5 is shown on the right.
Figure 7: APS of all background events above in a single CTA pointed observation over (for an instrumental response optimised for 50  of observation). We also show the APS for the observation time distributed over four observations, pointed offset from the central pointing position (“wobble observation mode”). The dotted line shows correspondingly . The vertical dashed line denotes the cut applicable for the quarter-sky survey APS (figure 5).
Figure 8: APS of all background events above towards eight different positions in the sky and a total observation time of (for an instrumental response optimised for 50  of observation). We also show the APS for these observations combined with the events from a quarter-sky survey.

The authors of [40] also investigated whether a single or some few deep-exposure observations with IACTs can be used for an APS analysis. For their simplified Gaussian model of the CTA field of view (see figure 1), they found a slight increase in sensitivity when investing a fixed total observation time in a single deep field instead of distributing the time over observations of several fields. However, we argue that small field observation data may hardly compete with large-area survey data for a realistic, non-Gaussian CTA instrumental model: Figure 7 shows the APS of the residual background above in a deep exposure () towards a single spot in the sky, where . The figure reveals dominant lobes in the spectrum which are not present in the multipole transformation of a simple Gaussian function [40]. The dashed vertical line again shows the corresponding of our standard survey setup. Adopting a small set of several deep-exposure data sets does not much improve the situation, as we show in figure 8. Here, are distributed (in wobble observation mode) over eight pointings on the sky. The figure also shows that combining deep-field observations with survey data does not attenuate the multipole lobes. For uncorrelated data sets, the joint APS of added data sets is the sum of the individual spectra, and the oscillations are preserved.

In principle, for known mask shapes, the spectral leakage can be eliminated from the spectrum by a rigorous calculation of from Eq. 9, as done by, e.g., the PolSpice [54, 55] or Master [56] algorithms. However, even provided a perfect knowledge of the masking window in a given energy interval, the unmasked spectrum can only be reconstructed at the expense of noise amplification [27]: For a mask being orders of magnitude larger than the physical signal in space, any information about the signal is likely to be buried in the shot noise and the systematic uncertainty about the mask. However, a rigorous investigation of unmasking small field-of-view -ray data is still to be done, e.g., whether the APS lobes can be suppressed by a suitable apodisation or tapering of the data. Such an advanced study would be particularly helpful to assess the usage of data from deep observations of dark spots (as, e.g., foreseen for the search for dark matter annihilation in dwarf spheroidal galaxies) for an APS analysis. For the remainder of this paper, we will focus on the APS analysis according to section 2, applied to a wide-field CTA extragalactic survey with .

4 CTA sensitivity to anisotropies in the DGRB

4.1 Sensitivity analysis description

We calculate the sensitivity to small-scale anisotropies in a CTA large-area survey data set at the 95% confidence level (C.L.) detection threshold. This C.L. was also used by [40]. Therefore, we proceed as follows: We first generate mock skymaps, each containing point sources with equal flux level, , randomly distributed on the quarter-sky survey field. With [47], such maps have an intrinsic anisotropy of constant over all multipoles, where is the intensity arising from these fluxes. Consequently, these sources generate a fixed fluctuation APS of . The flux of all the sources is modelled to follow the DGRB spectrum, Eq. 18, . We use the ctobsim ctool to simulate skymaps of the -ray and background events in the survey field and we increase the normalisation of until we detect an anisotropy signal, , in the mock data at the C.L. with the likelihood test described in section 2. According to Eq. 21, a C.L. detection corresponds to a value of our test statistic. We call the number of -ray events obtained from at this sensitivity threshold . We repeat the calculation 25 times, varying the source positions on the sphere.111111Note that by this procedure, we do not vary the intrinsic anisotropy of a constant intensity, but vary the intensity of a constant fluctuation anisotropy. By this, the total mean number of events, , is not constant in our search for the intensity threshold. However, as and moreover in the scanned range of anisotropy levels (starting from ), this difference in the approach is negligible.

We perform our analysis in the four energy bins between 30 and 100 GeV, 100 and 300 GeV, 300 GeV and 1 TeV, and between 1 and 100 TeV; additionally we report the integrated sensitivity using all events above 30 GeV. As we collect all events in these energy bins, disregarding their individual energies, our results are virtually independent of the exact shape of the used DGRB spectrum. The only minor dependence of our analysis to the hypothesised DGRB spectral shape enters via the modelling of the average angular resolution in each energy bin, Eq. 20. However, we consider this dependence negligible for the rescaling of the results between different DGRB spectra in the next subsection.

Figure 9: Illustration of the APS likelihood fit in the energy bin . Left: Global spectrum for an extreme fluctuation of the DGRB. It can be seen how the DGRB events (orange curve) add power on the top of the shot noise from the cosmic-ray background (blue), which however is obfuscated in the data by the sky mask at (black curve). At , the instrumental angular resolution suppresses the intrinsic signal power to the overall shot noise level. We evaluate the spectrum at . The green curve shows the likelihood fit to the data. It can be seen that the spectrum of simulated -ray events alone (orange curve) is well estimated and recovered by the fit. Right: Closeup of the interval for a signal APS at the 95% C.L. sensitivity threshold, with the average shot noise power already subtracted and the beam suppression undone.

Figure 9 illustrates the analysis and the equations from section 2. In the left panel, we show an APS containing -ray and background events, however, for an extreme anisotropy level for illustration purpose only. To scrutinise our analysis procedure, we also show a “pure” spectrum from -ray events only with no cosmic-ray background (orange curve). This information is of course not available in a later analysis of real data. In the right panel, we show the APS after subtracting the shot noise and after unfolding (“signal APS”). The blue band indicates the expected unbinned error, for the background-only hypothesis. This variance band is based on the scaling empirically derived from a set of MC simulations.121212See [57] for a comprehensive study of the APS variance of non-Gaussian fields. The black crosses in the right panel show the signal APS with the error, Eq. 13, after binning and maximising the likelihood. Finally, the green curve shows the likelihood fit to the binned and in both panels well agrees with the anisotropy of the injected -ray events (orange curves).

4.2 Results

In table 3, we present the 95% C.L. sensitivity to small-scale anisotropies in a CTA large-area survey,


In the top panel, we show the anisotropies in the total data set (including all background and DGRB events) to which CTA is sensitive to, independent of the origin of these anisotropies. While we have carefully excluded expected artefacts from the analysis in § 3.3, in principle, if any instrumental or physical origin generates anisotropies of these magnitude and on scales below , they are seen by our analysis. We will later discuss in section 5 how this sensitivity can be interpreted with respect to small-scale cosmic-ray anisotropies.

Sensitivity to anisotropies in all events of the survey
Energy “Measured” total Fluctuation sensitivity,
interval events,
Sensitivity to anisotropies in a presumed DGRB intensity (top panel rescaled)
Energy Presumed events Sensitivity to Presumed events Sensitivity to
interval , Eq. 18 , Eq. 24
Table 3: Top panel: 95% C.L. sensitivity to any anisotropies in all events of a CTA extragalactic survey data set. We report the median sensitivities and, because of our analysis’ limited statistic, a symmetric 68% C.I. Bottom panel: Sensitivities to small-scale anisotropies in the DGRB for an exponentially suppressed DGRB, Eq. 18 (left columns) or a broken power-law extrapolation Eq. 24 (right columns).

In a second step, we assume that all anisotropies are generated by the -rays of the DGRB. Taking the results from the top panel of table 3 and some additional knowledge about the average DGRB intensity, Eq. 11 can be used to draw a statement about the relative fluctuation level in the DGRB detectable by our analysis. This is presented in the lower panel of table 3. Here, we show this sensitivity assuming an exponentially suppressed DGRB intensity, Eq. 18, on the left. While [10] have excluded a pure power-law scaling of the DGRB in the VHE regime at the level, they cannot distinguish between an exponential cut-off of the DGRB (suggested by the absorption of distant -rays on the EBL) and a continuation of the spectrum with a steeper power-law slope (which would hint to a different source class at the highest energies). Therefore, we additionally investigate our results for a broken power-law DGRB,


with , , , and obtained from a fit to the Fermi-LAT data from [10, model B]. The APS sensitivity for such a broken power-law DGRB intensity is shown in the right columns of table 3. It can be seen that a power-law behaviour of the DGRB is needed in the TeV regime for generating detectable anisotropies, above . The CTA sensitivity assuming a broken power-law extrapolation of the DGRB is finally also displayed in figure 10. We stress that an interpretation of a detected anisotropy power in terms of the DGRB relies on complementary measurements of its diffuse intensity and the corresponding expected number of CTA events.

Our results of the CTA sensitivity can be compared to the previous study from [40] using . The authors of [40] fixed , and expressed their sensitivity in terms of the detectable fraction of dark matter induced anisotropies. With a detectable , we find our result comparable to their most conservative assumptions about the CTA performance ().131313Independent of the presumed DGRB spectrum, the CTA configuration , an observation time of , and background rate from [40] yields the sensitivity to in all events, most similar to our analysis ( between 300  and 1 , see upper panel of table 3).

5 Discussion

DGRB anisotropies measured by the Fermi-Lat:

Figure 10: Sensitivity of the CTA extragalactic survey to anisotropies in the DGRB compared to the measurement by the Fermi-LAT [28]. The figure shows the relative anisotropy levels, , with respect to the broken power-law DGRB intensity, Eq. 24. The highest energy bin for the CTA sensitivity shows the sensitivity integrated up to .
Figure 11: Comparison of the binned APS in the three overlapping bins between the Fermi-LAT analysis [28] (masking 3FGL sources) and our CTA sensitivity calculation at the 95% C.L. As in figure 10, we show the fluctuation APS with respect to the DGRB intensity of Eq. 24. Note that the CTA points and the fit in the upper right panel are the same as shown in figure 9 (right). Also note the larger vertical axis scale in the lower right panel.

Small-scale anisotropies in the DGRB have already been detected in the data from the Fermi-LAT. Taking 81 months of reprocessed Pass7 data between 500  and 500 , the authors of [28] have recently updated the earlier result from [27], based on 22 months of Fermi-LAT data. In addition to covering a larger energy range and better precision compared to the first measurement, [28] conclude on two power-law populations of point sources dominating the DGRB anisotropies below and above 2 , respectively. To compare with our results, in figure 10, we show the Fermi-LAT measurement of the auto-correlation APS in terms of the fluctuation APS, , relative to the average DGRB intensity of Eq. 24.141414Note that in figure 10, the relative location of the measured data points and the CTA sensitivity in terms of the fluctuation APS is independent of the assumed . The authors of [28] studied the APS for two slightly different data sets, excluding all 3FGL sources (orange crosses in figure 10) and only excluding the 2FGL sources (blue crosses). As expected, the DGRB power is reduced when masking the power originating from additionally resolved sources in the 3FGL. From the comparison between the Fermi-LAT measurement and our results, it becomes evident that the existing detection of DGRB anisotropies from [28] already weakly excludes the small-scale DGRB anisotropies accessible with CTA. However, the Fermi-LAT measurement is based on a slightly lower multipole range (). Also, the statistical uncertainty of the Fermi-LAT measurement in its highest energy bins is large; and in the regime of saturation of the Fermi-LAT detector, systematic biases may additionally increase the uncertainty. This is illustrated in figure 11 by a more close comparison between the three overlapping energy bins of the Fermi-LAT measurement [28] and our CTA sensitivity. Notably, a potential oversubtraction of the mask in the lowest -bins in the Fermi-LAT measurement may artificially reduce the detected signal. We also remark that the energy APS from [28] shows low-significant features which are not yet fully understood, e.g., a peak of the angular power in the single energy bin between 50 and 95.27  when using Pass8 data (which is not present for the default analysis relying on Pass7 data). Therefore, we conclude that the reach of small-scale -ray anisotropies with CTA is not ultimately excluded by Fermi-LAT results, in particular at the highest energies above 300 , and a complementary analysis with CTA may help to reduce instrument-related systematic uncertainties.

Expected VHE anisotropy levels of various source classes:

Several estimates have been drawn to assess the expected anisotropy level from unresolved members of various source classes contributing to the DGRB [19, 31, 49, 29, 36, 58]. While this has been already extensively discussed by [27, 12], we present a short overview of these estimates together with our calculation of the CTA sensitivity in table 4. For the APS signal from Galactic DM substructure, we rely on our previous work [49], where we have comprehensively bracketed the expected intensity and APS characteristics from Galactic DM. While we have presented in figure 7 of [49] the dimensional intensity APS for a particular particle physics model, we report in table 4 with the average emission from all Galactic DM (smooth halo and subhalos) in our definition of the CTA extragalactic survey field.

Using these estimates, the second to last column in table 4 indicates the fraction of intensity of the respective source class to the DGRB to create anisotropies detectable with CTA. This is confronted with a coarse estimate of the predicted DGRB fractions (last column of table 4). It is worth noting that [17] claim that in the VHE range, , the DGRB is dominated by the emission from unresolved blazars with , although their anisotropy imprint may be too low to be detectable with CTA. A signal from Galactic DM substructure would only bear the chance of a detection for the most optimistic clustering model from [49] and the full DGRB (at the peak of the DM-induced -ray spectrum) originating from DM annihilation. This would apply to annihilation cross sections of for DM particle masses and annihilation into bottom quarks. Although this suggests that the indirect search for DM signals with CTA via -ray anisotropies is not very competitive, we stress that this conclusion only applies for the standard paradigm of weakly interacting massive dark matter particles in a CDM Universe and only undergoing DM self-interactions. Moreover, a dedicated study tailored to spatially extended Galactic and extragalactic DM structures is excluded from this work. Overall, we finally conclude that an APS analysis with CTA is most promising to constrain the contribution of misaligned AGN (as already pointed out by [31]) and possibly also unresolved blazars to the DGRB in the VHE regime.

Source class Predicted Detectable Predicted
Millisecond pulsars [19] 38% [20]
Misaligned AGN () [31] [14]
Galactic DM substructure (HIGH) [49, this work] depending on
Unresolved blazars [29] [17]
Galactic DM substructure (LOW) [49, this work] depending on
Extragalactic DM structure [36] depending on
Star-forming galaxies [58] [13]
Table 4: Predicted contributions to the DGRB of various unresolved source classes. If not explicitly stated, all quoted values for and are assumed constant in the range. For DM annihilation, the predicted intensity depends on the velocity-averaged annihilation cross-section, .

Are small-scale cosmic-ray anisotropies excluded?

Throughout this paper, we have assumed a perfectly isotropised CTA residual background of electrons and misclassified hadrons. Due to the large relative amount of these background events, however, even tiny amplitudes of in the background could dominate over any -ray anisotropies. In fact, cosmic-ray anisotropies at scales are since long known [59, 60], with hadronic cosmic-ray anisotropies at the level of up to [61, 62]. However, even when extrapolating this level beyond , given that less than 1% of all hadrons survive the CTA -hadron separation cuts, remaining hadronic cosmic-ray anisotropies in the CTA background are smaller than and thus negligible.

The Fermi-LAT collaboration has recently published a search for anisotropies in the electron intensity based on almost seven years of data [63]. No electron anisotropy has been found, with an energy-dependent upper limit on the dipole anisotropy of . As their analysis suggests a sensitivity to at and ,151515See their figure 11 in the supplementary material to [63]. it remains to comprehensively study whether CTA’s sensitivity to small-scale anisotropies may supersede Fermi-LAT’s ability to measure electron anisotropies at .

However, although low-amplitude small-scale electron anisotropies may be within CTA’s sensitivity range, we judge their existence to be highly unlikely. Only primary cosmic electrons, accelerated in point-like sources and with sufficiently large Larmor-radii, are able to generate small-scale anisotropies. Because of their fast synchrotron cooling, such primary electrons must originate from nearby () sources, and only a few supernova remnants constitute plausible candidates [64]. Also, the scenario of DM annihilation in nearby subhalos generating electron anisotropies is considered to be already ruled out [65, 66].

6 Conclusions

In this work, we have studied whether CTA is able to probe small-scale anisotropies in the DGRB, based on the latest knowledge about the future instrument. Using standard background rejection cuts and a binned likelihood analysis, we have found a sensitivity to -ray anisotropies with relative amplitudes of above 30 GeV and scales smaller than . Beyond , such anisotropies may only be detectable if the DGRB is not exponentially suppressed at sub-TeV energies. While small-scale anisotropies are most likely solely generated by -rays, our results can also be interpreted in terms of any anisotropies present in the data, to which our results show a sensitivity to .

We have confronted this sensitivity to previous detections of DGRB anisotropies in the Fermi-LAT data. CTA may be in reach to detect the anisotropy levels found at low significance by the Fermi-LAT, in particular at -ray energies above 300 . A complementary analysis of CTA data may also help to rule out systematic biases affecting the current and future [67] Fermi-LAT measurements. Discussing various classes of unresolved astrophysical -ray emitters, a CTA analysis for VHE -ray anisotropies will be primarily able to probe current models for the population of unresolved blazars and misaligned AGN.

Investigating various anticipated CTA data sets, we have found the CTA extragalactic survey to be most suitable for this kind of analysis; and the analysis presented in this paper proposes an additional way of studying the survey data. On the other hand, we have stated that the analysis of single-FOV data for a non-Gaussian acceptance of the instrument is heavily impeded by artefacts in the multipole transformation. Discussing various sources of pollution to an ideal survey data set, we have obtained that a homogeneous exposure is most crucial for the success for the search of small-scale anisotropies. We therefore recommend a survey tiling as small as or even a slew survey, where data is taken continuously. Whether the sensitivity of CTA to -ray anisotropies can be further improved by such a slew survey (at the potential expense of a worse angular resolution), using divergent telescope pointing, dedicated analysis cuts, or a sophisticated modelling of the background APS, is finally left for future work.


This work was conducted in the context of the CTA Extragalactic and Fundamental Physics Working Groups and has gone through internal review by the CTA Consortium; we thank the internal reviewers Josep Martí and Germán Gómez-Vargas for their valuable remarks. Also, we thank Elisa Pueschel for her comments that helped to substantially improve the quality of the manuscript. This research has made use of the CTA instrument response functions provided by the CTA Consortium and Observatory, see http://www.cta-observatory.org/science/cta-performance (version prod3b-v1) for more details. This work was supported by the Research Training Group 1504, “Mass, Spectrum, Symmetry”, of the German Research Foundation (DFG) and the Max Planck Society (MPG). Some of the results in this article have been derived using the HEALPix package [52].


  • [1] W. B. Atwood, A. A. Abdo, M. Ackermann, W. Althouse, B. Anderson, M. Axelsson et al., The Large Area Telescope on the Fermi Gamma-Ray Space Telescope Mission, ApJ 697 (June, 2009) 1071–1102, [0902.1089].
  • [2] F. Acero, M. Ackermann, M. Ajello, A. Albert, W. B. Atwood, M. Axelsson et al., Fermi Large Area Telescope Third Source Catalog, ApJS 218 (June, 2015) 23, [1501.02003].
  • [3] S. P. Wakely and D. Horan, TeVCat: An online catalog for Very High Energy Gamma-Ray Astronomy, International Cosmic Ray Conference 3 (2008) 1341–1344.
  • [4] W. L. Kraushaar, G. W. Clark, G. P. Garmire, R. Borken, P. Higbie, V. Leong et al., High-Energy Cosmic Gamma-Ray Observations from the OSO-3 Satellite, ApJ 177 (Nov., 1972) 341.
  • [5] C. E. Fichtel, R. C. Hartman, D. A. Kniffen, D. J. Thompson, G. F. Bignami, H. Ögelman et al., High Energy Galactic Gamma Radiation Observed by the SAS-2 Satellite, International Cosmic Ray Conference 1 (1975) 29.
  • [6] M. Ackermann, M. Ajello, W. B. Atwood, L. Baldini, J. Ballet, G. Barbiellini et al., Fermi-LAT Observations of the Diffuse -Ray Emission: Implications for Cosmic Rays and the Interstellar Medium, ApJ 750 (May, 2012) 3, [1202.4039].
  • [7] M. Su, T. R. Slatyer and D. P. Finkbeiner, Giant Gamma-ray Bubbles from Fermi-LAT: Active Galactic Nucleus Activity or Bipolar Galactic Wind?, ApJ 724 (Dec., 2010) 1044–1082, [1005.5480].
  • [8] J.-M. Casandjian, I. Grenier and for the Fermi Large Area Telescope Collaboration, High Energy Gamma-Ray Emission from the Loop I region, ArXiv e-prints (Dec., 2009) , [0912.3478].
  • [9] A. A. Abdo, M. Ackermann, M. Ajello, W. B. Atwood, L. Baldini, J. Ballet et al., Spectrum of the Isotropic Diffuse Gamma-Ray Emission Derived from First-Year Fermi Large Area Telescope Data, Phys. Rev. Lett. 104 (Mar., 2010) 101101, [1002.3603].
  • [10] M. Ackermann, M. Ajello, A. Albert, W. B. Atwood, L. Baldini, J. Ballet et al., The Spectrum of Isotropic Diffuse Gamma-Ray Emission between 100 MeV and 820 GeV, ApJ 799 (Jan., 2015) 86, [1410.3696].
  • [11] C. D. Dermer, The Extragalactic -ray Background, in The First GLAST Symposium (S. Ritz, P. Michelson and C. A. Meegan, eds.), vol. 921 of American Institute of Physics Conference Series, pp. 122–126, July, 2007. 0704.2888. DOI.
  • [12] M. Fornasa and M. A. Sánchez-Conde, The nature of the Diffuse Gamma-Ray Background, Phys. Rep. 598 (Oct., 2015) 1–58, [1502.02866].
  • [13] I. Tamborra, S. Ando and K. Murase, Star-forming galaxies as the origin of diffuse high-energy backgrounds: gamma-ray and neutrino connections, and implications for starburst history, J. Cosmology Astropart. Phys. 09 (Sept., 2014) 043, [1404.1189].
  • [14] M. Di Mauro, F. Calore, F. Donato, M. Ajello and L. Latronico, Diffuse -Ray Emission from Misaligned Active Galactic Nuclei, ApJ 780 (Jan., 2014) 161, [1304.0908].
  • [15]