PSFMap#
- class gammapy.irf.PSFMap(psf_map, exposure_map=None)[source]#
Bases:
gammapy.irf.core.IRFMap
Class containing the Map of PSFs and allowing to interact with it.
- Parameters
Examples
from astropy.coordinates import SkyCoord from gammapy.maps import WcsGeom, MapAxis from gammapy.data import Observation from gammapy.irf import load_cta_irfs from gammapy.makers import MapDatasetMaker # Define observation pointing = SkyCoord("0d", "0d") filename = "$GAMMAPY_DATA/cta-1dc/caldb/data/cta/1dc/bcf/South_z20_50h/irf_file.fits" irfs = load_cta_irfs(filename) obs = Observation.create(pointing=pointing, irfs=irfs, livetime="1h") # Create WcsGeom # Define energy axis. Note that the name is fixed. energy_axis = MapAxis.from_energy_bounds("0.1 TeV", "10 TeV", nbin=3, name="energy_true") # Define rad axis. Again note the axis name rad_axis = MapAxis.from_bounds(0, 0.5, nbin=100, name="rad", unit="deg") geom = WcsGeom.create( binsz=0.25, width="5 deg", skydir=pointing, axes=[rad_axis, energy_axis] ) maker = MapDatasetMaker() psf = maker.make_psf(geom=geom, observation=obs) # Get a PSF kernel at the center of the image geom=exposure_geom.upsample(factor=10).drop("rad") psf_kernel = psf_map.get_psf_kernel(geom=geom)
Attributes Summary
Mask safe for the map
Methods Summary
containment
(rad, energy_true[, position])Containment at given coords
containment_radius
(fraction, energy_true[, ...])Containment at given coords
containment_radius_map
(energy_true[, fraction])Containment radius map.
copy
()Copy IRF map
cutout
(position, width[, mode])Cutout IRF map.
downsample
(factor[, axis_name, weights])Downsample the spatial dimension by a given factor.
from_gauss
(energy_axis_true[, rad_axis, ...])Create all -sky PSF map from Gaussian width.
from_geom
(geom)Create psf map from geom.
from_hdulist
(hdulist[, hdu, hdu_bands, ...])Create from
HDUList
.get_psf_kernel
(geom[, position, max_radius, ...])Returns a PSF kernel at the given position.
Normalize PSF map
peek
([figsize])Quick-look summary plots.
plot_containment_radius_vs_energy
([ax, fraction])Plot containment fraction as a function of energy.
plot_psf_vs_rad
([ax, energy_true])Plot PSF vs radius.
read
(filename[, format, hdu])Read an IRF_map from file and create corresponding object"
sample_coord
(map_coord[, random_state])Apply PSF corrections on the coordinates of a set of simulated events.
slice_by_idx
(slices)Slice sub dataset.
stack
(other[, weights, nan_to_num])Stack IRF map with another one in place.
to_hdulist
([format])Convert to
HDUList
.to_image
([spectrum, keepdims])Reduce to a 2-D map after weighing with the associated exposure and a spectrum
to_region_nd_map
(region)Extract IRFMap in a given region or position
write
(filename[, overwrite, format])Write IRF map to fits
Attributes Documentation
- energy_name#
- mask_safe_image#
Mask safe for the map
- psf_map#
- required_axes = ['rad', 'energy_true']#
- tag = 'psf_map'#
Methods Documentation
- copy()#
Copy IRF map
- cutout(position, width, mode='trim')#
Cutout IRF map.
- Parameters
- Returns
- cutout
IRFMap
Cutout IRF map.
- cutout
- downsample(factor, axis_name=None, weights=None)#
Downsample the spatial dimension by a given factor.
- classmethod from_gauss(energy_axis_true, rad_axis=None, sigma=<Quantity 0.1 deg>, geom=None)[source]#
Create all -sky PSF map from Gaussian width.
This is used for testing and examples.
The width can be the same for all energies or be an array with one value per energy node. It does not depend on position.
- classmethod from_geom(geom)[source]#
Create psf map from geom.
- Parameters
- geom
Geom
PSF map geometry.
- geom
- Returns
- psf_map
PSFMap
Point spread function map.
- psf_map
- classmethod from_hdulist(hdulist, hdu=None, hdu_bands=None, exposure_hdu=None, exposure_hdu_bands=None, format='gadf')#
Create from
HDUList
.- Parameters
- hdulist
HDUList
HDU list.
- hdustr
Name or index of the HDU with the IRF map.
- hdu_bandsstr
Name or index of the HDU with the IRF map BANDS table.
- exposure_hdustr
Name or index of the HDU with the exposure map data.
- exposure_hdu_bandsstr
Name or index of the HDU with the exposure map BANDS table.
- format{“gadf”, “gtpsf”}
File format
- hdulist
- Returns
- irf_map
IRFMap
IRF map.
- irf_map
- get_psf_kernel(geom, position=None, max_radius=None, containment=0.999, factor=4)[source]#
Returns a PSF kernel at the given position.
The PSF is returned in the form a WcsNDMap defined by the input Geom.
- Parameters
- geom
Geom
Target geometry to use
- position
SkyCoord
Target position. Should be a single coordinate. By default the center position is used.
- max_radius
Angle
maximum angular size of the kernel map
- containmentfloat
Containment fraction to use as size of the kernel. The max. radius across all energies is used. The radius can be overwritten using the
max_radius
argument.- factorint
oversampling factor to compute the PSF
- geom
- Returns
- kernel
PSFKernel
the resulting kernel
- kernel
- plot_containment_radius_vs_energy(ax=None, fraction=(0.68, 0.95), **kwargs)[source]#
Plot containment fraction as a function of energy.
The method plots the containment radius at the center of the map.
- plot_psf_vs_rad(ax=None, energy_true=<Quantity [ 0.1, 1., 10. ] TeV>, **kwargs)[source]#
Plot PSF vs radius.
The method plots the profile at the center of the map.
- classmethod read(filename, format='gadf', hdu=None)#
Read an IRF_map from file and create corresponding object”
- Parameters
- filenamestr or
Path
File name
- format{“gadf”, “gtpsf”}
File format
- hdustr or int
HDU location
- filenamestr or
- Returns
- irf_map
PSFMap
,EDispMap
orEDispKernelMap
IRF map
- irf_map
- sample_coord(map_coord, random_state=0)[source]#
Apply PSF corrections on the coordinates of a set of simulated events.
- Parameters
- map_coord
MapCoord
object. Sequence of coordinates and energies of sampled events.
- random_state{int, ‘random-seed’, ‘global-rng’,
RandomState
} Defines random number generator initialisation. Passed to
get_random_state
.
- map_coord
- Returns
- corr_coord
MapCoord
object. Sequence of PSF-corrected coordinates of the input map_coord map.
- corr_coord
- slice_by_idx(slices)#
Slice sub dataset.
The slicing only applies to the maps that define the corresponding axes.
- stack(other, weights=None, nan_to_num=True)#
Stack IRF map with another one in place.
- to_hdulist(format='gadf')#
Convert to
HDUList
.- Parameters
- format{“gadf”, “gtpsf”}
File format
- Returns
- hdu_list
HDUList
HDU list.
- hdu_list
- to_image(spectrum=None, keepdims=True)[source]#
Reduce to a 2-D map after weighing with the associated exposure and a spectrum
- Parameters
- spectrum
SpectralModel
, optional Spectral model to compute the weights. Default is power-law with spectral index of 2.
- keepdimsbool, optional
If True, the energy axis is kept with one bin. If False, the axis is removed
- spectrum
- Returns
- to_region_nd_map(region)#
Extract IRFMap in a given region or position
If a region is given a mean IRF is computed, if a position is given the IRF is interpolated.
- write(filename, overwrite=False, format='gadf')#
Write IRF map to fits
- Parameters
- filenamestr or
Path
Filename to write to
- overwritebool
Whether to overwrite
- format{“gadf”, “gtpsf”}
File format
- filenamestr or