Fermi-LAT Gamma-Ray data
analysis of the Supernova Remnant
Cassiopeia A
Author: Orestis Pavlou
Project Supervisor: dr. Jacco Vink
Daily Supervisors: Rachel Simoni, dr. Fabio Zandanel
February 2018
A B S T R A C T
In this work, we performed data analysis of Gamma-ray observa-tional data of the Supernova Remnant Cassiopeia A, provided by the Fermi-LAT satellite mission between August 2008 and March 2016, in the energy ranges between 100 MeV - 100 GeV.
The analysis was done using the Fermi Science Tools (as well as other user contributed tools) to extract the Spectral Energy Distribu-tion of Cassiopeia A (Cas A). A comparison is made with the most recent official releases of data for this Supernova Remnant (SNR) by the Fermi-LAT collaboration in 2013 (1), to show the improvement in statistics.
Model fitting was done using a Python fitting package called Naima, to perform a data comparison between the leptonic and hadronic emission models.
We also performed a theoretical examination on the enhancement factor used by the Fermi-LAT collaboration and present the effect a different nuclear composition of the Remnant and the surrounding environment would have on the total estimated explosion kinetic en-ergy.
D E D I C AT I O N
To my family and friends, for their love and support through everything.
A C K N O W L E D G E M E N T S
First and foremost, I would like to thank my thesis advisor, dr. Jacco Vink of the Anton Pannekoek Institute at the University of Amster-dam, a leading member of the Astroparticle Physics GRAPPA pro-gram. His guidance and input, as well as his mentorship, were of utmost importance to the inception, development and continual evo-lution of this work.
I would also like to thank my daily supervisors, Rachel Simoni and dr. Fabio Zandanel, for the unwavering daily support and guidance throughout the development of this work. Their passionate dedica-tion while assisting me on this project has been invaluable and they have provided me with an experience I will always carry with me in the future.
I would also like to acknowledge the importance of the input that the entire research group of dr. Jacco Vink has given me, through regular meetings where the exchange of ideas and research develop-ments have helped me throughout my entire research project.
I also am very thankful to dr. David Berge, as the second reader and referee of the present thesis, and I am grateful for his additional input and assistance.
Finally, I would like to express my unending gratitude and love towards my entire family and friends, whose support, caring and encouragement throughout my entire studies have been the corner-stones of what I was able to achieve the past few years.
Author Orestis Pavlou
C O N T E N T S
1 i n t r o d u c t i o n 15 1.1 Cassiopeia A 16 1.2 Cosmic Rays 17
1.3 Gamma-ray production 21 1.4 Diffusive Shock Acceleration 23 1.5 SNRs as CR accelerators 25 1.6 Fermi - Large Area Telescope 27
1.7 Fermi-LAT Instrument Response Functions 30 1.7.1 Performance Plots 31
1.7.2 Sensitivity Plots 33
1.7.3 Low Energy Performance Plots 35 1.7.4 Pass 8 versus Pass 7 37
2 a na ly s i s m e t h o d 39 2.1 Data Selection 39 2.2 Likelihood statistics 41 2.3 Likelihood analysis 42
2.4 Model and Residuals Maps 44 2.5 Spectral Energy Distribution 46 3 d i s c u s s i o n 49
3.1 Analysis Results 49 3.1.1 Errors 49
3.1.2 Diffuse background and residuals 51 3.1.3 Model Fitting 53 3.1.4 Hadronic Model 54 3.1.5 Leptonic Model 57 3.1.6 Conclusions 59 Appendices 63 a e n h a n c e m e n t f a c t o r 65
a.1 1st estimation of the enhancement factor 65 a.2 2nd estimation of the enhancement factor 66
L I S T O F F I G U R E S
Figure 2 Chandra image of Cas A (9), taken between 2004 and 2005 using the Advanced CCD Im-age Spectrometer (ACIS). Energy ranges: Red = 0.5−1.5 keV, Green = 1.5−2.5 keV, Blue = 4.0−6.0 keV. 17
Figure 3 A sample of telescopes operating in all wave-lengths by 2013 (11). 19
Figure 4 Updated Cosmic Ray Spectrum (2016), (12). For a detailed explanation of the features, see text. 20 Figure 5 Above: Compton Scattering, resulting in a high
energy electron. Below: Inverse Compton Scat-tering, resulting in a high energy photon. (13) 21 Figure 6 Bremsstrahlung: Radiation (photon) is emitted
after the deflection of a moving electron by a stationary positively charged particle. (14) 22 Figure 7 Synchrotron Radiation: Radiation (in this case
radio waves) is produced by the helix orbit of a charged particle around a magnetic field line. (15) 22
Figure 8 Pion decay: A neutral pion decays into to gamma-rays (16). 22
Figure 10 Non-linear shock acceleration: Instead of 2 re-gions (shocked and unshocked), we now have 3 regions by adding the ”sub-shock” between the two. 25
Figure 11 Concept of a gamma ray producing a pair that passes through the LAT. 29
Figure 12 The modules that constitute the LAT detec-tor. 29
Figure 13 Effective area as a function of energy for inci-dent photons with θ = 0 (called ’normal inci-dence photons’). 31
Figure 14 Acceptance of P8R2 SOURCE V6, as a func-tion of energy for the total amount of events as well as the FRONT and BACK event types. 32 Figure 16 Acceptance weighted energy resolution as a
function of energy, showing the difference in quality between each EDISP event type. 33 Figure 17 5-year exposure map (projected as a 10-year map), not including the modified observation strategy during Dec 2013-Dec 2014. 34
Figure 18 Differential sensitivity for the projected 10-year map in 4 separate bins per decade of energy between the energies of 31.6 MeV - 1 TeV. 34 Figure 20 Effective area versus energy, down to 10 MeV. 36 Figure 21 Acceptance versus energy, down to 10 MeV. 36 Figure 22 Energy resolution versus energy, down to 10
MeV, for FRONT/BACK and TOTAL event types. 36 Figure 25 Counts Map created using energies of 100 MeV
- 100 GeV in a 10-degree radius around Cas A. The numbers on each axis are the number of pixels of the image (100×100). The colorbar on the right side shows the photon count for each pixel. 40
Figure 26 Likelihood fit plot of all the sources in the source file. The top figure is the plot in colour, with the colours labeled for each source on the right hand side. The bottom figure demonstrates the errors for each of the 30 Energy bins, in com-parison with the best fit line. 43
Figure 27 Model Map for the area around Cas A. This is a heat map, with Cas A in the center of the figure and the Perseus arm emission notice-able in the lower-right part. The numbers on each axis are the number of pixels of the im-age (100x100). The colorbar on the right side shows the photon count for each pixel. 44 Figure 29 SED (E2×Differential Flux in each energy bin)
for our analysis (blue data points) and for the Fermi 2013 paper (red data points). Our errors our computed with 1σ standard deviation. 48 Figure 30 Comparison between the acceptance of the
lat-est version of Pass 7 to the version of Pass 8 used in the present work. Improvements of ≈25% are noticeable at energies ranging from of a few hundred MeV up to a few GeV. 50 Figure 31 Comparison between the broadband sensitiv-ity of the latest version of Pass 7 to the version of Pass 8 used in the present work. Improve-ments of ≈ 50% are noticeable at energies of a few hundred MeV to ≈ 30% at GeV ener-gies. 51
Figure 32 Composition map of the Perseus arm taken from Ungerechts et al. (2000) (31). There are two areas (D & E) that show an excess of CO and13CO cloud density and also several gaps in cloud density around Cas A. 52
List of Figures 13
Figure 33 Hadronic Model Fit: E2xDifferential Flux for
each energy bin for a pion decay radiative model. 55 Figure 34 Hadronic Model Fit: E2xDifferential Flux for
each energy bin for a pion decay radiative model with the minimum proton energy set at 10 GeV. An improvement in agreement between data and model is observed. 56
Figure 35 Leptonic Model Fit: The E2x Differential Flux for each energy bin for a leptonic radiative model which includes Inverse Compton scattering, Syn-chrotron and Bremsstrahlung radiations. 58 Figure 36 Inverse Compton Model Fit: The E2x
Differ-ential Flux for each energy bin for a separate Inverse Compton radiative model, with a spec-tral index of 2.6. 58
1
I N T R O D U C T I O N
In 2013 the Fermi-LAT collaboration published an intensive study of the GeV emission of Cassiopeia A: ”the detection of the GeV break in the gamma-ray spectrum of the SuperNova Remnant Cassiopeia A”. This study used 5 years of data and the reconstruction package Pass 7 for analysis. The goal of the present work is to pursue a similar analysis with 8 years of data and using an updated reconstruction package Pass 8. We discuss our results in the perspective of former studies and provide a deep study of the spectrum using models generated by the Naima research tools.
We will begin our introduction to the current project by providing a historical overview of the Supernova remnant Cassiopeia A in section 1.1, followed by an introduction to the detection of cosmic rays in section 1.2, their production and emission processes in section 1.3, as well as a short discussion on Diffusive Shock Acceleration in section 1.4 and Supernova remnants as cosmic ray accelerators in section 1.5. We conclude our introduction chapter by presenting a few details about the Fermi satellite mission, the Fermi-LAT instrument and its Instrument Response Functions in sections 1.6 and 1.7.
In section 2 we present our analysis procedure in depth, starting from our data selection and subsequently showing our likelihood statistics and analysis, the creation of our maps and the Spectral En-ergy Distribution (SED) we generated. A second analysis was also performed, using the user contributed tool called FermiPy.
In section 3 we discuss our analysis results by comparing our out-put to that of previous papers published by the Fermi-LAT Collabo-ration, to point out improvements due to the use of 8 years of Fermi data of Cas A and Pass 8. Finally, we present the procedure and results of our model fitting using the Naima tools, followed by a detailed comparison between the suggested hadronic and leptonic emission models and their efficiency, as well as a short discussion on the enhancement factor and its effect on the proposed composition of Cassiopeia A.
1.1 c a s s i o p e i a a
Cassiopeia A is one of the youngest and most unique known Super-nova Remnants (SNRs). It was first recorded around 350 years ago and has been a very intriguing celestial object, rigorously and regu-larly studied during the past decades. Cas A lies ≈ 11kly (3.4kpc) away, with galactic coordinates: l= 111.7◦, b = −2.1◦ (2) and is con-sidered to be the brightest source of radio waves outside the solar system (3). It is also a source of non-thermal emissions of relativis-tic parrelativis-ticles ranging from radio and X-ray up to TeV gamma rays, originating from the forward and reverse shocks (4), (5).
The original star of Cas A is estimated to have been around 15 -20 solar masses. This star went supernova around 1680, as records of a star being in that location were found only up to that year. The Astronomer Royal John Flamsteed is thought to be the first to discover it, as he was the first to catalogue a star near the position of Cas A on August 16, 1680 (6). It is also the only confirmed Galactic Supernova to be visible by the naked eye to date! The relatively young age of Cas A and its short distance from us have been helping us understand the anatomy of a supernova explosion, as it can be observed by many instruments and at many relatively early stages of its expansion. Cas A was first observed among a list of the earliest discrete radio sources by radio astronomers in Cambridge, in 1947 (7) and optically detected in 1950 (8).
In 2006, Cas A was used to study how cosmic rays are accelerated by a supernova explosion by Chandra, in X-rays (9). These studies provided the first hints that electrons, as well as protons and ions, can be accelerated up to energies where they escape from the SNR bubble, thus offering some of the first evidence that SNRs can be the key source of the most energetic Galactic cosmic rays.
Scattering inside the shock front and disentangling by strong mag-netic fields around the remnant generates these highly energetic cos-mic rays, through some physical processes that will be discussed and overviewed in sections 1.2 - 1.5 of this introduction. Figure 1 is a deep Chandra image from 2006 (9), one of the first highly detailed images of Cas A, which shows different parts of the SNR in different X-ray energy bands. In red and green, the lowest energy X-rays are represented, which are emitted from a big part of the ejecta of the SN explosion. At the outer edge, particles from the surrounding interstel-lar medium are accelerated to very high energies, after they cross the shock front back and forth several times, emitting the non-thermal X-rays shown in blue in this figure (9). This high energy emission is caused by the expansion of the outer shell of the SNR called the shock front.
1.2 cosmic rays 17
Figure 2.: Chandra image of Cas A (9), taken between 2004 and 2005 using the Advanced CCD Image Spectrometer (ACIS). En-ergy ranges: Red = 0.5−1.5 keV, Green = 1.5−2.5 keV, Blue = 4.0−6.0 keV.
1.2 c o s m i c r ay s
Since their first measurement, in 1912 by Victor Hess, cosmic rays have been occupying the forefront of astrophysical studies. Using balloons, Hess detected a ”radiation of very great penetrating power” entering our atmosphere ”from above”. It was later discovered that these cosmic rays are charged particles, which penetrate the atmo-sphere, interacting with particles in it and generate particle showers. Galactic cosmic rays consist of electrons, protons and nuclei, with velocities ranging from 0.45c - 0.996c (where c is the speed of light). Whereas low energy cosmic rays are abundant and easily detected by balloons and satellites, very high energy cosmic rays are very rare, as shown in Figure 3 in page 19, which translates to a low number of detections of particles per year. In addition to their rarity, these high energy cosmic rays are harder to detect.
While the nature and sources of low energy cosmic rays is well un-derstood, high energy cosmic rays pose questions about their origin since their discovery. Distinctly significant is the flux of cosmic rays at energies around 1015 eV. Although the cosmic ray flux decreases drastically for higher energies, around this particular area the slope
of the flux decrease changes. Thus there is a small bump in the spec-trum, which is called the ”knee”. Around 1018 eV, the slope shape changes again and this is called the ”ankle”. The leading theory ex-plaining these phenomena is that in this energy range, between the ”knee” and the ”ankle”, the dominant sources responsible for the cos-mic rays switch from Galactic to extra-Galactic. Below the ”knee” the dominant sources are thought to be of galactic origin, and more specifically SNRs which are proven to be CR sources, while above the ”ankle” the spectrum is entirely dominated by the emission of extra-galactic sources, for example Active Galactic Nuclei (AGN) and Gamma-ray Bursts (GRBs). The reasoning behind this is that only extra-galactic sources have been detected to create CRs of those en-ergies and nothing locally (in our Galaxy) has ever been observed to create CRs with energies above the ”ankle” (1018eV). Thus the area in the spectrum between 1015 and 1018 eV is where the galactic sources, which are more abundant producers of CRs that reach us in ener-gies below 1015eV, start losing their dominance and the extra-galactic sources begin to dominate.
Another interesting feature of the cosmic ray spectrum, is the de-tections at ultra-high energies, around 1019- 1020eV. Special relativity predicts a limit for sources emitting at these energies around 1019 eV, because particles with such high energies would interact with the Cosmic Microwave Background (CMB) photons and could not travel further than 150 million light-years from the source. Therefore detec-tions of CR particles above these energies means that they couldn’t be produced farther than 150 million lys away from us. This is called the ”Greisen-Zatsepin-Kuzmin cutoff” or ”GZK cutoff” for short. Since CRs have been detected multiple times (although very rarely, around 1 particle per km2per century), there is difficulty justifying their extra-Galactic origin, because of this limit in distance, although there are hundreds of Galaxies within this limit.
To determine the sources of creation and emission of these CR par-ticles, charged protons and nuclei are poor messengers, as they are deviated by the galactic magnetic field. Neutral secondary particles produced at the CR acceleration site, as gamma-rays, neutrinos, or high energy X-rays, are used to examine CR astronomy. A combina-tion of satellite telescope detectors and ground-based telescopes are used to examine these secondaries, with the former mainly detecting X-rays and gamma-rays, whereas ground-based telescope array tectors, which consist of Cherenkov telescopes, are able to make de-tections on gamma-ray particles of energies around 1 TeV (10). There are certain limitations encountered by ground-based telescopes com-pared to space-based ones, such as being limited by the large isotropic background from cosmic-ray showers, created when CRs meet atmo-spheric particles. Cherenkov radiation observatories can avoid the latter. On the other hand, space-based detectors are much more
ex-1.2 cosmic rays 19
pensive to launch and maintain, much more difficult to repair and the higher energy threshold of the ground-based detectors means that lower energy detections are eliminated by cascades high up in the atmosphere, and thus don’t interfere as much with proper de-tections. Satellite telescopes also offer a much smaller effective area of simultaneous observations compared to ground-based Cherenkov telescopes. Figure 2 shows roughly the energy dominant ranges for many space-based and ground-based detectors operating as of 2013 (11).
Figure 3.: A sample of telescopes operating in all wavelengths by 2013 (11).
SNRs have been observed to accelerate cosmic rays to very high en-ergies via a process called ”Diffusive Shock Acceleration” or ”Fermi Acceleration”. This process is described in detail in section 1.4 of this chapter. The basic idea behind this process is that particles in the surrounding medium around SNRs are swept by the expanding front of the shocked medium. These particles are accelerated by the kinetic energy of the shock deposited in them, until they escape the shock front and re-enter it at a later time. This back and forth in and out of the shock front (called crossings) results in some of these particles gaining enough energy to escape from the system outwards. These high energy particles are what we detect and study as cosmic rays.
Figure 3 shows a plot of the cosmic ray spectrum versus the energy, with detections from various experiments.
Figure 4.: Updated Cosmic Ray Spectrum (2016), (12). For a detailed explanation of the features, see text.
1.3 gamma-ray production 21
1.3 g a m m a-ray production
The LAT instruments detect and reconstruct the paths of high energy incident particles and their products inside the detector. These parti-cles are results of radiative processes happening in the CR accelera-tion sites. The main radiative processes considered to be responsible for the production of the incident CR particles are:
1. Inverse Compton Scattering (IC)
2. Synchrotron Radiation / Bremsstrahlung 3. Pion Decay
Inverse Compton Scattering is the reverse process of Compton Scat-tering, as suggested by its name. Compton Scattering is the leptonic process during which a high energy photon interacts with a lower energy electron producing a high energy electron. During Inverse Compton scattering, a population of high energy electrons deposits its energy into a population of lower energy photons, which escape with high energies (see Figure 6).
Figure 5.: Above: Compton Scattering, resulting in a high energy elec-tron. Below: Inverse Compton Scattering, resulting in a high energy photon. (13)
Synchrotron Radiation is a radiative process during which charged particles are emitted from a source, after interacting via
Bremsstrahlung with the magnetic field. Bremsstrahlung is the interac-tion between two charged particles, in which a moving charged parti-cle interacts via its charge with a stationary charged partiparti-cle resulting in the change of direction and momentum of the moving charged par-ticle. This momentum change from the charged particle, produces a photon with this excess momentum. In Synchrotron Radiation, the incident moving charged particle interacts with a magnetic field line in such a way that it can follow a helix path around it, constantly emitting photons during this movement.
Both these processes are shown Figures 7 and 8 (14), (15).
Figure 6.: Bremsstrahlung: Radiation (photon) is emitted after the deflection of a moving electron by a stationary positively charged particle. (14)
Figure 7.: Synchrotron Radiation: Radiation (in this case radio waves) is produced by the helix orbit of a charged particle around a magnetic field line. (15)
Pion Decay in astrophysics, usually refers to the physical process of the decay of a neutral pion (π0) into two photons. Pions consist of a
quark and an anti-quark (mesons). They are unstable, with a decay time for charged pions being in the order of nanoseconds (26×10−9s) and for the neutral pions being much shorter (8.4×10−17s). Charged
pions, π+,− decay into muons and muon neutrinos/anti-neutrinos, whereas neutral pions decay into two gamma-rays. Pions are the products of high energy collisions between hadrons (i.e. high energy protons or nuclei).
1.4 diffusive shock acceleration 23
1.4 d i f f u s i v e s h o c k a c c e l e r at i o n
The main idea of diffusive shock acceleration originates from Enrico Fermi, who, in 1949, suggested that particles can be accelerated to higher energies (or lose energy) after bouncing around moving clouds. If the cloud is moving towards the particle, which means it has kinetic energy to offer to the particle, then the particle will be accelerated to higher energies. Since the cloud is moving through a medium, then the likelihood of particles gaining energy is higher than them losing energy. Since the cloud, which in the case of SNRs is the cloud of the shock expanding, deposits huge amounts of its kinetic energy into the surrounding particles, then Fermi argued that these particles could theoretically be accelerated to very high energies. This process is also called ”Fermi Acceleration”.
Shocks can be encountered in many astrophysical sources, i.e. ac-cretion shocks caused by compact objects, shocks caused by solar wind and coronal mass ejections, by stellar winds and supernova rem-nants, as well as relativistic shocks created by Active Galactic Nuclei (AGN), Gamma-Ray Bursts (GRBs) or Galaxy clusters.
Using first order Fermi Acceleration and assuming the notation shown in the following images, in Figure 8, we can derive the shock equations.
(a) Regions (shocked and unshocked) around the shock, from the frame of the shock.
(b) Regions (shocked and unshocked) around the shock, from the frame of an outside observer (particle frame).
The shock equations, also known as Rankine-Hugoniot equations, that are derived using the above notation, are presented below:
Mass-flux conservation:
ρ1υ1 =ρ2υ2
Pressure equilibrium at the shock: P1+ρ1υ21 =P2+ρ2υ22
Enthalpy/Energy-flux:
Using the equations for the dimensionless quantities of the com-pression ratio: χ = ρ2
ρ1 =
υ1
υ2 and for the Mach number: M =
υ
cs
, cs =
q
γPρ, we can solve these conservation equations to get a
quadratic equation of dimensionless quantities: [( 2γ γ−1) 1 γM2 +1]χ 2− ( 2γ γ−1)(1+ 1 γM2)χ+ [( 2γ γ−1) −1] =0
This equation has one trivial solution and the other solution (shock so-lution) is: χ= (γ+1)M2
(γ−1)M2+2. For strong shocks, where the Mach number
goes towards infinity (M→∞), this solution becomes: χ→ γ+1
γ−1. For
relativistic gas a conventional γ= 43 is used, whereas, for monatomic gas γ = 53 is used. These give χ = 7 and χ = 4 respectively. The latter solution for χ will be used to extract the value of the particle spectrum q in the particle distribution equation: dN(E)/dE'E−q.
For χ=4 for strong shocks, we have an expected particle spectrum of a slope q= (1+ ( 3
χ−1)) =2 and this gives a particle distribution of
dN(E)/dE' E(1+(χ−31)) =E−q =E−2.
We can also calculate the mean free path of the particles in the shock, to determine if we can justifiably assume that the shocks are collisionless. Due to the low density of astrophysical shocks, the mean free path for a proton traveling in a shock that moves with a velocity of v ≈ 1000km/s and a column density of 1 particle per cubic centimeter, we find a mean free path: λp ≈1020n−p1 cm, which
is larger than the size of most supernovae detected. We can therefore assume collisionless shocks, meaning that the shocks don’t ”heat up” due to particle-particle collisions, but only from other effects, such as electro-magnetic waves. Collisionless shocks have been directly observed (17).
Today the most popular theory for particle acceleration in shocks, is non-linear shock acceleration. According to this theory, concentrated accelerated particles push against the plasma, compressing and pre-heating it. This results in a lower Mach number and a mix of relativis-tic and non-relativisrelativis-tic parrelativis-ticles in this additional region, called the sub-shock, that lies between the shocked and unshocked regions. This is considered to be a more realistic scenario than the linear (”step”-like) acceleration, for a few reasons. First of all, a lower Mach num-ber is present at the sub-shock due to lower compression there. This leads to a much higher compression overall, due to escape. Also, the spectrum is now curved and not a pure power-law. A steeper spec-trum detected in young SNRs like Cas A, is another supporter for non-linear acceleration, as it suggests lower compression ratios at the sub-shock than linear shock acceleration (χ < 4). Ejecta observed close to the shock front, imply higher overall compression. This cou-pled with the fact that the temperatures observed in some SNRs are lower than expected by linear shock acceleration (18), (19), make the case for non-linear shock acceleration even stronger. Figure 11
illus-1.5 snrs as cr accelerators 25
trates the 3 regions of non-linear shock acceleration: the shocked, unshocked and sub-shock regions.
Figure 10.: Non-linear shock acceleration: Instead of 2 regions (shocked and unshocked), we now have 3 regions by adding the ”sub-shock” between the two.
1.5 s n r s a s c r a c c e l e r at o r s
As mentioned earlier, Supernova Remnants are CR acceleration sources to energies up to the ”knee” (≈1015 eV). As discussed in the Gamma-ray Production section, CRs can be accelerated to produce Gamma-Gamma-ray radiation with both leptonic and hadronic processes, namely via IC scattering, Synchrotron via Bremsstrahlung interaction with a strong magnetic field and via pion decay. Although distinguishing between these emission processes has been difficult, direct evidence for the ac-celeration of CR protons from SNRs, with the observation of the pion decay feature, have been obtained (20).
Particles are gaining energy by scattering elastically each time they enter the shock. Performing these ”crossings” in and out of the shock, enables the particles to accelerate to very high energies after n-crossings, given by the relation: E= E0(1+2α)n. As far as the
lep-tonic emission processes are considered, the radio index for SNRs has been measured to be around: α ≈0.5(0.3−0.8)(21). The radio signal indicates that SNRs emit Synchrotron radiation. The radio index (α) correlates to the particle index (q), with the relation: q = 2×α+1.
This means that the detected radio index of ≈ 0.5 corresponds to a particle index of q ≈ 2.0, which is close to the detected cosmic ray particle index of q=2.7 (22).
For synchrotron radiation, after taking into account synchrotron losses, the total energy of the particles and the magnetic field is
esti-mated by the following relation: Wtotal =ηVREEmax
min κEN(E)dE+V B2
2µ0, where η = correction for protons,
κ = ( 1
1−B2 E636(t)t)
−q+2 = a factor that accounts for synchrotron losses,
N(E)dE = κE−q = electron energy distribution, B = magnetic field
and µ0 = magnetic constant. After some calculations it can be shown
that: Wmag= 34Wparticles.
A very important step towards explaining the cosmic ray spectrum, using galactic SNRs as sources, came in the 1990’s in the form of X-ray data. Up until the 1980’s the astrophysical community was not at all convinced that SNRs could be responsible for accelerating particles up to the ”knee”, based solely on radio data of synchrotron emission. SNRs are estimated to deposit 5−10% of their total energy into CRs, and this is the percentage that CRs needed to be found to be accelerated to, around SNRs, if their case could be supported.
Using estimations for the magnetic field of B≈5µG and turbulence
η≈ 1, researchers found it unlikely for SNRs with V ≈ 5000km/s to
accelerate particles up to the required energies. For this to happen, basically either the magnetic field should be orders of magnitude higher than this, and/or the turbulence around the SNR should be very high. X-ray observations (first in 1995) implied the presence of electrons around SNRs with energies up to 10−100 TeV. Also, the fact that some young SNRs are emitters of such high energy elec-trons, suggests high turbulence (1 < η < 20), high magnetic fields
(B ≈ 10−100µG) and that electrons don’t need a lot of time to be accelerated to high energies (this can happen in a few decades-centuries). This also implies particle acceleration to energies of the order of 1014 eV. These data then present a very good argument
sup-porting some SNRs as accelerators of CRs up to 1015as it seems a lot more plausible, with the observed higher magnetic fields and turbu-lence.
1.6 fermi - large area telescope 27
1.6 f e r m i - large area telescope
One of the most modern tools obtaining observations on galactic gamma-rays is the Fermi Gamma-ray Space Telescope (FGST) satel-lite. This satellite was launched on June 11th, 2008, to an initial alti-tude of ≈ 565 km and has been performing observations since June 24th, 2008. The Fermi team claims that it is able to obtain measure-ments with a sensitivity factor 30 times better that of any other earlier gamma ray missions.
The missions lifetime was initially approved for 5 years and later extended to 10 years, in August 2013. This means the satellite and instruments will remain operational until 2018. The schedule of the mission was separated into three main phases: Phase 0 which repre-sents the initial in-orbit check of 2 months, then Phase 1 which consists of the first year ’verification’ period during which the entire sky was first surveyed and finally Phase 2 which is the main body of the mis-sion and consists of years of observations, until its eventual mismis-sion ending. The gamma-ray burst data were made available after Phase 1, whereas the LAT data and the corresponding Science Tools were made available after Phase 2 began (14 months after launch).
One of the main instruments of the Fermi satellite is the Large Area Telescope (LAT), which performs high sensitivity long term gamma-ray observations of the entire sky, studying high energy astrophysical phenomena, in the energy range between 20 MeV and>300 GeV (23). The other instrument of the Fermi satellite is the Gamma-Ray Burst Monitor (GBM), which covers a lower energy range than the LAT, between 10 keV - 25 MeV. This was also the instrument that recorded the electro-magnetic counterpart of the gravitational wave from the two merging Neutron stars in August 2017, published by the LIGO Scientific Collaboration and VIRGO Collaboration in October 2017 (24). Both of these instruments have very wide fields of view (FOVs). For more than 80% of its lifetime, the LAT is in the scheduled sur-vey mode, but the instruments can either be commanded to inves-tigate a particular part of the sky or, based on data obtained from the instruments themselves, request autonomous re-pointing of the spacecraft to point to a certain target. In addition, the observing mode can be switched to observe short events in or near the FOV of the LAT. Afterwards the scheduled survey mode is resumed. A window of five hours is currently provided for such automated re-points. The pointing itself has an accuracy of< 2 degrees (with a 1σ significance) and a relative pointing error of<10 arcsec.
In this essay we focus only on the LAT instrument, as the selected energy range of the obtained data is between 100 MeV - 300 GeV. The LAT has multiple uses, such as localizing point sources, extended sources’ spatial mapping and spectra, measuring the diffuse gamma-ray background (DGRB) as well as the flux of cosmic gamma-ray electrons.
It is also being used to send notifications about gamma-ray transient events quickly and create catalogues of gamma-ray sources.
The LAT’s detector uses pair conversion: When gamma-rays enter the detector, they interact with tungsten and produce an electron-positron pair (γ → e−+e+). Whereas the energy range of 1−30 MeV is observed by Compton scattering detectors, which depend on the Compton scattering of a photon on a charged particle in the de-tector, energies above 30 MeV are dominated by pair production. To achieve this, pair production detectors are layered with converters, which consist of heavy metals. These heavy metals are the target for the incident gamma-ray photon. The photon interacts with the heavy metal and produces the initial positron-electron pair. The electrons of this pair then ionize the gas in the chamber, which in turn leads to the triggering which gives an electrical signal. Interpreting this signal permits a reconstruction of the path of the pair through the layers of the detector. We know the x-y positions where the pair particles hit each layer and adding the layers together provides the 3rd dimension (z). Thus a three-dimensional picture of the path the pair travelled through can be obtained, which can help point to the direction the incident gamma-ray photon hit the detector. This, in addition to the amplitude and properties of the signal itself, provides the position in the sky the gamma-ray photon came from, as well as it’s initial energy when it hit the detector (25).
The silicon strip detectors in the LAT track the pair, while it travels with the direction of the incident gamma ray. This is because the gamma ray’s energy is many times higher than the energy of the pair produced. The pair finally enters the lowest part of the LAT, which contains a calorimeter made of 1536 Caesium Iodide - CsI(Tl) crystals, arranged as 8 rows of 96 crystals, which provide an energy measurement of the pair.
Background events are rejected by the LAT, using an anti-coincidence detector (ACD), which consists of 89 plastic scintillator tiles and is placed around the LAT. The ACD is made out of scintillator tiles, that detect these background events and issue a ”veto signal”. Pulse height signals are produced by the pair passing through the tracker and the calorimeter and by combining the x-y coordinates of each of the silicon strip detectors that are hit, the particle trajectory and energy losses are reconstructed and determined.
A detailed schematic of the structure of the inside of each tracker module is shown in Figure 10.
The LAT consists of 16 tracker (TKR) modules placed next to each other on a 4x4 plane, with 16 calorimeter (CAL) modules underneath them, and all are surrounded by an ACD shield. The TKR modules are made of 18 XY tracker planes, which are tungsten plates used to track the path of the positron-electron pairs created by the gamma-ray hitting the detector within it’s field of view and have an argamma-ray of
1.6 fermi - large area telescope 29
Figure 11.: Concept of a gamma ray producing a pair that passes through the LAT.
silicon-strip tracking detectors (SSDs) which detect charged particles. The SSDs are then used to track the initial direction of the pair. Pairs at low energy present a challenge as they result in multiple scatterings in the conversion plain, limiting the accuracy of the calculation of the incident direction. Also, the identification of the particles that result in each interaction is key to distinguishing between pair scatterings and incident cosmic rays.
The LAT’s Data Acquisition System (DAQ), which is placed under-neath the calorimeter modules, is used as a filter that reduces the background events and transmits to the ground. The DAQ is also in charge of controlling, commanding, monitoring and housekeeping of the modules using software that can be uploaded for updating. The average bit rate of the DAQ’s processing is 1.2 Mbps. Figure 11 shows the general structure of the LAT, presenting the place of each module.
Figure 12.: The modules that constitute the LAT detector. The LAT analyzes the events created by incident particles, reducing the raw trigger rate (around 10kHz), down to≈400 events per second
(≈ 400Hz). Only ≈ 2−5Hz of these represent actual photons. The LAT stores data packets of variable length for these events, to the spacecraft’s Solid State Recorder, which are then transmitted to the ground.
1.7 f e r m i-lat instrument response functions
In the following paragraphs we will detail the performance of the Instrument Response Functions (IRFs) of the Fermi-LAT while also presenting some performance plots on the effective area, the accep-tance, the Point Spread Function (PSF), the energy resolution and the sensitivity. We will also show comparisons between the Pass 7 and Pass 8 performance to clarify the improvements between the two re-leases.
At the time of the launch of the satellite, the initial calculations of the instrument underestimated the importance of the effects from the interactions between background sources and gamma-rays. These are called ”ghost events”. These events should be identified and rejected by the instrument before calculating the incident particle’s track. A ghost particle passing through the detector triggering a significant amount of trackers, at the same time the incident particle passes through, will interfere with the calculation process. These ghost events can lead to miscalculating the incident particle’s direction or energy, or even be responsible for valid detections being rejected by the Monte Carlo simulation’s event classification algorithms. These ghost events turned out to be much more important than previously thought, which lead the team to produce updated versions of IRFs. Pass 6 or P6 v3 was the name of the first package to be released by the Fermi team, shortly after the launch, in 2010 (26). This package took in-orbit effects into account, to correct for the pile-up or coincident accidental effects. It was later updated by the version P6V11. A more effective and important update was Pass 7, which was made public in August 2011 (27), was the package used by the Fermi Collaboration in their 2013 paper (1), which was the publication that prompted the present work.
The most up-to-date performance package is the Pass 8 analysis framework (28) which, like Pass 7 before it, takes into account and im-proves upon measurements made in-flight, something lacking from the initial setup at launch. Additionally, the aforementioned ghost events are handled differently by this release: Ghost events that co-incide with proper triggering events in the detector, are kept and analyzed, resulting in many more events being recovered that were previously lost due to this effect.
1.7 fermi-lat instrument response functions 31
1.7.1 Performance Plots
The detections of gamma-rays by the instrument are divided into separate event classes which are defined by sets of event selections and have unique sets of corresponding response functions. As was already mentioned, the LAT has 16 tracker (TKR) modules, 12 lay-ers of which have 3% radiation length tungsten convertlay-ers, referred to as the FRONT section of the module, and the rest 4 layers con-sist of 18% radiation length tungsten converters, referred to as the BACK section. The FRONT section has approximately 2 times bet-ter point spread function (PSF) than the BACK section. The LAT point spread function is defined to be a function of the energy of the in-cident photon and its inclination angle. The PSF is represented by: x = Sδ p
p(E), where δp = 2sin
−1(p0−p
2 ) is the scaled-angular deviation
and Sp(E) =
q
[c0(100MeVE )−β] +c21 is called the scale factor. For the
latter, c0 and c1 are parameters in radians which are given distinct
values for each event type, whereas the β parameter is given a set value for all event types of β=0.8.
Two new forms of categorization introduced by Pass 8 are the PSF event types and the EDISP event types. The PSF event types indicate the quality of the reconstructed direction, by dividing the data into quartiles, the lowest quality of which is PSF0 and highest is PSF3. The Energy Dispersion (EDISP) is defined as the same as PSF but is used to represents the quality of the energy reconstruction.
In the following pages, we present the plots provided by the Fermi Collaboration, created using the LAT science tools. The plots consist of 4 points per energy decade, connected for clarity, as stated by the team. The plots are self-explanatory with information provided in the captions, but clarifications will be made on certain plots, when deemed necessary. Figure 12 shows the effective area of the LAT in-strument, which is the area with acceptable resolution that can be used for calculations, as a function of the incident photon energy.
Figure 13.: Effective area as a function of energy for incident photons with θ =0 (called ’normal incidence photons’).
Acceptance here represents the effective area integrated over the solid angle. It is noticeable that acceptance is ”slower to turn on” as stated by the team, which means it becomes as significant at higher energies than the effective area, something that is caused by the de-pendence of the field of view on the energy. Effective acceptance is derived after taking into account certain details of the strategy used to take observation, such as downtimes of the instrument and the inability to take data at the South Atlantic Anomaly.
Figure 14.: Acceptance of P8R2 SOURCE V6, as a function of energy for the total amount of events as well as the FRONT and BACK event types.
Using Monte Carlo simulations, the PSF is derived, including a dependence on the off-axis angle. The following plots in Figure 14 show the containment angles, at 1σ and 2σ significance respectively, for the different event types of the PSF.
(a) PSF at 68% containment for each of the PSF event types.
(b) PSF at 95% containment for each of the PSF event types.
1.7 fermi-lat instrument response functions 33
Figure 15 displays the energy resolution (acceptance weighted), at 1σ significance, for each energy dispersion (EDISP) event type.
Figure 16.: Acceptance weighted energy resolution as a function of energy, showing the difference in quality between each EDISP event type.
1.7.2 Sensitivity Plots
Point source sensitivity is also improved with the Pass 8 update. Fig-ure 16 shows the 10-year integral flux sensitivity for an isolated point source for P8R2 SOURCE V6. Integral flux sensitivity is the mini-mum flux above 100 MeV, that after 10 years of observation obtains the 5σ detection with a power law spectrum index = 2. The 5-year exposure map, obtained from August 2008 until August 2013, is mul-tiplied by 2 to generate a 10-year projected survey-mode exposure map. The loss of observing time caused by South Atlantic Anomaly (SAA) passages of the satellite and the downtime of the instruments (except from a modification in observing strategy that took place dur-ing December 2013 - December 2014) are taken into account. Durdur-ing that year, beginning from December 4th 2013 until December 4th 2014, the Fermi-LAT instrument operated in a ”modified observing mode that roughly doubled the rate of increase of exposure in the galactic center relative to normal survey mode” (29). This observational time is not taken into account. Point source effects are ignored and only diffuse backgrounds are considered for the calculation and creation of this map (Figure 16). The calculation uses gll iem v06.fits for
the galactic IEM template and iso P8R2 SOURCE V6 v06.txt for the isotropic template, which are explained in more detailed in the Anal-ysis Method chapter.
Figure 17.: 5-year exposure map (projected as a 10-year map), not including the modified observation strategy during Dec 2013-Dec 2014.
Figure 17 shows the differential sensitivity of P8R2 SOURCE V6 for ten years between 31.6 MeV and 1 TeV. A point source with a power-law spectrum with index 2 is assumed, with a uniform background around it. The four different curves represent four different locations of the Galaxy (in Galactic coordinates): Galactic Center, Intermediate Latitudes, North Galactic Pole, and North Celestial Pole. The effect of the point sources is ignored, so only diffuse backgrounds are consid-ered. A requirement was set for a minimum of 10 photons per energy bin.
Figure 18.: Differential sensitivity for the projected 10-year map in 4 separate bins per decade of energy between the energies of 31.6 MeV - 1 TeV.
1.7 fermi-lat instrument response functions 35
The next set of plots displays the broadband sensitivity to power-law sources of P8R2 SOURCE V6. The broadband sensitivity is the maximum flux of a power law source, at the threshold of LAT de-tections. The plot on the left shows this for a set of spectral indices for different power law spectral models. The plot on the right shows the broadband sensitivity for the four different Galactic locations as described above (similarly to the differential sensitivity plot).
(a) Broadband sensitivity for a set of different power law spectral mod-els at the detection threshold of the LAT.
(b) Broadband sensitivity for the dif-ferent Galactic locations: Galac-tic Center, Intermediate Latitudes, North Galactic Pole, and North Ce-lestial Pole.
1.7.3 Low Energy Performance Plots
In order to determine the performance of data at low energies, we should point out the improvements made by Pass 8 over Pass 7, which struggled with performing at energies under 100 MeV. As the team indicates, it is not recommended to use data below 30 MeV, as the systematic uncertainties of the LAT’s IRF become considerably larger.
Below, we present three relevant plots, in Figures 19, 20 and 21: (1) The effective area (Log scale) versus energy (down to 10 MeV) for normal incidence photons (θ=0).
(2) The acceptance (Log scale) versus energy (down to 10 MeV). (3) The energy resolution (acceptance weighted) versus energy (down to 10 MeV), for FRONT/BACK and TOTAL event types.
Figure 20.: Effective area versus energy, down to 10 MeV.
Figure 21.: Acceptance versus energy, down to 10 MeV.
Figure 22.: Energy resolution versus energy, down to 10 MeV, for FRONT/BACK and TOTAL event types.
1.7 fermi-lat instrument response functions 37
1.7.4 Pass 8 versus Pass 7
In this short subsection, we present the comparison plots between Pass 8 (P8R2 SOURCE V6) and Pass 7 (P7REP SOURCE V15), in all the different IRF categories discussed above. The first noticeable dif-ference is the slight increase in energy range. The plots are (in order of appearance, from left to right):
(1) Acceptance (2) Effective Area
(3) PSF at 68% containment (acc. weighted)
(4) Energy dispersion (EDISP) at 68% containment (acc. weighted)
As a last comparison between Pass 8 and Pass 7, the plots that fol-low show differential (left) and broadband (right) sensitivities. This is specifically for the Galactic Pole, for a 10-year observation in survey mode.
The comparison plots show multiple improvements of the Pass 8 analysis over the Pass 7 analysis. More specifically, the acceptance is improved by 25% at low and intermediate energy ranges (100 MeV - 100 GeV) and up to 50% improvement at high energies (TeV). The effective area is also improved up to 50% at low energies, by 20% at intermediate energies and up to 100% at high energies. A sim-ilar trend in improvement is achieved in energy dispersion. These improvements result in an improvement in the sensitivity of detec-tion of the differential flux (E2× Flux). The comparison ratio in the bottom of the last two figures shows improvements of 100% for low and high energies and around 30% for intermediate energies. The resulting effects of these significant improvements will be displayed in the following chapters of the present work, where we will directly compare our results with the results presented in the most recent Fermi-LAT Collaboration release of 2013 (1), which used the Pass 7 analysis framework.
2
A N A LY S I S M E T H O D
Using data collected by the Fermi-LAT instrument and provided by the Fermi Science Team website, specifically in the Fermi SSC Data Server Library (23), we performed a data analysis of the area around Cas A. In the sections that follow, we present a detailed explanation of the analysis, from the selection of data to the processes that were followed to generate the results of this work. As our aim is to com-pare our results with previous papers published by the Fermi-LAT collaboration, we followed a similar procedure to perform our data analysis.
After the selection of data, we present our likelihood statistics and analysis procedure, our generated model and residual maps and the obtained Spectral Energy Distribution (SED) plot, in comparison to the previously generated distribution provided by the Fermi-LAT col-laboration (1).
2.1 d ata s e l e c t i o n
For the analysis, data was taken from the Fermi-LAT satellite for al-most 8 years of observations (91 months), from August 5th 2008 until March 3rd 2016. The mission elapsed time is from 239587201 until 478656004 with a radius of 20◦×20◦ degrees around the Cas A su-pernova remnant (equitorial coordinates 350.85, 58.815). The energies considered are from 100 MeV up to 300 GeV, with a data size of 2.57 GB. The tools used for the first analysis were the Fermi-LAT Science Tools provided by the Fermi website and more specifically the gtlike tool was used to process the data, whereas the user contributed pack-age called FermiPy was used for a second analysis. The same data set was used for both analyses. The third Fermi catalogue (3FGL) was taken, combined with the P8SOURCE class event selection and P8 V6 IRFs. The radius of interest (ROI) was taken to be 10◦.
The parameters for sources above ROI = 8◦ were left fixed (frozen), whereas the parameters for sources below this ROI were left free.
For this analysis we followed the tutorial provided in the Fermi website using Binned Analysis. Our first step is filtering the data, by using the tool ”gtselect”, in the suggested energy range of 100 MeV
100 GeV and a radius around Cas A of 100. We also exclude events
with a zenith angle greater than 900. These filtered data were used to create a 2-dimensional counts map (CMAP) of 100×100 pixels with a bin size of 0.1 (10 degrees radius / 0.1 = 100 pixels size) in an Aitoff projection. A counts map is a 2-dimensional projection of the sum of photons detected in each position which helps us distinguish po-tential sources of emission, with stronger candidates having a higher count of photons than their surroundings. The generated counts map is presented in Figure 24.
Figure 25.: Counts Map created using energies of 100 MeV - 100 GeV in a 10-degree radius around Cas A. The numbers on each axis are the number of pixels of the image (100×100). The colorbar on the right side shows the photon count for each pixel.
Subsequently, we generate an exposure map, divided in 20 energy bins, as well as a cube map (ccube), divided in 30 energy bins. We use these two to create a binned exposure map, again in 30 energy bins. Using the ”gtlike” tool we run a script using the 3FGL LAT cat-alog file to create the source file, initially with free parameters for all sources. After manually editing the initial source file to ”freeze” the parameters of sources outside a radius of 8 degrees around Cas A, a method suggested by the Fermi Science Tools team on their website, so that we can produce a converging log-likelihood plot. This likeli-hood method is explained in detail in the following sections and the aforementioned plot is displayed in Figure 25.
2.2 likelihood statistics 41
2.2 l i k e l i h o o d s tat i s t i c s
At this point, before proceeding any further, it is suitable to explain in detail the process behind the log-likelihood estimation as well as what this represents.
The likelihood (L) represents the probability that indicates the sig-nificance of the photon counts in our maps given the input gamma-ray data set model. This likelihood analysis (or log-likelihood in this case) aims to provide us with data that can answer key questions like the significance of a detection, as well as an accurate measurement of the flux of the detected sources, an estimation of the spectral index, the sources’ locations in the sky and the errors of these values. To do this, a demonstration of some key parts of the statistics around the likelihood estimation is required.
For a set of data X = xi with model parametersΘ= θi, the
proba-bility rule for independent events is: P(X|Θ) =∏iP(xi|Θ) = L(Θ|X),
where L represents the likelihood. To estimate the parameters, it is better to maximize the likelihood and for this it is more suitable to work with the log-likelihood:
lnL(Θ) =lnL(Θ|X) =
∑
i
lnP(xi|Θ) (1)
For each parameter θ, we have a probability:
P(xi|Θ) = 1 √ 2πσi e− (xi−θ)2 2σ2i (2) Thus, the variable part of the log-likelihood is equal to:
lnL(θ) = −
∑
(xi−θ)2
2σ2 i
(3) The Gaussian error of the variable θ is given by:
σθ =
1 q
∑(1/σi2)
(4)
The Test Statistic for each measurement of a detection is given by: TS=2(lnLd−lnLn) ≈χ2(d) (5)
where Lddenotes the likelihood of an actual detection and Lndenotes
the likelihood of a null detection. The significance is: σ=√TS. With these tools in our disposal, we can now proceed to explain in detail the rest of our analysis process, as well as discuss its results and their significance.
2.3 l i k e l i h o o d a na ly s i s
The next step was to generate a source map, using a list of all the sources around Cas A from the latest 3FGL catalog as well as the diffuse background emission provided by the Fermi-LAT collabora-tion with the files: ”iso P8R2 SOURCE V6 v06” for the Instrument Response Functions (IRFs) and ”gll iem v06” for the diffuse back-ground emission. A list of 100 separate sources was generated, one of which was Cas A (3FGL catalog name: ”3FGL J2323.4+5849”).
Using this source map and the previous exposure map, we per-formed a binned analysis with the ”P8R2 SOURCE V6” IRFs to cre-ate a binned exposure map. With the tool "pyLikelihood" and the binned exposure map at our disposal, we perform a likelihood anal-ysis. As instructed by the Fermi Science website, we used a tolerance of 10−8. Using the source file that has the parameters of sources out-side of a radius of 8 degrees from Cas A as ”frozen”, the likelihood calculation converged. Using the previous source file, we received a non-convergence and we were prompted to freeze the aforemen-tioned sources’ parameters. A list of all the sources included in the analysis (pulled from the source file) with their corresponding test statistic (TS) values was also printed and investigated, to see which of the sources had a more dominant effect.
The likelihood fit plot is displayed in the Figure 25, for all the sources in the source file. The figure in colour makes it easier to distinguish between the dominant and non-dominant sources in our model, although this had to be done in combination with having to manually go through source file to identify the source names that are more dominant and compare them to the plot, because the number of sources is too high to be distinguishable solely by colour sepa-ration. The second plot is the same as the first one, in black and white, and only included to indicate the proper axis labels of the first (coloured) plot: Counts/Bin v Energy. The last plot showcases the errors the 30 energy bins that were used and each bin is compared to the best fit line. The errors become quite high at higher energies and so does the deviation of the bins from this line at higher ener-gies. Apart from the diffuse background emission, we identified 12 other dominant sources with a TS > 20 in the likelihood plot, other than Cas A. 7 of these had a TS > 100, the highest of which was ”3FGL J0002.6+6218”, with a TS = 748.37. For reference, Cas A was
2.3 likelihood analysis 43
Figure 26.: Likelihood fit plot of all the sources in the source file. The top figure is the plot in colour, with the colours labeled for each source on the right hand side. The bottom figure demonstrates the errors for each of the 30 Energy bins, in comparison with the best fit line.
2.4 m o d e l a n d r e s i d ua l s m a p s
The next step is creating the model map. Having obtained the log-likelihood data for Cas A, we make use of the ”like.logLike.writeXml” command to create a new source file, which we use as our source model. Combining this source model with the combined source maps we generated earlier for all the sources around Cas A, we are finally able to create the model map. This model map is presented in Figure 26.
Figure 27.: Model Map for the area around Cas A. This is a heat map, with Cas A in the center of the figure and the Perseus arm emission noticeable in the lower-right part. The num-bers on each axis are the number of pixels of the image (100x100). The colorbar on the right side shows the pho-ton count for each pixel.
By subtracting our generated model map from the counts map and then dividing this subtraction over our model map, we can create what is called a ”Residual Map”. Residual maps are a very useful way to gauge how representative is our model map, made by using all the sources and the diffuse background emission, compared to the raw counts map. It is a tool to compare the two and see if there are any under-corrections or over-corrections, by visualizing this comparison in a heat map. In the case of under-correction, we will notice spots on the residuals map that are brighter than their surroundings. The opposite is true for over-correction: we will notice darker spots dom-inating in some areas, and those are the areas we over-estimated the
2.4 model and residuals maps 45
emission, thus over-correcting it. This means that the ”smoother” or, in better words, the more uniform our map is, the more suited is our model map to the observations. The residual map created using the ”farith” tool from the Fermi Tools website is displayed below, in Figure 27.
Something that is barely noticeable with the first residual map is that there are areas of over-correction (1st image). This effect is more pronounced when we smooth the image using ds9 (2nd image). Specifically, there are two areas of over-correction: one at l = 110.15 and b = 0.24 and one at l = 111.46 and b = −3.81. The first one seems to be more significant than the second. Although we will an-alyze this result in a later section of the ”Discussion” chapter, we can point out that a very similar over-correction effect was observed by and discussed in the latest paper from the Fermi-LAT collaboration (1) at exactly the same area close to Cas A and it is related to an over-estimation of the material present in the region by the current NASA diffuse background maps.
(a) Residual Map generated in our analysis, by subtracting the model map from the counts map and then dividing the result over the model map. The over-corrections are not clearly visible in this un-smoothed map.
(b) Same Residual Map, smoothed, with a different heat map. The effect of over-correction is now clearly visible at two regions around Cas A.
By studying the output Test Statistic (TS) values for each of the sources in our output source file, we can directly gauge the domi-nance of each individual source in the ROI. Comparing the TS value of Cas A to the highest TS values of our residual sources provides a clue whether or not our results are contaminated. The Fermi website indicates that TS values that are< 20 should be considered insignifi-cant enough, so that we can accept our results as not contaminated.
As stated at the end of the likelihood analysis section, the TS value for Cas A is: 3293.7841696178075. We also obtained the TS values for
the surrounding sources. The function ”like.Ts(source)” helps us print these values, so that we can go over them manually.
2.5 s p e c t r a l e n e r g y d i s t r i b u t i o n
We progressed using a highly suggested (by the Fermi website) user-contributed set of tools, called ”fermipy”. Fermipy is an external tool that uses simplified commands to perform some of the commands that the Fermi Science Tools package itself can perform. It is way more streamlined and straight-forward to use, but the support is not as up-to-date as the Fermi Science Tools package, naturally, as it is user-contributed.
Using the GTAnalysis tool from fermipy, we can run the analysis to get the Spectral Energy Distribution tables. For this, we first edit our own configuration file with our preferences. More specifically, we edit it to take an ROI of 15 degrees around Cas A, with energies rang-ing between 100 MeV - 100 GeV and we set our event types and time range, as well as any other indicator that is required, to be exactly the same as in our previous Fermi Tools analysis. It is also noted that we use exactly the same background diffuse models as well.
The straight-forward gta.sed commands let us retrieve and print out tables for the flux, differential flux and E2×differential flux in each
bin, energy bin centers, the 1σ errors and 2σ upper limits for each of the fluxes, as well as the TS values for each bin. These tables are provided below. It should be noted that we use the 2σ upper lim-its as upper limlim-its for the values of the E2×Differential Flux which are many orders-of-magnitude lower than their corresponding upper limits. Thus, the first 6 values as well as the last value of the 2σ upper limits table are used instead of their corresponding E2×Differential Flux values, as seen in Figure 28 which is a Spectral Energy Distribu-tion (SED) comparison plot between our results and the Fermi-LAT Collaboration 2013 paper data (1). Additionally, the values for en-ergy bins < 100 MeV are not used, as they are unreliable, based on the statistics of the instruments discussed in section 1.7.
2.5 spectral energy distribution 47
Flux in each bin:
2.40604217e−08 2.71572724e−09 5.93094713e−22 2.37082499e−16 2.14867418e−15 1.40644921e−21 1.03709659e−09 1.37236685e−09 1.37413991e−09 1.08383677e−09 8.69581802e−10 7.66576086e−10 6.60988702e−10 4.40890446e−10 2.85268553e−10 1.56106876e−10 1.38515997e−10 1.35085368e−10 6.92401525e−11 3.64631667e−11 4.32562443e−11 2.73975930e−11 1.50672102e−11 1.11069364e−11
Energy bin centers:
115.47819847 153.99265261 205.35250265 273.84196343 365.17412725 486.96752517 649.38163158 865.96432336 1154.78198469 1539.92652606 2053.52502646 2738.41963426 3651.74127255 4869.67525166 6493.81631576 8659.6432336 11547.81984689 15399.26526059 20535.25026457 27384.19634264 36517.41272548 48696.75251659 64938.16315762 86596.43233601
E2×Differential Flux in each bin:
9.62009104e−06 1.44797844e−06 4.21696503e−19 2.24789201e−13 2.71673021e−12 2.37137371e−18 2.33182071e−06 4.11477569e−06 5.49423082e−06 5.77882879e−06 6.18281696e−06 7.26827272e−06 8.35737680e−06 7.43372747e−06 6.41401315e−06 4.68056177e−06 5.53829239e−06 7.20251641e−06 4.92304678e−06 3.45724638e−06 5.46921198e−06 4.61942964e−06 3.38773003e−06 3.33019934e−06
1σ error of E2×Differential Flux in each bin:
9.01800091e−07 8.30044616e−07 9.57931778e−08 9.06682516e−08 2.56553738e−07 1.99424521e−07 5.57273083e−07 5.39180911e−07 5.35652592e−07 5.25997834e−07 5.56140310e−07 6.17043968e−07 6.89796081e−07 7.38056673e−07 7.65367712e−07 7.30126341e−07 9.14827808e−07 1.13563295e−06 1.09967904e−06 1.08669169e−06 1.55495681e−06 1.60972932e−06 1.55028770e−06 1.90774999e−06
2σ upper limits of E2×Differential Flux in each bin:
1.11126606e−05 2.81976260e−06 2.53101265e−07 2.39489535e−07 5.86602868e−07 4.71566259e−07 3.26397216e−06 5.01939266e−06 6.39649685e−06 6.66889806e−06 7.12685249e−06 8.31939053e−06 9.53347925e−06 8.70464718e−06 7.74188343e−06 5.97407246e−06 7.16599314e−06 9.22068337e−06 6.94413314e−06 5.52238555e−06 8.39395573e−06 7.74673472e−06 6.56154128e−06 7.39134964e−06
Figure 28 is a direct comparison of the SED between our analysis (labeled as CasA Analysis) vs the Fermi-LAT Collaboration 2013 paper analysis. On the top graph our only these two are presented, whereas on the bottom graph the MAGIC + VERITAS SED data points are also included, to provide a perspective on higher energies.
Figure 29.: SED (E2× Differential Flux in each energy bin) for our analysis (blue data points) and for the Fermi 2013 paper (red data points). Our errors our computed with 1σ standard deviation.
What is immediately obvious is the much lower errors of our analy-sis compared to the 2013 Fermi-LAT Collaboration’s paper one, which is mainly due to two reasons: (a) The higher amount of data over a longer period of time, as we have used almost 8 years of observa-tions in our analysis compared to the≈5 years of observations used in the 2013 Fermi-LAT Collaboration paper analysis. (b) Our anal-ysis framework was Pass8, which offers a significant improvement over the previously used Pass7 analysis framework in the 2013 paper by the Fermi-LAT collaboration. The latter was discussed in length in section 1.7 Fermi-LAT instrument response functions. Another thing that is immediately obvious are the much lower upper-limits of our analysis on the lower energies, in the range between 100 MeV - 400 MeV. These results will be discussed and analyzed further in the final chapter of the present work.
3
D I S C U S S I O N
3.1 a na ly s i s r e s u lt s 3.1.1 Errors
As mentioned earlier, the striking result from the analysis results pre-sented in the previous section, is the discrepancy in SED errors be-tween our analysis and the Fermi 2013 paper presented ones. As the errors are estimated to be roughly 1 order of magnitude (or 10 times) smaller than those of the 2013 paper by Fermi, we need to perform some investigation into what causes this vast improvement. As mentioned in an earlier section, we suspect two culprits to be be-hind this: The observation time difference and the major statistical improvements caused by the implementation of Pass 8 over Pass 7.
Beginning with the difference in observation time, our observations are taken from August 5th 2008 until March 3rd 2016. This is a period of ≈91 months, compared to the observation period of≈ 44 months that is taken in the 2013 Fermi-LAT Collaboration paper. To consider the potential effect this has on the flux errors, we can perform some ”back of the envelope” calculations. Let’s consider FBerr as our flux
er-rors and FAerr as the flux errors of the 2013 Fermi paper. We know that the change in errors scales with the square root of the ratio between the two corresponding observation times: tB as our observation time
and tA as the observation time of the 2013 Fermi paper:
FBerr ≈α× (r tA
tB
) ×FAerr (6)
where α is a constant.
Therefore we can extrapolate from this that the improvement in errors we would expect, purely from an observation time standpoint, is approximately: q tA tB = q 44 91 ≈ 0.695. Thus, a 70% improvement
can be assumed to come solely from the larger observation time. The other, and perhaps most important, cause of this dramatic im-provement in errors, is the transition from Pass 7 to Pass 8. Looking back at the Fermi-LAT Instrument Response Functions section of this work and more precisely at the Pass 8 versus Pass 7 subsection, we can investigate further into the matter. As shown in the aforementioned
section, we see an improvement in Acceptance. This improvement is around 25% both for energies around a few GeV and for energies around a few hundreds of MeV. It is also worth pointing out that there is an improvement in flux sensitivities. Focusing on the ratio be-tween Pass 8 and Pass 7, we can see that, for broadband sensitivities, there is an improvement in flux sensitivity ranging from 50% around a few hundred MeV down to ≈ 30% around a few GeV. Both these comparison figures are presented again below, in Figures 29 and 30, to ease the reading process.
It is very important to note that these comparison plots are between the latest version of Pass 7 and the version of Pass 8 we used. How-ever, the 2013 paper of the Fermi-LAT collaboration (1) didn’t use the latest version of Pass 7 (which is version 15), but an earlier Pass 7 version, Pass7v6. Although there are no figures provided by the Fermi-LAT collaboration to compare Pass8v2 to Pass7v6, the Pass7v6 figures display even lower acceptance and broadband & differential sensitivities. This means that the improvement due to the version of Pass 8 used for this work is even more pronounced and significant.
Figure 30.: Comparison between the acceptance of the latest version of Pass 7 to the version of Pass 8 used in the present work. Improvements of≈25% are noticeable at energies ranging from of a few hundred MeV up to a few GeV.
