PSFMap

class gammapy.cube.PSFMap(psf_map, exposure_map=None)[source]

Bases: object

Class containing the Map of PSFs and allowing to interact with it.

Parameters
psf_mapMap

the input PSF Map. Should be a Map with 2 non spatial axes. rad and true energy axes should be given in this specific order.

exposure_mapMap

Associated exposure map. Needs to have a consistent map geometry.

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.cube 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')

Methods Summary

containment_radius_map(self, energy[, fraction])

Containment radius map.

copy(self)

Copy PSFMap

cutout(self, position, width[, mode])

Cutout map psf map.

from_energy_dependent_table_psf(table_psf)

Create PSF map from table PSF object.

from_geom(geom)

Create psf map from geom.

from_hdulist(hdulist[, psf_hdu, …])

Convert to HDUList.

get_energy_dependent_table_psf(self, position)

Get energy-dependent PSF at a given position.

get_psf_kernel(self, position, geom[, …])

Returns a PSF kernel at the given position.

read(filename, \*\*kwargs)

Read a psf_map from file and create a PSFMap object

sample_coord(self, map_coord[, random_state])

Apply PSF corrections on the coordinates of a set of simulated events.

stack(self, other[, weights])

Stack PSFMap with another one in place.

to_hdulist(self[, psf_hdu, psf_hdubands, …])

Convert to HDUList.

to_image(self[, spectrum, keepdims])

Reduce to a 2-D map after weighing with the associated exposure and a spectrum

write(self, filename[, overwrite])

Write to fits

Methods Documentation

containment_radius_map(self, energy, fraction=0.68)[source]

Containment radius map.

Parameters
energyQuantity

Scalar energy at which to compute the containment radius

fractionfloat

the containment fraction (range: 0 to 1)

Returns
containment_radius_mapMap

Containment radius map

copy(self)[source]

Copy PSFMap

cutout(self, position, width, mode='trim')[source]

Cutout map psf map.

Parameters
positionSkyCoord

Center position of the cutout region.

widthtuple of Angle

Angular sizes of the region in (lon, lat) in that specific order. If only one value is passed, a square region is extracted.

mode{‘trim’, ‘partial’, ‘strict’}

Mode option for Cutout2D, for details see Cutout2D.

Returns
cutoutPSFMap

Cutout psf map.

classmethod from_energy_dependent_table_psf(table_psf)[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_psfEnergyDependentTablePSF

Table PSF

Returns
psf_mapPSFMap

Point spread function map.

classmethod from_geom(geom)[source]

Create psf map from geom.

Parameters
geomGeom

PSF map geometry.

Returns
psf_mapPSFMap

Point spread function map.

classmethod from_hdulist(hdulist, psf_hdu='PSFMAP', psf_hdubands='BANDSPSF', exposure_hdu='EXPMAP', exposure_hdubands='BANDSEXP')[source]

Convert to HDUList.

Parameters
psf_hdustr

Name or index of the HDU with the psf_map data.

psf_hdubandsstr

Name or index of the HDU with the psf_map BANDS table.

exposure_hdustr

Name or index of the HDU with the exposure_map data.

exposure_hdubandsstr

Name or index of the HDU with the exposure_map BANDS table.

get_energy_dependent_table_psf(self, position)[source]

Get energy-dependent PSF at a given position.

Parameters
positionSkyCoord

the target position. Should be a single coordinates

Returns
psf_tableEnergyDependentTablePSF

the table PSF

get_psf_kernel(self, 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.

Parameters
positionSkyCoord

the target position. Should be a single coordinate

geomGeom

the target geometry to use

max_radiusAngle

maximum angular size of the kernel map

factorint

oversampling factor to compute the PSF

Returns
kernelPSFKernel

the resulting kernel

classmethod read(filename, **kwargs)[source]

Read a psf_map from file and create a PSFMap object

sample_coord(self, map_coord, random_state=0)[source]

Apply PSF corrections on the coordinates of a set of simulated events.

Parameters
map_coordMapCoord 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.

Returns
corr_coordMapCoord object.

Sequence of PSF-corrected coordinates of the input map_coord map.

stack(self, other, weights=None)[source]

Stack PSFMap with another one in place.

Parameters
otherPSFMap

the psfmap to be stacked with this one.

to_hdulist(self, psf_hdu='PSFMAP', psf_hdubands='BANDSPSF', exposure_hdu='EXPMAP', exposure_hdubands='BANDSEXP')[source]

Convert to HDUList.

Parameters
psf_hdustr

Name or index of the HDU with the psf_map data.

psf_hdubandsstr

Name or index of the HDU with the psf_map BANDS table.

exposure_hdustr

Name or index of the HDU with the exposure_map data.

exposure_hdubandsstr

Name or index of the HDU with the exposure_map BANDS table.

Returns
hdu_listHDUList
to_image(self, spectrum=None, keepdims=True)[source]

Reduce to a 2-D map after weighing with the associated exposure and a spectrum

Parameters
spectrumSpectralModel, 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

Returns
psf_outPSFMap

PSFMap with the energy axis summed over

write(self, filename, overwrite=False, **kwargs)[source]

Write to fits