Sensitivity Tools
The fermipy-flux-sensitivity
script calculates the LAT flux
threshold for a gamma-ray source in bins of energy (differential
sensitivity) and integrated over the full LAT energy range (integral
sensitivity). The source flux threshold is the flux at which the
median TS of a source (twice the likelihood ratio of the best-fit model with and
without the source) equals a certain value. Primary inputs to this script
are the livetime cube (output of gtltcube
) and the model cube for the galactic diffuse
background. The obs_time_yr
option can be used to rescale the
livetime cube to a shorter or longer observation time.
$ fermipy-flux-sensitivity --glon=30 --glat=30 --output=lat_sensitivity.fits \
--ltcube=ltcube.fits --galdiff=gll_iem_v06.fits --event_class=P8R2_SOURCE_V6 \
--ts_thresh=25.0 --min_counts=10.0
If no livetime cube is provided then the sensitivity will be computed assuming an “ideal” survey-mode operation with uniform exposure over the whole sky and no Earth obscuration or deadtime. By default the flux sensitivity will be calculated for a TS threshold of 25 and at least 3 counts.
A map of sensitivity with WCS or HEALPix pixelization can be generated
by setting the map_type
argument to either wcs
or hpx
:
# Generate a WCS sensitivity map of 50 x 50 deg centered at (glon,glat) = (30,30)
$ fermipy-flux-sensitivity --glon=30 --glat=30 --output=lat_sensitivity_map.fits \
--ltcube=ltcube.fits --galdiff=gll_iem_v06.fits --event_class=P8R2_SOURCE_V6 \
--map_type=wcs --wcs_npix=100 --wcs_cdelt=0.5 --wcs_proj=AIT
# Generate a HPX sensitivity map of nside=16
$ fermipy-flux-sensitivity --output=lat_sensitivity_map.fits \
--ltcube=ltcube.fits --galdiff=gll_iem_v06.fits --event_class=P8R2_SOURCE_V6 \
--map_type=hpx --hpx_nside=16
The integral and differential sensitivity maps will be written to the
MAP_INT_FLUX
and MAP_DIFF_FLUX
extensions respectively.
By default the flux sensitivity will be computed for a point-source
morphology. The assumed source morphology can be changed with the
spatial_model
and spatial_size
parameters:
It is possible to choose among PowerLaw, LogParabola and
PLSuperExpCutoff SED shapes using the option sedshape
.
# Generate the sensitivity to a source with a 2D gaussian morphology
# and a 68% containment radius of 1 deg located at longitude 30deg and
# latitude 30 deg and with a PLSuperExpCutoff SED with index 2.0 and
# cutoff energy 10 GeV
$ fermipy-flux-sensitivity --output=lat_sensitivity_map.fits \
--ltcube=ltcube.fits --galdiff=gll_iem_v06.fits --event_class=P8R2_SOURCE_V6 \
--spatial_model=RadialGaussian --spatial_size=1.0 --glon=30 --glat=30
--sedshape=PLSuperExpCutoff --index=2.0 --cutoff=1e4
# Generate the sensitivity map in healpix with nside 128 of a point source with
# LogParabola SED and with spectral index 2.0 and curvature index beta=0.50
# between 1 and 10 GeV
$ fermipy-flux-sensitivity --output=lat_sensitivity_map.fits \
--ltcube=ltcube.fits --galdiff=gll_iem_v06.fits --event_class=P8R2_SOURCE_V6 \
--spatial_model=PointSource --sedshape=LogParabola --index=2.0 --beta=0.50 \
--hpx_nside=128 --map_type=hpx --emin=1000 --emax=10000
The output FITS file set with the output
option contains the
following tables. Note that MAP
tables are only generated when
the map_type
argument is set.
DIFF_FLUX
: Differential flux sensitivity for a gamma-ray source at sky positition set byglon
andglat
.INT_FLUX
: Integral flux sensitivity evaluated for PowerLaw sources with spectral indices between 1.0 and 5.0 at sky positition set byglon
andglat
. Columns starting withebin
contain the source amplitude vs. energy bin.MAP_DIFF_FLUX
: Sky cube with differential flux threshold vs. sky position and energy.MAP_DIFF_NPRED
: Sky cube with counts amplitude (NPred) of a source at the detection threshold vs. position and energy.MAP_INT_FLUX
: Sky map with integral flux threshold vs. sky position. Integral sensitivity will be computed for a PowerLaw source with index equal to theindex
parameter.MAP_INT_NPRED
: Sky map with counts amplitude (NPred) of a source at the detection threshold vs. sky position.
The output file can be read using the Table
module:
from astropy.table import Table
tab = Table.read('lat_sensitivity.fits','DIFF_FLUX')
print(tab['e_min'], tab['e_max'], tab['flux'])
tab = Table.read('lat_sensitivity.fits','INT_FLUX')
print(tab['index'], tab['flux'])