EDispMap#

class gammapy.irf.EDispMap(edisp_map, exposure_map=None)[source]#

Bases: gammapy.irf.core.IRFMap

Energy dispersion map.

Parameters
edisp_mapMap

The input Energy Dispersion Map. Should be a Map with 2 non-spatial axes. migra and true energy axes should be given in this specific order.

exposure_mapMap, optional

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

Examples

# Energy dispersion map for CTA data
import numpy as np
from astropy import units as u
from astropy.coordinates import SkyCoord
from gammapy.maps import WcsGeom, MapAxis
from gammapy.irf import EnergyDispersion2D, EffectiveAreaTable2D
from gammapy.makers.utils import make_edisp_map, make_map_exposure_true_energy

# Define energy dispersion map geometry
energy_axis_true = MapAxis.from_edges(np.logspace(-1, 1, 10), unit="TeV", name="energy_true")
migra_axis = MapAxis.from_edges(np.linspace(0, 3, 100), name="migra")
pointing = SkyCoord(0, 0, unit="deg")
geom = WcsGeom.create(
        binsz=0.25 * u.deg,
        width=10 * u.deg,
        skydir=pointing,
        axes=[migra_axis, energy_axis_true],
)

# Extract EnergyDispersion2D from CTA 1DC IRF
filename = "$GAMMAPY_DATA/cta-1dc/caldb/data/cta/1dc/bcf/South_z20_50h/irf_file.fits"
edisp2D = EnergyDispersion2D.read(filename, hdu="ENERGY DISPERSION")
aeff2d = EffectiveAreaTable2D.read(filename, hdu="EFFECTIVE AREA")

# Create the exposure map
exposure_geom = geom.squash(axis_name="migra")
exposure_map = make_map_exposure_true_energy(pointing, "1 h", aeff2d, exposure_geom)

# Create the EDispMap for the specified pointing
edisp_map = make_edisp_map(edisp2D, pointing, geom, exposure_map)

# Get an Energy Dispersion (1D) at any position in the image
pos = SkyCoord(2.0, 2.5, unit="deg")
energy_axis = MapAxis.from_energy_bounds(0.1, 10, 5, unit="TeV", name="energy")
edisp = edisp_map.get_edisp_kernel(energy_axis, position=pos)

# Write map to disk
edisp_map.write("edisp_map.fits")

Attributes Summary

edisp_map

mask_safe_image

Mask safe for the map.

required_axes

tag

Methods Summary

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_diagonal_response(energy_axis_true[, ...])

Create an all-sky EDisp map with diagonal response.

from_geom(geom)

Create energy dispersion map from geometry.

from_hdulist(hdulist[, hdu, hdu_bands, ...])

Create from HDUList.

get_edisp_kernel(energy_axis[, position])

Get energy dispersion at a given position.

normalize()

Normalize PSF map.

peek([figsize])

Quick-look summary plots.

read(filename[, format, hdu, checksum])

Read an IRF_map from file and create corresponding object.

sample_coord(map_coord[, random_state, ...])

Apply the energy dispersion 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_edisp_kernel_map(energy_axis)

Convert to map with energy dispersion kernels.

to_hdulist([format])

Convert to HDUList.

to_region_nd_map(region)

Extract IRFMap in a given region or position.

write(filename[, overwrite, format, checksum])

Write IRF map to fits.

Attributes Documentation

edisp_map#
mask_safe_image#

Mask safe for the map.

required_axes = ['migra', 'energy_true']#
tag = 'edisp_map'#

Methods Documentation

copy()#

Copy IRF map.

cutout(position, width, mode='trim')#

Cutout IRF 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’}, optional

Mode option for Cutout2D, for details see Cutout2D. Default is “trim”.

Returns
cutoutIRFMap

Cutout IRF map.

downsample(factor, axis_name=None, weights=None)#

Downsample the spatial dimension by a given factor.

Parameters
factorint

Downsampling factor.

axis_namestr

Axis to downsample. By default, spatial axes are downsampled.

weightsMap, optional

Map with weights downsampling. Default is None.

Returns
mapIRFMap

Downsampled IRF map.

classmethod from_diagonal_response(energy_axis_true, migra_axis=None)[source]#

Create an all-sky EDisp map with diagonal response.

Parameters
energy_axis_trueMapAxis

True energy axis.

migra_axisMapAxis, optional

Migra axis. Default is None.

Returns
edisp_mapEDispMap

Energy dispersion map.

classmethod from_geom(geom)[source]#

Create energy dispersion map from geometry.

By default, a diagonal energy dispersion matrix is created.

Parameters
geomGeom

Energy dispersion map geometry.

Returns
edisp_mapEDispMap

Energy dispersion map.

classmethod from_hdulist(hdulist, hdu=None, hdu_bands=None, exposure_hdu=None, exposure_hdu_bands=None, format='gadf')#

Create from HDUList.

Parameters
hdulistHDUList

HDU list.

hdustr, optional

Name or index of the HDU with the IRF map. Default is None.

hdu_bandsstr, optional

Name or index of the HDU with the IRF map BANDS table. Default is None.

exposure_hdustr, optional

Name or index of the HDU with the exposure map data. Default is None.

exposure_hdu_bandsstr, optional

Name or index of the HDU with the exposure map BANDS table. Default is None.

format{“gadf”, “gtpsf”}, optional

File format. Default is “gadf”.

Returns
irf_mapIRFMap

IRF map.

get_edisp_kernel(energy_axis, position=None)[source]#

Get energy dispersion at a given position.

Parameters
energy_axisMapAxis

Reconstructed energy axis.

positionSkyCoord

The target position. Should be a single coordinates.

Returns
edispEnergyDispersion

The energy dispersion (i.e. rmf object).

normalize()[source]#

Normalize PSF map.

peek(figsize=(15, 5))[source]#

Quick-look summary plots.

Plots corresponding to the center of the map.

Parameters
figsizetuple

Size of figure.

classmethod read(filename, format='gadf', hdu=None, checksum=False)#

Read an IRF_map from file and create corresponding object.

Parameters
filenamestr or Path

File name.

format{“gadf”, “gtpsf”}, optional

File format. Default is “gadf”.

hdustr or int

HDU location. Default is None.

checksumbool

If True checks both DATASUM and CHECKSUM cards in the file headers. Default is False.

Returns
irf_mapPSFMap, EDispMap or EDispKernelMap

IRF map.

sample_coord(map_coord, random_state=0, chunk_size=10000)[source]#

Apply the energy dispersion corrections on the coordinates of a set of simulated events.

Parameters
map_coordMapCoord

Sequence of coordinates and energies of sampled events.

random_state{int, ‘random-seed’, ‘global-rng’, RandomState}, optional

Defines random number generator initialisation. Passed to get_random_state. Default is 0.

chunk_sizeint

If set, this will slice the input MapCoord into smaller chunks of chunk_size elements. Default is 10000.

Returns
MapCoord.

Sequence of energy dispersion corrected coordinates of the input map_coord map.

slice_by_idx(slices)#

Slice sub dataset.

The slicing only applies to the maps that define the corresponding axes.

Parameters
slicesdict

Dictionary 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 dictionary are kept unchanged.

Returns
map_outIRFMap

Sliced IRF map object.

stack(other, weights=None, nan_to_num=True)#

Stack IRF map with another one in place.

Parameters
otherIRFMap

IRF map to be stacked with this one.

weightsMap, optional

Map with stacking weights. Default is None.

nan_to_num: bool, optional

Non-finite values are replaced by zero if True. Default is True.

to_edisp_kernel_map(energy_axis)[source]#

Convert to map with energy dispersion kernels.

Parameters
energy_axisMapAxis

Reconstructed energy axis.

Returns
edispEDispKernelMap

Energy dispersion kernel map.

to_hdulist(format='gadf')#

Convert to HDUList.

Parameters
format{“gadf”, “gtpsf”}, optional

File format. Default is “gadf”.

Returns
hdu_listHDUList

HDU list.

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.

Parameters
regionSkyRegion or SkyCoord

Region or position where to get the map.

Returns
irfIRFMap

IRF map with region geometry.

write(filename, overwrite=False, format='gadf', checksum=False)#

Write IRF map to fits.

Parameters
filenamestr or Path

Filename to write to.

overwritebool, optional

Overwrite existing file. Default is False.

format{“gadf”, “gtpsf”}, optional

File format. Default is “gadf”.

checksumbool, optional

When True adds both DATASUM and CHECKSUM cards to the headers written to the file. Default is False.