It was found in Example 3 that a non-zero value of the GCE template is preferred by a fit in the galactic center. In this example we will test the point source interpretation of this excess by including, in addition to the Poissonian templates considered, non-Poissonian point source templates of various morphologies.
Here we simply perform the run, a detailed analysis of the results can be found in the next example. A Python script version of this notebook can be found as Example8_Galactic_Center_nonPoissonian.py, which can be run faster on multiple processors with MPI (see example in Example8_Galactic_Center_Batch.batch.
NB: Even with nlive=100, this notebook can take up to an hour to complete. This highlights that for realistic non-Poissonian runs, running on multiple cores becomes necessary. We show an explicit application of this in Example 10.
NB: This example makes use of the Fermi Data, which needs to already be installed. See Example 1 for details.
In [1]:
# Import relevant modules
%matplotlib inline
%load_ext autoreload
%autoreload 2
import numpy as np
from NPTFit import nptfit # module for performing scan
from NPTFit import create_mask as cm # module for creating the mask
from NPTFit import dnds_analysis # module for analysing the output
from NPTFit import psf_correction as pc # module for determining the PSF correction
We first need to
NPTF from npfit.pyThese are done identically to Example 3, and we refer to that notebook for details.
In [2]:
n = nptfit.NPTF(tag='GCE_Example')
In [3]:
fermi_data = np.load('fermi_data/fermidata_counts.npy').astype(np.int32)
fermi_exposure = np.load('fermi_data/fermidata_exposure.npy')
n.load_data(fermi_data, fermi_exposure)
In [4]:
pscmask=np.array(np.load('fermi_data/fermidata_pscmask.npy'), dtype=bool)
analysis_mask = cm.make_mask_total(band_mask = True, band_mask_range = 2,
mask_ring = True, inner = 0, outer = 30,
custom_mask = pscmask)
n.load_mask(analysis_mask)
In [5]:
dif = np.load('fermi_data/template_dif.npy')
iso = np.load('fermi_data/template_iso.npy')
bub = np.load('fermi_data/template_bub.npy')
gce = np.load('fermi_data/template_gce.npy')
dsk = np.load('fermi_data/template_dsk.npy')
n.add_template(dif, 'dif')
n.add_template(iso, 'iso')
n.add_template(bub, 'bub')
n.add_template(gce, 'gce')
n.add_template(dsk, 'dsk')
# Remove the exposure correction for PS templates
rescale = fermi_exposure/np.mean(fermi_exposure)
n.add_template(gce/rescale, 'gce_np', units='PS')
n.add_template(dsk/rescale, 'dsk_np', units='PS')
In [6]:
n.add_poiss_model('dif', '$A_\mathrm{dif}$', fixed=True, fixed_norm=12.85)
n.add_poiss_model('iso', '$A_\mathrm{iso}$', [0,2], False)
n.add_poiss_model('gce', '$A_\mathrm{gce}$', [0,2], False)
n.add_poiss_model('bub', '$A_\mathrm{bub}$', [0,2], False)
This time we add a non-Poissonian template correlated with the Galactic Center Excess and also one spatially distributed as a thin disk. The latter is designed to account for the unresolved point sources attributed to the disk of the Milky Way (known sources in the 3FGL are masked).
In [7]:
n.add_non_poiss_model('gce_np',
['$A_\mathrm{gce}^\mathrm{ps}$','$n_1^\mathrm{gce}$','$n_2^\mathrm{gce}$','$S_b^{(1), \mathrm{gce}}$'],
[[-6,1],[2.05,30],[-2,1.95],[0.05,40]],
[True,False,False,False])
n.add_non_poiss_model('dsk_np',
['$A_\mathrm{dsk}^\mathrm{ps}$','$n_1^\mathrm{dsk}$','$n_2^\mathrm{dsk}$','$S_b^{(1), \mathrm{dsk}}$'],
[[-6,1],[2.05,30],[-2,1.95],[0.05,40]],
[True,False,False,False])
In [8]:
pc_inst = pc.PSFCorrection(psf_sigma_deg=0.1812)
f_ary, df_rho_div_f_ary = pc_inst.f_ary, pc_inst.df_rho_div_f_ary
In [9]:
n.configure_for_scan(f_ary, df_rho_div_f_ary, nexp=1)
As noted above, we take a small value of nlive simply to ensure the run finishes in a reasonable time on a single core.
In [10]:
n.perform_scan(nlive=100)
This can take up to an hour to run. The output of this run will be analyzed in detail in the next example.
In [11]:
from IPython.display import Image
Image(url = "https://imgs.xkcd.com/comics/compiling.png")
Out[11]: