AGN in the Extragalactic Gamma Ray Background : A Possible Detection Method for Unresolved AGN Based on Intrinsic Variability in Sources Observed by Fermi-LAT

(1)

MSc Astronomy and Astrophysics

GRAPPA track

AGN in the Extragalactic Gamma Ray Background

A Possible Detection Method for Unresolved AGN

Based on Intrinsic Variability in Sources Observed by Fermi-LAT

Master Thesis

David Homan

5670268

60 ECTS

July 2012 - July 2014

Supervisor:

Daily Supervisor:

Examiner:

(2)

Abstract

The precise nature of the extragalactic gamma ray background (EGRB) is still unclear. It is expected that gamma ray emission by unresolved Active Galactic Nuclei (AGN) forms an important contribution. AGN are among the most lumi-nous sources in the Universe and many AGN show a variability in their gamma ray emission. The method presented in this thesis will use this variability to attempt to discern the signal of unresolved AGN from the constant EGRB. The Fourier power spectrum, characterising the variability of a signal, is calculated for a sample of 64 bright AGN and 1000 background regions. The spectra are based on lightcurves measured by the Fermi Large Area Telescope (LAT), cov-ering 268 weeks of observations. The flux data is divided into weekly bins. The difference in variability between the AGN signal and that of the background is limited to the low frequency range of the spectra and is summarised in the ‘variability score,’ a sum over the frequency bins 1-8. This test statistic is then used to simulate the detection of unresolved AGN. The variability score of a region containing a simulated AGN is compared with that of the empty back-ground. If there is a 5σ certainty that the signal belongs to an AGN, a detection is registered. The simulation provides a decidedly negative result regarding the viability of this method: no detections were made. Although it is possible to refine the comparison with the background, it appears that AGN variability is not distinct enough to be used as a detection tool.

(3)

1. Introduction

The research presented in this thesis is an investigation into the possibility of using Fourier analysis in the reduction of the diffuse gamma ray background (DGRB). The focus will be on one particular component of the DGRB, the con-tribution of unresolved Active Galactic Nuclei (AGN). As most AGN display some intrinsic variability, this property could be used to detect individual AGN whose flux lies below current detection limits. The idea follows from a simi-lar detection method for the contribution of unresolved pulsars to the DGRB, suggested by Geringer-Sameth and Koushiappas (from here on referred to as GSK) (Geringer-Sameth and Koushiappas 2012). The procedure outlined here will take a more detailed approach in determining of the relevant test statistic than proposed by GSK. The data used in this research is collected by the Large Area Telescope (LAT) aboard the Fermi Gamma-ray Space Telescope.

The development of the detection method is largely reflected in the outline of this thesis. The Introduction will discuss the astrophysical properties of the ob-jects of study and will attempt to position the topic of this thesis within the field of high energy astrophysics. The first object to be discussed is the gamma ray background itself: after all resolved sources have been removed from the gamma ray sky, there is a remaining flux, which is referred to as the Diffuse Gamma Ray Background. This background contains a large galactic component of dif-fuse emission from the Milky Way. After this component has also been removed, the remainder is the so-called Isotropic Galactic Background Radiation (IGRB), which contains a possible irreducible local and a large galactic component.

It is generally believed that AGN form a large component of this isotropic background as they are among the most powerful sources of radiation in the Universe and are strongly represented among gamma ray sources. There are most likely also contributions from other undersolved sources, such as Gamma Ray Bursts, normal star forming galaxies and star burst galaxies. A possi-ble diffuse component to the background is the emission due to Dark Matter self-annihilation or decay. The exploration of the gamma ray background there-fore provides possible insights for both the study of Dark Matter and of AGN. Section 1.1 will discuss the possible contributions to the IGRB as well as the methods that have been used to study the entire gamma ray background.

Because the AGN contribution is the topic of research in this thesis, the Introduction will conclude with a discussion of AGN. AGN astronomy is a con-tinually developing field and there is still much discussion about the emission process, the interpretation of AGN energy spectra and the powering mecha-nisms. Section 1.2 will cover the prevailing AGN emission model and its corre-spondence to blazar spectra, the classification and possible evolution of different blazar types and the variability of AGN.

The development of the detection method itself will begin in section 2, with the description of the Fourier techniques necessary to properly analyse the vari-ability of a given signal. Specifically, this will include the definition of the Fourier Transform, the associated Discrete Fourier Transform (DFT) which must be

(5)

ap-1 INTRODUCTION 1.1 The Gamma Ray Background

plied when dealing with binned data and a discussion of the DFT response to astrophysical signals. The section will conclude with the specific Fast Fourier Transform (FFT) algorithm that was used to analyse the Fermi data.

A brief description of the Large Area Telescope will follow in section 3. The LAT is the main instrument aboard the Fermi Gamma-ray Space Telescope. It is a high energy gamma ray telescope with the broad scientific objective of furthering our understanding of the gamma ray sky. The satellite was launched on June 11, 2008 and the LAT has been gathering data for over five years. With a field of view of ∼2 sr and an angular resolution of ∼0.8◦ at 1 GeV, the LAT provides the most detailed observations of the gamma ray sky to date. The sec-tion will discuss the LAT’s method of detecsec-tion and known response funcsec-tions to incoming signals.

Section 4 will describe the detection method. The process starts with the analysis of known AGN to establish a general AGN Fourier power spectrum. In comparing this spectrum with that of the background, the most efficient test statistic to distinguish AGN signal from DGRB can be determined. This test statistic will be labelled the ‘variability score’ (similiar to the ‘periodicity score’ suggested by GSK). The calculated Fourier spectra will subsequently be used in a simulation of the gamma ray sky, which contains a given population of unresolved AGN in the background. The purpose of the simulation is to es-tablish the effectiveness of using the variability score as a detection method. A possible detection is simulated by calculating the variability score for a region of the simulated sky containing an AGN and subsequently checking whether this variability score is sufficiently distinct from that of the background. The detection level is the total number of detections in a simulated sky. The results of the simulation will conclude section 4.

From these results it is clear that the detection method is not sensitive enough to register any new AGN. The proper establishment of the background level from which the simulated AGN signals must be differentiated is the most problematic aspect of this method. In the presented simulation the number of background regions with a high variability score is too large. Possible improve-ments to the method, particularly as far as the analysis of the background is concerned, are discussed in section 5.

1.1 The Gamma Ray Background

Observations of the gamma ray sky have evolved over several generations of ex-periments: from the OSO 3 mission (1968), which first observed cosmic gamma rays, to SAS-2 (1972), COS-B (1975-’82), the Compton Gamma-Ray Observa-tory (1991-2000) and the current Fermi Space Telescope (2008-). It was the SAS-2 that first detected the presence of a diffuse background, a flux of photons unassociated with any resolved source. The increased angular resolution and energy range of the Energetic Gamma-Ray Experiment Telescope (EGRET) aboard the Compton observatory, allowed for a much improved analysis of the background, both of its spatial distribution and of its energy spectrum. The

(6)

1 INTRODUCTION 1.1 The Gamma Ray Background

data collected by the LAT continues this track of increasing refinement, provid-ing the most detailed information about the gamma ray background to date (as exemplified in figure 1).

The first step in establishing the background is to remove the point sources from the gamma ray sky. The most recent list of resolved LAT sources is the Fermi Second Source Catalog (Abdo et al. 2012). The resulting radiation, the DGRB, can be decomposed into two distinct contributions: the galactic diffuse emission, which varies at large angular scales, and the Isotropic Gamma Ray Background. The latter component is comprised of the Extragalactic Gamma Ray Background (EGRB) and possibly an isotropic flux of galactic origin. The IGRB is effectively the ‘remainder’ after known components have been removed from the gamma ray sky and there are currently large uncertainties about both its origin and intensity. Making a proper measurement of the IGRB or the EGRB is therefore a complicated process.

The uncertainties in the calculation of the IGRB stem largely from the sepa-ration of the diffuse galactic component from the rest of the DGRB. The galac-tic component, which is subtracted from the DGRB, is established theoregalac-tically, through the use of analytical models. The galactic diffuse emission is believed to originate in the interaction of cosmic rays with particles in the Interstellar Medium (ISM) and with the Interstellar Radiation Field (ISRF). In fact, the study of the DGRB provides an opportunity to study cosmic rays in regions other than Earth’s immediate vicinity. The isotropic distribution of cosmic rays results in the diffuse nature of the radiation, whereas the galactic distribution of gas and dust in the ISM causes the fluctuation at large angular scales. The varied nature of cosmic rays (the term covers both atomic nuclei and electrons) and their possible interactions within the galaxy, mean that the diffuse galactic component itself also has a number of different components.

The model used by the Fermi Collaboration to determine the diffuse galactic emission is GALPROP (Strong, Moskalenko, and Ptuskin 2007), which simulates both cosmic ray propagation in the galaxy and the gamma ray flux resulting from interactions. The input of the model is a given cosmic ray flux and cos-mic ray spectrum (derived from coscos-mic rays at Earth), gas column densities of neutral and molecular hydrogen in the Milky Way, based on multiwavelength surveys, the structure of the galactic magnetic field and the photon density of the ISRF. The cosmic ray interactions that result in gamma ray production are collisions of protons and He nuclei with gas in the ISM resulting in pion (π0₎

de-cay, bremsstrahlung of electrons interacting with the gas in the ISM and inverse Compton scattering of ISRF photons on cosmic ray electrons. GALPROP also includes the production and interactions of secondary cosmic rays, which con-tribute to the diffuse emission as well. Synchrotron emission due to the galactic magnetic field occurs at energies below that of gamma rays.

The general procedure to find the IGRB energy spectrum, or Spectral En-ergy Density (SED), from an all sky map of the DGRB begins with masking the galactic plane, by restricting data to regions with |b| > 10◦, or even higher galactic latitudes. The contribution of the diffuse galactic emission at higher latitudes is provided by GALPROP and subtracted from the DGRB map. These

(7)

Figure 1: The EGRB energy spectrum (SED) as detected by EGRET (blue and red) and the LAT (black). The LAT spectrum has a comparable index of Γ ≈ 2.4, but a lower integrated flux, approximately 1.03 × 10−5ph cm−2s−1sr−1, for photon energies higher than 100 MeV. Figure and data from Abdo et al. 2010d.

steps are repeated for different energy ranges. The result is a spectrum for the remaining IGRB, shown in figures 1 and 2. Figure 1 displays the spectra derived both from EGRET and from LAT data; the superiority of LAT data is quite clear. Figure 2 shows the composition of the total gamma ray flux as determined by the Fermi Collaboration. The complexity of modelling the diffuse emission still leaves relatively large uncertainties. The Fermi Collaboration notes that particularly IC processes can occur up to high latitudes, depending on the size of the galactic halo, without a significant gradient in the intensity of the gamma ray flux. This would contaminate the determined IGRB (Strong, Moskalenko, and Reimer 2004; Abdo et al. 2010d).

The nature of the IGRB has been determined with even less certainty than that of the galactic gamma ray background. The energy spectrum is relatively smooth and featureless. Both Strong et al., working with EGRET data, and the Fermi Collaboration note that fitting a simple power law to the spectrum results in a spectral index Γ ≈ 2.4, very similar to the spectral index associated with blazars (discussed in section 1.2). As mentioned earlier, the IGRB is generally believed to have a large extragalactic component, possibly constituting 100% of the remaining background. In fact the terms EGRB and IGRB are sometimes used interchangeably; as this thesis focusses on the detection of extragalactic components to the background, the term EGRB will be used from here on out. Like the components of the gamma ray sky, the EGRB is itself the sum of many (possible) contributions.

Of the possible components to the EGRB, unresolved blazars are the most likely candidates. The similarity of their average spectrum to that of the EGRB,

(8)

Figure 2: The composition of the total gamma ray flux as determined by the Fermi Collaboration. The galactic diffuse component is calculated using GALPROP, the sources are listed in the Fermi 2nd _{Source Catalog and the CR background is an}

irreducible local background for the LAT. The model matching the DGRB data points is a combination of a modelled IGRB and GALPROP. Figure from Abdo et al. 2010d.

as well as their probably relatively isotropic distribution across the sky mean that blazars can make up a significant section of the diffuse extragalactic flux. Estimates for this contribution range from 15% to 100% of the EGRB. Other known source classes that could have significant contributions are non-blazar AGN, starburst galaxies, normal star-forming galaxies, extended emission from galaxy clusters and millisecond pulsars. Of these classes the non-blazar AGN and starburst galaxies could have the same relevance as the blazars.

Note that the contributions not only vary in relative strength, but also over the energy range, as different source classes have different energy spectra. The contribution of normal star-forming galaxies has been estimated to be as high as one third of the total EGRB for energies around 1 GeV (Pavlidou and Fields 2002). The uncertainties on the limits of the individual components are so large, however, that the presence or absence of a diffuse component, possibly due to dark matter annihilation, cannot be confirmed. The investigation into the EGRB can therefore provide valuable insights for the study of dark matter as well. The main issue is to find reliable contribution levels from the different source classes.

In order to properly establish the different contributors and their relative strength, the EGRB has been investigated using several different approaches. The first method is to use the EGRB energy spectrum, as described above. The resemblance of the background spectrum to that of blazars is a clear indication

(9)

of their presence in the EGRB, but the lack of features makes it difficult to pinpoint other specific components. More information is therefore required. A second approach is based on obtaining a solid hold on the populations of the different source classes, for fluxes below the LAT’s detection limit.

The observed population of a given source class can be extrapolated to lower fluxes. This method has become much more reliable since the LAT was used to detect and identify many new blazars and observe previously detected objects in greater detail. In a comprehensive approach, the source count distribution as well as the distribution of spectral indices among the objects can be established for all known objects of a certain class. The Fermi Collaboration performed this study for high galactic latitudes and found that for different classes of blazar the source count distribution, meaning the number of objects at a certain flux, was best fitted by a broken power law, with a steeper power law for higher fluxes. In this approach, the distribution of spectral indices of the objects was established as well, and the total flux of unresolved AGN-like point sources was estimated by integrating a count distribution dependent on both flux and spectral index over the entire range of both variables. The contribution of blazars as estimated by the Fermi Collaboration is show in figure 3 (Abdo et al. 2010c).

Another method to approximate the population of a known source class be-low the flux detection limit is to use photon counting statistics, a technique that has been applied in radio and X-ray astronomy as well (Malyshev and Hogg 2011). The main idea of this approach is that at least some component of the EGRB is due to unresolved point sources, generally associated with AGN, and is therefore not a truly diffuse component. The contribution of many indi-vidual sources will affect the photon distribution in the background, even when the sources cannot be detected directly. The distribution of photons across the background is modelled as a function of the distribution of AGN and an average AGN energy spectrum. This model is fitted to the LAT measurements of the DGRB (again for high latitudes). The model parameters providing the best fit give the most likely source count distribution for AGN like point sources. Both the Fermi Collaboration and Malyshev & Hogg find a maximum blazar contribution to the EGRB of ∼ 25%-40%, of course depending strongly on the modelled galactic component. The results of the two methods to determine the source populations are consistent with the notion that the EGRB has a number of different origins.

Other than the energy spectrum and the estimation of possible source populations, a third approach in analysis of the EGRB is provided by spatial information about the gamma ray sky. Whilst the EGRB is isotropic at large scales, small angular scale anisotropies can be indicative of an underlying struc-ture to the gamma ray emission. The high angular resolution of the newer generation of gamma ray telescopes makes this method particularly viable. The angular power spectrum is calculated by expanding the spatial distribution of the detected fluxes in spherical harmonics, where the power is given by the

(10)

Figure 3: The contributions of the two classes of blazar Flat Spectrum Radio Quasars (FSRQs) and BL-Lac objects (section 1.2) to the EGRB (labelled here simply as “Diffuse background”). The total flux is calculated by integrating over the source count distribution, as indicated in the text. The total blazar contribution appears insufficient to explain the whole EGRB. Figure from Abdo et al. 2010c.

square of the coefficient alm in

I(θ, φ) =X

lm

almYlm (1)

Here I represents the intensity of the radiation as a function of direction on the sky and Ym

l are the spherical harmonics. This procedure can be considered as

the spatial (spherical) equivalent of the Fourier transform in the time domain, which is discussed in section 2. The coefficients weigh the strength of intensity variations at different angular separations and therefore provide a measure of the ‘bias’ or clustering of unresolved point sources or of the diffuse dark matter. The angular spectrum can be used to infer the distributions of known source classes. Different source classes have different spatial distributions, which can be derived from observations of resolved objects of a given source type. A benefit of this method is that the spatial distribution can be derived from observations in other wavebands than gamma rays, creating a larger statistical base for the comparison of different contributors to the background. This method can also be used to detect a possible dark matter annihilation signature, as the angular spectrum associated with the distribution of dark matter differs significantly from that of unresolved point sources. In a two component simulation (blazars and DM) the distribution of dark matter is modelled in halos and subhalos, following well known, computationally established density distributions, such as the Navarro-Frenk-Wright (NFW) profile. In specific energy intervals where the DM annihilation component is approximately one third of the EGRB, which

(11)

1 INTRODUCTION 1.2 Active Galactic Nuclei

is definitely possible in the likely scenario that blazars make up a maximum of ∼40%, the dark matter component should be detectable through the use of the angular spectrum. This signature has not yet been detected (Ando et al. 2007). In a study based on a little under two years of LAT data, the Fermi Col-laboration calculated the angular power spectrum for the isotropic background and fitted it to two different model backgrounds. The total angular power for high multipoles (155 ≤ l ≤ 504) shows a clear excess over the Poisson spectrum, which is the angular spectrum associated with the shot noise of an otherwise perfectly isotropic distribution. The excess was not as high as had been pre-viously modelled for different classes of unresolved point sources. The relative independence from l of the angular power in the given multipole range indi-cates that the portion of the EGRB associated with this excess is likely due to unclustered point sources below the flux detection limit. Using the spatial distributions established for known sources (in multiple studies and in different wavebands) the maximal contributions of different classes of point sources to the EGRB were established. Based on the angular spectrum, the contribution of blazars was limited to ∼20%. Other classes such as normal starforming galaxies and dark matter distribution could still contribute up to 100% of the EGRB, whereas the maximal millisecond component is almost negligible (Ackerman et al. 2012a).

The three main approaches to EGRB analysis outlined above, study of its energy spectrum, extrapolation from established source class populations and the angular power spectrum, all make use of the LAT’s high angular and en-ergy resolution. The difference in approach suggested by GSK was to apply the timing data registered for all LAT events as well, thereby opening up a new di-mension for study. The method presented in this thesis will use the timing data to target a specific subset of the source class candidates for the EGRB: AGN, both blazars and ‘regular’ AGN. By calculating the Fourier power spectrum the inherent variability of these sources can be quantified. The Introduction will continue with a brief discussion of AGN and their known emission properties.

1.2 Active Galactic Nuclei

1.2.1 AGN Unification and Spectral Classification

Active Galactic Nuclei are among the most powerful sources of radiation in the Universe, with luminosities from 109 _{L to 10}14 _{L. They are characterised} by emission over a broad spectral range, from radio up to gamma rays and TeV energy photons, by emission from a relatively small region and by relativistic internal motions. It is believed that AGN are powered by a supermassive black hole (106_-109

M) located at the centre of a host galaxy. The population of AGN has many subtypes and variations, which can have strikingly different appearances. One class is formed by radio galaxies, which are marked by two symmetric lobes extending from the centre of the galaxy, approximately orthog-onal to the galactic plane. A narrower substructure, referred to as the jet, can be seen closer to the galactic centre for most of these objects. The lobes can

(12)

stretch up to several kpc and emit non-thermal radiation, peaking in the radio band. Seyfert galaxies are optically bright spiral galaxies that show emission lines in their spectra, in contrast to the usual galactic spectrum comprised of stellar spectra containing absorption lines.

An apparently different class altogether is that of the quasars, quasi stellar objects, which were first detected as point sources in the 1960’s. These objects have such complex and distinctly non-thermal spectra as well as such high red-shifts that they could not possibly be identified as stars. The realisation that all these seemingly disparate objects were different representations of a single or small group of source types was a gradual process. Although a great effort has been made to better understand these sources, the precise AGN unification scheme is still under debate (Sparke and Gallagher 2007).

The main aspect of AGN unification, in particular the unification of quasars and radio-galaxies, is the notion that the differences between the objects are caused by a difference in viewing angle. This idea was developed and refined in the 1980’s and ‘90’s (Urry and Padovani 1995). As mentioned before, AGN are considered to be located at the centre of a host galaxy. In the case of quasars the emission is so powerful (and the redshift so large) that the AGN outshines its host. The basic, and clearly simplified, model of an AGN shows a highly anisotropic emission structure, resulting in large differences in ob-served radiation for different viewing angles. This model is shown in figure 4.

Figure 4: Basic AGN model: a black hole surrounded by an accretion disk and torus, with relativistic jets emitted at right angles to the accretion disk plane. The line emit-ting regions are clouds of hot gas, where the observed difference in line width is due to a difference in rotational velocity. Fig-ure from Urry and Padovani 1995.

The central structure of an AGN is formed by a central black hole, with an accretion disk of galactic gas surround-ing it. Around this inner region is a larger band, or torus, of obscuring gas and dust. Two relativistic jets point away from the black hole, in opposite directions and at approximately a right angle to the accretion disk. The broad and narrow line regions are clouds of hot gas that produce thermal radiation. The gas closer to the black hole moves at significantly higher speeds, ∼107 m s−1_{, as it is deeper in the gravitational}

potential, resulting in broader emission lines (in fact: the width of the lines is used to determine the gas rotational ve-locity). All different regions contribute to the overall radiation, but the jet is by far the most luminous component, when viewed in the direction along the jet axis. The jet luminosity is based on highly energetic charged particles, which are ejected at relativistic speeds by the black hole in a process that is

(13)

not well understood. These particles, often modelled as electrons, lose their energy through synchrotron radiation, due to the presence of a strong magnetic field, and through Compton scattering of available photons.

The intensity of the radiation observed by someone in the direction of the jet axis is caused by a combination of several relativistic effects. The first is relativistic beaming, or the aberration of light: for a stationary observer the emission of a moving source, radiating isotropically in its own rest frame, is compressed into the direction of motion. In the case of a high Lorentz fac-tor (Γbulk) of the bulk motion of the jet (which is the case for powerful AGN)

this implies that approximately half of all the radiation is emitted in a cone centered on the jet axis with opening angle Γ−1_bulk. The second factor is a rela-tivistic Doppler correction, due to the constant motion of the radiating matter along the jet. The final factor is simply the effect of the Lorentz transformation of the jet frame to that of a stationary observer. Combined these effects provide a ‘Doppler boost’ factor to the observed radiation energy of (γ(1 − βcos θ))4_.

Here γ and β are both defined by the bulk motion and θ is the angle between the jet axis (the direction of motion of the radiating paricles) and the line of sight (Ghisellini et al. 1998).

The result of this distribution of the emission is that the intensity of the source and the energies that are observed are highly dependent on the angle of observation. Observers looking straight into the jet will detect intense, high energy radiation. Due to the beaming effect a small increase in θ will cause the observer to note a sharp decrease in luminosity; on average a change from 0◦ to 10◦ will result in a flux reduction by a factor ∼250. For θ ≈ π/2 (i.e. galaxies that are viewed edge-on) the jet luminosity is minimal and the broad line region and disk are obscured by the torus, which means these AGN are classified as narrow line objects. AGN with one of their jets directed at Earth are known as blazars. In general blazars are characterised by high luminosity, flat radio spectra associated with non-thermal radiation, high optical polarisa-tion and variability at several timescales (Abdo et al. 2011). Blazars’ luminosity and their hard spectra make them prime candidates for Fermi-LAT as well as important contributors to the EGRB.

The unification scheme suggested by Urry and Padovani (1995) applies to ra-dio loud galaxies and blazars. All known blazars are rara-dio loud, which is defined by a relation of radio flux to optical flux of 10:1, and no jets have been detected in radio quiet galaxies. The scheme is based on a subdivision of both blazars and radio galaxies into two classes. The division of blazars is between BL Lac objects (named after the prototypical BL Lacertae) and Flat Spectrum Radio Quasars (FSRQs), where the basis of separation is that FSRQs show emission lines in their spectra, whereas BL Lacs do not. The subtypes of radio galaxies are the narrow line classes (Fanaroff-Riley classification) FRI and FRII. FRI galaxies have radio lobes that decrease in luminosity with distance from the nu-cleus, whereas FRII galaxies have lobes with bright, clearly defined outer edges. FRII galaxies also have ‘blobs’ of increased radiation at pc scales, visible inside the lobes. The unification works as follows: BL Lacs are FRI galaxies with the jet pointed towards Earth and similarly for FSRQs and FRII galaxies. This

(14)

sorting is based on matching power in the radio spectrum, similarity in host galaxy and in the case of FSRQs on similarity in narrow line emission (Urry and Padovani 1995). The distinction based on host galaxy is not very strong, as all blazars have been detected in giant ellipticals.

The effect of the anisotropic emission in AGN is that observations in the gamma ray part of the spectrum are mainly constrained to blazars, as the high-est energy output is formed almost entirely in the jet and the emission is strongly beamed. The high resolution observations of the LAT have therefore resulted in a large increase in the number of detected and confirmed blazars in recent years, providing a wealth of data about the populations and specific characteristics of the different blazars types. The definition used by the Fermi Collaboration to differentiate between BL Lacs and FSRQs is that the equivalent line width of the strongest emission line in the spectrum must be lower than 5 ˚A for a BL Lac blazar. The confirmation of a new blazar is made by comparing the coordi-nates of a promising candidate, based on its spectrum, with the coordicoordi-nates of galaxies detected at other wavelengths. In this manner the LAT AGN catalog contains ∼900 AGN, of which ∼45% are BL Lacs, ∼35% are FSRQs and ∼18% are of currently undetermined type. It should be noted that the relative dearth of observations at other wavelengths for the Southern Hemisphere are in part the cause of the large fraction of the unconfirmed blazar class, as the distinction between FSRQs and BL Lacs requires optical spectra to analyse the emission lines (Abdo et al. 2011).

The surveys of AGN by the Fermi Collaboration provide some clear prop-erties of the blazar populations. In general FSRQs have a higher luminosity than BL Lacs: quasars are intrinsically more powerful sources. The stronger gamma ray luminosity matches the difference in radio flux detected for the par-ent populations of FRI and FRII galaxies. The two categories of blazars are also different in their energy spectra. On average FSRQs have softer spectra, with a spectral index of 2.42, and their spectra are fitted better with a broken power law. The BL Lacs in the sample have a flatter average spectrum, with an index of ∼2.1, but also have a much larger distribution of the indices among the individual objects. This distribution is indicative of a subdivision of the BL Lac class, that will be discussed further on. Here, it is important to note that the spectral indices provide a clear distinction between FSRQs and Bl Lacs, as can be seen in figure 6, in addition to the differences in emission lines and luminosity (Abdo et al. 2011).

Only a very small percentage (<1%) of the AGN observed by the LAT are so-called misaligned AGN (MAGN), in which the angle between the jet and the line of sight is larger than for a blazar. These objects have significantly lower intrinsic fluxes and show steeper energy spectra than blazars. The reduced out-put in the highest energy domain is to be expected due to the strong decrease in beaming. This suggests an interpretation of the Steep Spectrum Radio Quasars (SSRQs) as broad line radio galaxies (in which the inner region around the black hole is not obscured by the torus) that have such such a high luminosity that the AGN outshines its host galaxy (Abdo et al. 2010a). The relatively low flux makes MAGN significantly less likely to be resolved by the LAT. However, their

(15)

contribution to the EGRB could still be significant, as discussed in section 1.1, particularly as the number of these objects in the sky will be higher than that of the blazars.

1.2.2 Modelling the Blazar Spectrum

The classification system of AGN, strengthened by the observational data from the LAT, provides the framework for an effort to understand the composition of AGN signals. One of the key topics in the analysis of the AGN emission process is an attempt to recreate the detected energy spectrum, over the en-tire frequency range. Given the intense luminosity of the jet, AGN spectra, in particular those of blazars, are dominated by the processes that take place in this component. However, other emission regions such as the disk, torus and line regions also contribute to the overall output, both through interaction with the jet and directly. The process of identifying the different components of the energy spectrum is focussed on blazar spectra, as these provide the clearest view of the jet. Currently, one of the most favored emission models is the one first put forth by Ghisellini et al., which describes the ‘double hump’ nature found in the blazar spectra and organises blazars in the ‘blazar sequence,’ which could provide insights into AGN evolution as well (Ghisellini et al. 1998).

The two main emission processes in this general model are synchrotron emis-sion and inverse Compton scattering. These two processes cause the two ‘humps’ in the spectrum. The low energy peak is related to synchrotron (peak frequency: νS

peak) and the high energy peak belong to the IC process (νpeakIC ). In both cases it

is often assumed that the population of radiating sources is limited to electrons, although hadronic models can also reproduce the desired spectrum (Ghisellini 2013; B¨ottcher et al. 2013). Electrons will be considered the source population for the following description and a short consideration of hadronic models will be provided further on in the text. The mechanism powering the jet acceler-ates the electrons up to relativistic speeds, after which the particles lose this energy radiatively or thourgh interactions with existing photons. The presence of a magnetic field in the jet causes the synchrotron radiation, which peaks at radio wavelengths, in the jet rest frame. The inverse Compton scattering can be dominated by two populations of photons: those resulting from the synchrotron emission, referred to as synchrotron self Compton (SSC) or external photons (EC) provided by the line emitting regions, disk or the CMB. In general the EC process is more efficient, allowing the electrons to cool quicker (Ghisellini 2013). An example of a typical spectrum (note: in the jet rest frame) is provided by the FSRQ 4C+0102, shown in figure 5, in which all individual fitted spectral components are identified.

The categorisation of the blazar sequence is based on the observation that on average ν_peakS and ν_peakIC occur at lower energies for FSRQs than for BL Lacs. The explanation for this phenomenon is that the more powerful FSRQs have a stronger magnetic field, allowing for faster cooling through synchrotron emis-sion. The electrons therefore have a lower average energy, resulting in a ‘redder’

(16)

Figure 5: Spectral energy density for the FSRQ 4C+0102, with the frequencies as observed in the jet rest frame. The grey datapoints are from the NED (NASA Ex-tragalactic Database), the other datapoints are from UVOT, XRT, BAT and LAT (in ascending max. energy). The grey area top right indicates the sensitivity of the LAT in the first three months of observations. The solid blue line is the sum of all spectral components, fitted to the data. The double hump structure of the synchrotron and IC components (summed together in the solid light green line) is clearly visible. Be-cause the source has a ‘red’ spectrum, the disk component is distinct from the other contributions. The label “corona” refers to the line emitting regions. Figure from Ghiselllini et al. 2011.

spectrum, in which both ν_peakS and ν_peakIC are shifted to lower energies. A second effect is that the relative strength of the IC component (referred to as “Comp-ton dominance”) is less for FSRQs. The efficiency of the cooling is increased by the presence of large numbers of external photons, from the strong line emitting regions present in FSRQs, which intensify electron energy loss via EC. The op-posite effect takes place in the less powerful BL Lacs: the electrons do not cool efficiently, causing them to have higher average energies. This shifts both νS

peak

and ν_peakIC to more ‘blue’ regions of the spectrum and increases the Compton dominance (Ghisellini 2013).

The process of classifying blazars by a peak emission frequency (ν_peakS ) re-quires a subdivision of the BL Lac class. Whereas all FSRQs will have ‘red’ spectra due to their high luminosity, the spread in power among BL Lacs is reflected in a wider distribution of νS

peak. The terminology associated with this

subdivision is Low, Intermediary and High Synchrotron Peak (LSP, ISP and HSP) blazars; all FSRQ are LSP, the BL Lacs can be different classes. Among BL Lacs, the highest power objects will be LSP BL Lacs running down to the

(17)

lowest power HSP BL Lacs. For HSP BL Lacs, ν_peakS will lie in the X-ray spec-trum (for a stationary observer aligned with the jet) and the LAT detected gamma ray spectrum will cover the falling slope of the IC component. LSP blazars will typically have ν_peakS in sub-mm range, and a steeper falling IC slope in the gamma spectrum. Interestingly, the spectral distribution for LSP blazars often leaves the thermal component due to the disk, generally peaking in soft X-rays, more clearly visible. A direct observation could allow for spectral fit that can constrain the mass of the central black hole (Ghisellini et al. 1998).

The classification set out in the previous paragraphs is the essence of the blazar sequence: the similarity between FSRQs and LSP BL Lacs suggests a close connection between the two. The structure of the classification is clear: from high power FSRQs with red spectra, to ever less luminous BL Lacs, ending in the blue spectrum HSP BL Lacs. The gradual transition between the different subclasses is suggestive of an evolutionary track. A possible mechanism for this evolution is a change in the accretion flow. This idea is supported by the fact that the divide between FSRQs and BL Lacs appears to follow a critical ratio between the Eddington luminosity (an accretion limit) and the luminosity of the blazar (Ghisellini 2013). The powerful FSRQ jets are thought to be powered by large, active accretion disks, which deplete over time. As the jet loses power, the cooling process becomes less efficient and the spectrum becomes more alike to LSP BL Lacs, moving down the line to HSP BL Lacs over time (Xie et al. 2006; Ghisellini et al. 1998).

The central ideas of the blazar sequence were developed in the 1990’s for the limited sample of blazars then observed, making a comparison of this model with the new LAT data particularly important. It turns out that the LAT sam-ple of AGN corresponds well with the ‘double hump’ spectral model and the associated blazar sequence. The Fermi Collaboration reports a clear correlation between νS

peak and the spectral index, as seen in figure 6. A high νpeakS implies a

low spectral index (i.e. a harder spectrum) and objects with a low νS

peak have a

higher index, corresponding to the anticipated spectra for an HSP BL Lac and FSRQ respectively. The average softer gamma ray spectra for FSRQs appear to indicate a lower νIC

peak for objects with strong emission lines, highlighting the

effectiveness of EC cooling by the more abundant external photons, over SSC. The theoretical spectral composition therefore coincides well with the observed blazar spectra (Abdo et al. 2011). The observed correlation between jet power and disk luminosity indicates that the accretion rate is a key component in the power of the jet, and subsequently the luminosity of the blazar (Ghiselllini et al. 2011).

A separate aspect of the LAT observations, particularly relevant to the evo-lutionary interpretation of the blazar sequence, is the distribution of the blazar classes, and their relevant spectral parameters, as a function of redshift. The relation between gamma ray luminosity and known redshift is given in figure 7. The split between FSRQ and BL Lac is clearly visible in this distribution as well: FSRQs are overwhelmingly found at high redshifts (up to z=3.2), whereas BL Lacs are exclusively found at lower redshifts. A bias in the detection is that

(18)

BL Lacs generally have lower luminosities, making it less likely to detect them at high redshifts. The AGN catalog lists approximately 55% of all BL Lacs with unknown redshift (Abdo et al. 2011). The main problem in finding this distri-bution is the lack of lines in the typical BL Lac spectrum. For these objects the redshift can be estimated from other parameters such as a reduced flux in higher wavebands due to HI absorption in the inter cluster medium, absorption lines by an intermediary object and fitting of the host galaxies. Using these methods it is possible to track the populations of all blazar types (Ajello et al. 2014).

The cosmic evolution of a population refers to the change in number density with increasing redshift. The FSRQs show a positive evolution, as could be expected from figure 7, with the number density peaking around z=1. On av-erage the BL Lacs appear to have a negative evolution, but taking the different subclasses into account, it turns out that only the HSP BL Lacs have a negative evolution, whereas ISP and LSP BL Lacs do increase in number with higher redshift. The peak of the BL Lac number density is around z≈0.6 and there are no BL Lacs with z>1.5. The general evolution of blazars is therefore one in which the luminosity increases and the average spectral index decreases with increasing redshift (Ajello et al. 2014).

Of particular interest to the proposed system of blazar evolution is the appar-ent inverse correlation between FSRQ abundance and HSP BL Lac abundance. At high redshift the powerful quasars are most prominent, possibly fueled by galaxy mergers. The expansion of the Universe decreases the number of mergers and the FSRQs are ‘burning’ through the available gas by means of the accretion disk. After z=1 the number of FSRQs starts to decline, whereas the number of BL Lacs is increasing. First LSP and ISP BL Lacs are more prominent, but the population of BL Lacs becomes increasingly dominated by HSP BL Lacs. In the LAT sample, ∼53% of the BL Lacs are of the HSP type. This scenario is certainly appealing due to its simple causal nature: high power, high accretion FSRQs are formed, run out of fuel and the lower power, low accretion HSP BL Lacs remain. As Ajello et al. 2014 stress, however, there is much reason for caution: the level of uncertainty in the BL Lac redshifts is still high and the complex nature of the objects in question means that there are many outliers in the distributions that do not fit the theoretical model well.

1.2.3 Alternative Emission Models and Variability

Following this phenomenological discussion of AGN, and blazars in particu-lar, it is informative to briefly consider alternative ways in which the radiative processes in the jet can be modelled, as well as the mechanisms powering the emission process. As mentioned earlier, the general approach to AGN emission is to assume electrons as the radiating particles. The model by Ghisellini et al. that successfully recreates the blazar energy spectrum is based on this as-sumption. It is, however, also possible to use both protons and electrons in the model, although these hadronic models are significantly more complex.

(19)

Figure 6: Relation between the spectral index and νpeakS for LAT blazars. Red indicates

FSRQs, green LSP BL Lacs, light blue ISP BL Lacs and dark blue HSP BL Lacs. Both the distinction between FSRQ and BL Lac as well as between the BL Lac subtypes is clearly visible. The high energy electrons, due to ineffective cooling, in the HSP BL Lac provide the ‘bluest’ νS

peakand flattest spectrum. Figure from Abdo et al. 2011.

Protons can radiate via proton synchrotron, if the magnetic field in the jet is strong enough to contain the Larmor radius inside the jet radius. Additional processes, applicable to protons only, are pion (π0_{) production and the}

produc-tion of ultra high energy photons (Eγ TeV), which both produce cascades

of decay particles. To combine these processes and their interaction with the leptonic component (particularly through secondary IC) time consuming Monte Carlo simulations are often necessary (B¨ottcher et al. 2013). Hadronic models can provide good fits to most blazar spectra, providing the same quality of fits as lepton-only models. They have the added benefit that they account for emis-sion and variability in the TeV range, which has been detected for ∼45 LAT blazars (Abdo et al. 2011).

Another aspect of modelling which has not been discussed up to this point is that of timing. The double hump model provides a good fit for the spectrum integrated over time, but variability is one of the key characteristics of blazar output. The recreation of the blazar lightcurve is therefore another important goal of AGN modelling. One of the time-dependent processes possible inside the jet are internal shock waves, that occur due to variations in the matter in-jection into the jet. The effect occurs when a higher than average density ‘blob’ of matter is injected into the jet, followed by a second blob moving at a higher velocity.

When the two blobs collide, a shock wave is created, which causes first order Fermi acceleration for the charged particles (here considered as electrons only) each time they run into the wavefront. As long as the collision between the blobs

(20)

Figure 7: Relation between gamma ray luminosity and redshift for LAT blazars. Red indicates FSRQs, green LSP BL Lacs, light blue ISP BL Lacs and dark blue HSP BL Lacs. The solid and dotted lines indicate estimated LAT detection limits based on blazar spectra with index 1.8 and 2.2 respectively. FSRQs and LSP BL Lacs are located at higher redshifts than most BL Lacs. There is a clear correlation between Lγand z, although a detection bias favouring brighter objects needs to

be taken into account. Figure from Abdo et al. 2011.

lasts, the radiation from the blazar is seen to increase in intensity. When the shock wave dissipates, the usual radiative processes cool the energised electrons. Because cooling processes with different inherent timescales contribute to dif-ferent regions of the energy spectrum, a time lag in variability between difdif-ferent wavebands is created. These internal shock processes are thought to contribute to the variability of blazars on timescales from seconds to hours (B¨ottcher and Dermer 2010).

Despite the large effort that has been made to improve the understanding of AGN emission, the underlying mechanism driving the entire process is still an issue of debate. As mentioned in the beginning of this section, the main agent in the process is thought to be the supermassive black hole at the centre of the system. The potential energy of the material in the accretion disk provides the main input for the mechanism. The rate at which this energy is released gravitationally is determined by the accretion efficiency of the disk. Although a correlation between accretion rate and blazar luminosity does appear to exist, the strength of this correlation and subsequently the precise way in which the disk and jets interact is still uncertain (Ghiselllini et al. 2011; Ajello et al. 2014). Another way in which energy can be extracted from the black hole and its surroundings is the Blandford-Znajek mechanism (Blandford and Znajek 1977). In this process the rotational energy of the black hole can be altered via

(21)

elec-1 INTRODUCTION 1.2 Active Galactic Nuclei

tromagnetic interactions with the magnetic field that is created by the rotating plasma of the accretion disk. In essence the black hole acts as an electromagnetic generator. This is possible because black holes effectively interact as resistors. The exact process depends on the geometry of the Kerr-metric and is too com-plex to discuss in detail here. However, the general principle can be understood as follows. Accreting plasma moves towards the black hole, dragging the mag-netic field lines due to the acrrection disk’s rotation along with it. The field lines become bunched and tangled together close to the black hole, due to its rotation. A voltage is generated because charges rotating alng with the black hole experience a Lorentz force due the present magnetic field. The associated electromotive force is proportional to the angular speed and size of the black hole and to the strength of the magnetic field. The voltage results in a current when begins to flow along the magnetic field lines. In this manner the rotational energy of the black hole is converted into kinetic energy of the plasma, which is ejected at the rotational poles (Hartle 2003). The Blandford-Znajek process and the direct accretion of matter are considered to be important processes in powering AGN.

The final aspect of AGN to be discussed is their variability, which will be of key importance for the detection method presented in this thesis. Before the availability of a large sample of observed gamma ray blazars, the study of AGN variability was mainly contained to X-ray data. The timescales of X-ray variability can be as low as minutes, which has led to the conclusion that vari-ability could be directly related to the black hole mass. The reasoning is that the compactness of the emitting region around the black hole implies that any fluctuations happen within a small space, which means that there is only a small crossing time for any photon or relativistic electron, this in turn means that the maximum fluctuation timescale is low. One important avenue of analysis in this field is to assume a linear scaling relation between black hole mass and X-ray variability. By comparing stellar mass black holes in X-ray binaries, such as Cygnus X-1, with AGN, the mass of the black hole in the AGN can be con-strained (Hayashida et al. 1998; Vaughan, Fabian, and Nandra 2002).

The approach is based on calculating the normalised power spectrum for the discrete Fourier transform (section 2) of the timing data, generally referred to as the Power Spectral Density (PSD). This spectrum is then fitted with a power law. Any break in the power law would be indicative of a characteristic timescale of the object, which is simply defined as a limit on variability im-posed by some aspect of the AGN. A break can be found in both the power laws fitted to stellar mass black holes and to AGN. In combination with the estimated mass of the stellar mass black hole, derived from the movement of the companion star in the binary system, the ‘variability timescales’ of the two object can be used to derive the AGN black hole mass. The main drawback to this approach is that X-ray observations are generally not continuous, which means that the calculated variabilities are based on a limited number of sepa-rate pointings (Hayashida et al. 1998; Vaughan, Fabian, and Nandra 2002).

In a survey of 104 AGN, based on XMM-Newton data, Gonz´ales-Mart´ın and Vaughan 2012 found that approximately 72% of the AGN show some degree of

(22)

variability, in at least one of the three applied X-ray energy ranges. The sample contained seven particularly bright blazars and all AGN are located at z<0.4. The timing data was binned in 100s wide bins, and the calculated power spectra showed an average spectral index of ∼2. Most of the X-ray variability was found to lie in the low frequency ‘red noise’ range, which the authors associate with perturbations in the innermost regions of the AGN system (i.e. in the accretion disk, near the black hole). No truly periodic components to the signal were found for any of the AGN, which is to be expected as AGN with such char-acteristic are exceedingly rare. The power law fits to the AGN spectra again showed a break for most sources, which allowed for a derivation of the black hole mass in comparison with stellar mass black holes. The results showed an average AGN black hole mass between 106

M and 108

M (Gonz´ales-Mart´ın and Vaughan 2012).

A study very similar to that in the X-ray band was performed by Nakagawa and Mori 2013 based on LAT gamma ray data. The analysis is based on blazar lightcurves from the Fermi monitored source list, provided on the LAT website (Fermi Collaboration 2014a), where the flux has been binned in daily bins. The sample consisted of 15 of the brightest blazars, for which the power spectrum was calculated. Only for four of the blazars a broken power law provided a better fit than a power law with a single spectral index. The main issue in the approach was that the presence of noise was difficult to account for, which could well mean that the flattening in the higher frequency part of the spec-trum, causing the break in the power law, was the result of a reduction to the flat white noise spectrum (Nakagawa and Mori 2013).

A more comprehensive review of AGN variability in the gamma ray spec-trum was provided by the Fermi Collaboration (Abdo et al. 2010b), and is based on a preliminary sample of the 106 brightest AGN (104 blazars) after 11 months of LAT observations. A number of techniques were used to determine the variability of the AGN. A χ2 _{test comparing the flux in a time bin with the}

average flux of the light curve gives a simple indication of whether a source is variable or not; approximately 64% of blazars in the sample are variable. The Discrete Autocorrelation Function (DACF) and Structure Function (SF) each measure the variability of the source in the time domain, comparing two points in the lightcurve separated by ∆t. The level of variability detected is commen-surate with that derived from the χ2 _{test and no periodic signal components}

are detected. Although the DACF and SF provide useful information, the au-thors state that characteristic variation time scales are better explored using the Fourier Power Spectrum.

The power spectrum is calculated for a subset of the best observed sources (9 FSRQs and 6 BL Lacs) for a time series where the bin width has been set to three days; this procedure approaches the weekly time binning used for the detection method presented in this thesis (Section 4). The spectra were fitted with power laws. There is a notable difference between the power spectra of FSRQs and BL Lacs: the FSRQs tend to have flatter spectra, average index ∼1.4, while the BL Lacs show an average spectral index of ∼1.7. There is no indication of a break in the power law for any of the spectra. An interesting

(23)

anomaly presents itself in the subclasses of the BL Lacs. There is a large spread in the distribution of the spectral indices among the subclasses, but the HSP BL Lacs appear to have the flattest spectra. This is not in accordance with the expected similarities between FSRQs and LSP BL Lacs, based on the blazar sequence. The size of the sample limits the possibility for strong conclusions however (Abdo et al. 2010b).

The slopes of the power spectra match up well with the detected slopes in optical and radio wavebands. The correlation with the X-ray slopes is less clear (although e.g. the HSP BL Lac Mrk 421 has an index op ∼1 in both gamma and hard X-rays), most likely as a result of the more intermittent type of obser-vation common for that waveband. Overall, the variability of AGN appears to fall into one of two “flavours”. The first is the so-called random walk processes: variations and turbulence in the accretion flow, resulting in complex time pro-files and detectable variability in the low frequency range. This is the variability on which the presented detection method will largely rely. The second flavour is that of a relatively constant flux with large flares at random intervals. The flares are strong increases in luminosity that can last from minutes to days. The internal shock model discussed above provides a possible explanation for this effect, when large irregularities in the accretion flow cause the formation of high density ‘blobs’ in the jet (Abdo et al. 2010b).

Concluding this section, it should be stressed that the investigation of AGN is relevant for a broad scientific field. It provides insights into very high energy processes, insights into galaxy formation and insights into large scale structure when blazars are applied as probes in cosmological studies. However, for this thesis the relation between AGN and the gamma ray background is the most important. The high luminosity and high energies of blazars in particular make AGN excellent candidates for components of the EGRB. As was covered in sec-tion 1.1, the exact contribusec-tion is still very much unclear, but is expected to be significant. The detection method presented here will make use of a broad sample of AGN (64), encompassing different spectral and variability types. The calculated spectra will be used individually, providing representations of their AGN types for the simulation. On average this will result in the simulation of a generic, unresolved AGN population. Any possible detection will indicate the viability of this detection method.

(24)

2. Fourier Analysis

One of the most effective mathematical tools for investigating the periodic na-ture of a signal is Fourier analysis. As this type of analysis is essential to the investigation of AGN variability, this section will cover both the mathematical underpinnings and the implementation of Fourier analysis on real astrophysical data. The discussion will start with the basic principles involved in Fourier anal-ysis and proceed with a derivation of the Discrete Fourier Transform (DFT), the primary tool used in the analysis of discretely binned data. The next step is the discussion of the effects of the noise in astrophysical observations as well as an approach to deal with statistical errors during a Fourier transformation. The final topic is the Fast Fourier Transform (FFT), the collection of algorithms with which numerical Fourier analysis is most commonly performed.

2.1 Fourier Series and Fourier Transform

The basic principle underlying Fourier analysis is the notion that any periodic function can be represented as the sum of trigonometric functions. Due to its orthogonality, the trigonometric system can be viewed as a basis for the vector space in which the function is defined. In this manner the analysis breaks a function down into its different periodic components. The most basic form is the Fourier series, which is valid for any periodic function over an interval the length of one period. The Fourier series is the sum over infinitely many possible trigonometric components, each multiplied with a factor that signifies the relevance of the contribution to the periodic function under consideration. As the sum is infinite, all possible periodic components are considered. In the case of a general function f (x) with period p = L_π the relation is given below, where cn is the weight factor for each contributing function. For brevity the

Fourier series is presented in complex form:

f (x) = ∞ X n=−∞ cne 2πinx p _c n= 1 p p/2 Z −p/2 f (x)e−2πinxp _dx ₍₂₎

The approach introduced by the Fourier series, valid for a function over a limited interval (as indicated by the integration limits in the definition of cn), can be extended to the Fourier integral. The Fourier integral contains all

values of n in the complex exponential, not just discrete values, and covers each trigonometric function over its entire range. The complex Fourier integral, for the same general function f (x), is defined as:

f (x) = ∞ Z −∞ ˆ f (w)e2πiwxdw f (w) =ˆ ∞ Z −∞ f (x)e−2πiwxdx (3)

In this case the ‘weight’ of the periodic components is determined by the func-tion ˆf (w). This function is referred to as the Fourier Transform of f (x). The

(25)

2 FOURIER ANALYSIS 2.2 The Discrete Fourier Transform

Fourier integral can therefore be interpreted as two sequential transformations: from f (x) to ˆf (w) and back to f (x). The ability to describe any function com-pletely as the sum of trigonometric functions has led to the most important application of Fourier analysis: the transformation of a periodic function to a function of its periodic components. In terms of bases, the transformation implies a change from the coordinate basis of f (x) to a basis consisting of the trigonometric functions. In other words: the Fourier transform lists the pro-jections of f (x) on the trigonometric basis (analogous to the set cn of Fourier

components in the discrete case). It should be noted that the Fourier transform is a linear transformation from the coordinate vector space to the ‘periodicity’ vector space and therefore is technically not a basis transformation. This dis-tinction does not apply to the discrete Fourier transform, which will be discussed further on.

In the the case of a time series, the Fourier transformation disassembles the signal into its different frequency components. An important tool in this analysis is the power spectrum, also referred to as the Power Spectral Density (PSD), which is defined as the absolute value of the (complex) Fourier trans-form squared. Because the Fourier transtrans-form often has complex components, it is more efficient to use the absolute value of the function as the unit of mea-surement. The power spectrum shows the relative strength of the frequency components. For continuous functions the Fourier transformation and power spectrum are given by (note that in this notation ω is the frequency, related to the angular frequency by ω =ωang

2π ): ˆ f (ω) ≡ Z ∞ −∞ f (t)e−2πiωtdt ⇒ P (ω) = | ˆf (ω)|2 (4)

If the Fourier transform is seen as the spectral density, counting the contribu-tions from each component of the spectrum, the power spectrum can be seen as a measure of the energy contained in the different components. It is therefore a measurement for the energy contained in the periodic nature of the system.

2.2 The Discrete Fourier Transform

In the case of astrophysical observations the data is generally binned, therefore the continuous Fourier transform must be modified to apply to discrete data. This is achieved in the form of the discrete Fourier transform (DFT). The DFT works as follows: if the binned data is considered as the sampling of a continuous signal, which can be periodic or not, the effect of the DFT is to find a complex trigonometric polynomial that interpolates the original signal between the data points. If the data consists of N data points, labelled nj, the following relation

between the data and the polynomial must hold:

nj= p(xj) = N −1

X

l=0

(26)

where cl is the Fourier coefficient for each bin. To find these coefficients, one

can follow a procedure very similar to the derivation of the ‘traditional’ Fourier coefficients, by multiplying the polynomial (5) with e−2πijk/N and summing over j from 0 to N − 1. A short calculation using the orthogonality of trigonometric functions results in the following coefficient:

ck= 1 N N −1 X j=0 nje−2πijk/N (6)

Note the similarity to the continuous Fourier transform. The integral has been replaced with a sum and there is an extra factor 1/N, but otherwise it is identical. In fact, the polynomial p(x) forms the Fourier transform from the frequency space to the time series, in the same way that the Fourier integral in (3) can be seen as the reverse of the Fourier transform f ( ˆw). The complex coefficients ck form the discrete Fourier components in the frequency space. Also note that

in this case the transform truly is a basis transformation, as it projects a series of complex polynomials (5), onto another series of complex polynomials (6). For simplicity the factor 1/N is usually dropped, as it can always be included in any type of normalisation. The exact definition of the unnormalised DFT component for frequency bin k is therefore (Kreyszig 2006; Bracewell 2000):

Ak = N −1

X

j=0

nje−2πijk/N (7)

With the DFT properly defined, it is possible to discuss several important as-pects of this transform. First is its frequency resolution: the level of difference in frequency which can be detected and the manner in which the frequencies are distributed over the Fourier bins. Second is the reaction of the DFT to a pure periodic signal; as it turns out, there is a possibility of some loss of signal strength in the case of signals at frequencies that do not precisely coincide with the binning, this process is referred to as scalloping. The effects of noise and the combination of signal and noise are discussed in the following subsection.

The resolution of the Discrete Fourier Transform is determined by the to-tal observation time and the size of the bins in the time domain. Consider a time series of N bins, each of width dt. The associated DFT will have the same amount of bins. The frequency associated with one particular bin, fk,

is k/(N·dt) = k/T. The equivalence to the continuous Fourier transform can be seen when one considers the term k/N in the definition of the DFT as the equivalent of the term ω in (4). The kth_{bin represents the frequency ω}

k

(mea-sured in units of dt−1) or k/N and has width 1/T. This means that the zeroth bin contains the constant component of the signal, sometimes referred to as the DC-component. The value in this bin is always real and equals the summed signal value over all time bins, as can be seen by inserting the value k = 0 into (7). In the case of photon counting data, such as provided by the LAT, the value of the zeroth bin equals the total number of photons in the observation.

(27)

frequency, with the associated bin k = N/2. At k = N/2 the frequency equals half the sampling rate. For frequencies higher than this value the resolution of the DFT is insufficient to distinguish the signal from lower frequency signals. This effect is called aliasing. The cause of aliasing is that variations in the sig-nal amplitude that occur between two different sampling times are not detected, which implies undersampling of the signal. When a continuous signal is under-sampled the registered fluctuations of the high frequency signal can resemble fluctuations at a lower frequency, as most of the high frequency fluctuations simply are not registered. It is theoretically possible to eliminate this source of error by passing the signal through a high band-pass filter before applying the DFT, this is however not an option for the LAT data.

In the case of the detection of individual photons, the continuous signal of the (astrophysical) source is already sampled by the division of the signal into discrete quanta. The ‘sampling rate’ cannot be controlled in the observation and therefore pre-imposes a limit on the detectable frequencies. The binning of the photons sets a further upper limit. An increased signal amplitude is in-dicated by shorter arrival times between photons. In the case of a single short flare in intensity, undersampling occurs when photons with a small difference in arrival times end up in the same bin as photons that arrived significantly earlier or later. If the frequency of photon arrival times increases, the probability of this frequency being erased by the width of the time bin is also increased. For a continuous signal component at a frequency higher that fmax, undersampling

can lead to aliasing. The effect is strongest for frequencies only slightly higher the Nyquist frequency. For higher frequencies the number of full cycles included in a single time bin is larger, which dilutes the effect of aliasing.

Other than the maximum resolution, the Nyquist frequency has one more important property, particularly relevant for astrophysical observations. For exclusively real (i.e. not complex) entries in the time series, as are common in astrophysical observations, the power spectrum is symmetric around the Nyquist frequency. This can be seen from the following relation:

AN −k= N −1 X j=0 nje−2πij(N −k)/N = N −1 X j=0 nje−2πij· e+2πijk/N = N −1 X j=0 nje+2πijk/N = N −1 X j=0 n∗_j(e−2πijk/N)∗= A∗_k

This property is useful in speeding up the calculation of power spectra, as it implies that one only needs to calculate the Fourier transform up to the Nyquist frequency, reducing the calculation time.

The second aspect to be discussed is the reaction of the DFT to a purely periodic signal (This discussion closely follows section 3.2.1 in Ransom, Eiken-berry, and Middleditch 2002). As will be covered further on in this section, the power in a single DFT frequency bin can be assessed through both an algebraic and a geometric approach. In the case of the response to a periodic function only, represented here by a simple trigonometric function, the geometric

AGN in the Extragalactic Gamma Ray Background : A Possible Detection Method for Unresolved AGN Based on Intrinsic Variability in Sources Observed by Fermi-LAT

MSc Astronomy and Astrophysics

GRAPPA track