PSFMap¶
-
class
gammapy.irf.
PSFMap
(psf_map, exposure_map=None)[source]¶ Bases:
gammapy.irf.irf_map.IRFMap
Class containing the Map of PSFs and allowing to interact with it.
- Parameters
Examples
import numpy as np from astropy import units as u from astropy.coordinates import SkyCoord from gammapy.maps import Map, WcsGeom, MapAxis from gammapy.irf import EnergyDependentMultiGaussPSF, EffectiveAreaTable2D from gammapy.makers.utils import make_psf_map, PSFMap, make_map_exposure_true_energy # Define energy axis. Note that the name is fixed. energy_axis = MapAxis.from_edges(np.logspace(-1., 1., 4), unit='TeV', name='energy') # Define rad axis. Again note the axis name rads = np.linspace(0., 0.5, 100) * u.deg rad_axis = MapAxis.from_edges(rads, unit='deg', name='theta') # Define parameters pointing = SkyCoord(0., 0., unit='deg') max_offset = 4 * u.deg # Create WcsGeom geom = WcsGeom.create(binsz=0.25*u.deg, width=10*u.deg, skydir=pointing, axes=[rad_axis, energy_axis]) # Extract EnergyDependentTablePSF from CTA 1DC IRF filename = '$GAMMAPY_DATA/cta-1dc/caldb/data/cta/1dc/bcf/South_z20_50h/irf_file.fits' psf = EnergyDependentMultiGaussPSF.read(filename, hdu='POINT SPREAD FUNCTION') psf3d = psf.to_psf3d(rads) aeff2d = EffectiveAreaTable2D.read(filename, hdu='EFFECTIVE AREA') # Create the exposure map exposure_geom = geom.to_image().to_cube([energy_axis]) exposure_map = make_map_exposure_true_energy(pointing, "1 h", aeff2d, exposure_geom) # create the PSFMap for the specified pointing psf_map = make_psf_map(psf3d, pointing, geom, max_offset, exposure_map) # Get an EnergyDependentTablePSF at any position in the image psf_table = psf_map.get_energy_dependent_table_psf(SkyCoord(2., 2.5, unit='deg')) # Write map to disk psf_map.write('psf_map.fits')
Attributes Summary
Methods Summary
containment_radius_map
(energy[, fraction])Containment radius map.
copy
()Copy IRF map
cutout
(position, width[, mode])Cutout IRF map.
downsample
(factor[, axis, weights])Downsample the spatial dimension by a given factor.
from_energy_dependent_table_psf
(table_psf[, …])Create PSF map from table PSF object.
from_gauss
(energy_axis_true[, rad_axis, sigma])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_energy_dependent_table_psf
([position])Get energy-dependent PSF at a given position.
get_psf_kernel
(position, geom[, max_radius, …])Returns a PSF kernel at the given position.
read
(filename[, 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])Stack IRF map with another one in place.
Convert to
HDUList
.to_image
([spectrum, keepdims])Reduce to a 2-D map after weighing with the associated exposure and a spectrum
write
(filename[, overwrite])Write IRF map to fits
Attributes Documentation
-
psf_map
¶
-
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=None, weights=None)¶ Downsample the spatial dimension by a given factor.
- Parameters
- factorint
Downsampling factor.
- axisstr
Which axis to downsample. By default spatial axes are downsampled.
- weights
Map
Map with weights downsampling.
- Returns
- map
IRFMap
Downsampled irf map.
- map
-
classmethod
from_energy_dependent_table_psf
(table_psf, geom=None)[source]¶ Create PSF map from table PSF object.
Helper function to create an allsky PSF map from table PSF, which does not depend on position.
- Parameters
- table_psf
EnergyDependentTablePSF
Table PSF
- table_psf
- Returns
- psf_map
PSFMap
Point spread function map.
- psf_map
-
classmethod
from_gauss
(energy_axis_true, rad_axis=None, sigma=<Quantity 0.1 deg>)[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)¶ 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.
- hdulist
- Returns
- irf_map
IRFMap
IRF map.
- irf_map
-
get_energy_dependent_table_psf
(position=None)[source]¶ Get energy-dependent PSF at a given position.
By default the PSF at the center of the map is returned.
- Parameters
- position
SkyCoord
the target position. Should be a single coordinates
- position
- Returns
- psf_table
EnergyDependentTablePSF
the table PSF
- psf_table
-
get_psf_kernel
(position, geom, max_radius=None, 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.
-
classmethod
read
(filename, hdu=None)¶ Read an IRF_map from file and create corresponding object
-
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.
- Parameters
- slicesdict
Dict of axes names and integers or
slice
object pairs. Contains one element for each non-spatial dimension. For integer indexing the corresponding axes is dropped from the map. Axes not specified in the dict are kept unchanged.
- Returns
- map_out
IRFMap
Sliced irf map object.
- map_out
-
stack
(other, weights=None)¶ Stack IRF map with another one in place.
- Parameters
- other
IRFMap
IRF map to be stacked with this one.
- weights
Map
Map with stacking weights.
- other
-
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
-
write
(filename, overwrite=False, **kwargs)¶ Write IRF map to fits