RegionGeom

class gammapy.maps.RegionGeom(region, axes=None, wcs=None, binsz_wcs='0.1 deg')[source]

Bases: gammapy.maps.Geom

Map geometry representing a region on the sky. The spatial component of the geometry is made up of a single pixel with an arbitrary shape and size. It can also have any number of non-spatial dimensions. This class represents the geometry for the RegionNDMap maps.

Parameters
regionSkyRegion

Region object.

axeslist of MapAxis

Non-spatial data axes.

wcsWCS

Optional wcs object to project the region if needed.

binsz_wcsfloat

Angular bin size of the underlying WcsGeom used to evaluate quantities in the region. Default size is 0.01 deg. This default value is adequate for the majority of use cases. If a wcs object is provided, the input of binsz_wcs is overridden.

Attributes Summary

as_energy_true

If the geom contains an energy axis rename it to energy true

axes

List of non-spatial axes.

axes_names

All axes names

binsz_wcs

Angular bin size of the underlying WcsGeom

center_coord

(astropy.coordinates.SkyCoord)

center_pix

Pixel values corresponding to the center of the region

center_skydir

Sky coordinate of the center of the region

data_shape

Shape of the Numpy data array matching this geometry.

data_shape_axes

Shape of data of the non-spatial axes and unit spatial axes.

frame

Coordinate system, either Galactic (“galactic”) or Equatorial (“icrs”).

has_energy_axis

Whether geom has an energy axis

is_allsky

is_flat

Whether the geom non spatial axes have length 1, i.e.

is_hpx

is_image

Whether the geom is an image without extra dimensions.

is_region

is_regular

npix

Number of spatial pixels

projection

region

SkyRegion object that defines the spatial component of the region geometry

wcs

WCS projection object.

width

Width of bounding box of the region.

Methods Summary

bin_volume()

If the RegionGeom has a non-spatial axis, it returns the volume of the region.

contains(coords)

Check if a given map coordinate is contained in the region.

contains_pix(pix)

Check if a given pixel coordinate is contained in the geometry.

contains_wcs_pix(pix)

Check if a given wcs pixel coordinate is contained in the region.

coord_to_idx(coords[, clip])

Convert map coordinates to pixel indices.

coord_to_pix(coords)

Convert map coordinates to pixel coordinates.

copy(**kwargs)

Copy and overwrite given attributes.

create(region, **kwargs)

Create region geometry.

crop()

Crop the geometry at the edges.

data_nbytes([dtype])

Estimate memory usage in megabytes of the Numpy data array matching this geometry depending on the given type.

downsample(factor, axis_name)

Downsample a non-spatial dimension of the region by a given factor.

drop(axis_name)

Drop an axis from the geom.

energy_mask([energy_min, energy_max, …])

Create a mask for a given energy range.

from_hdulist(hdulist[, format, hdu])

Read region table and convert it to region list.

from_regions(regions, **kwargs)

Create region geom from list of regions

get_coord([mode, frame, sparse, axis_name])

Get map coordinates from the geometry.

get_idx()

Get tuple of pixel indices for this geometry.

get_wcs_coord_and_weights([factor])

Get the array of spatial coordinates and corresponding weights

pad(pad_width, axis_name)

Pad the geometry at the edges.

pix_to_coord(pix)

Convert pixel coordinates to map coordinates.

pix_to_idx(pix[, clip])

Convert pixel coordinates to pixel indices.

plot_region([ax])

Plot region in the sky.

replace_axis(axis)

Replace axis with a new one.

resample_axis(axis)

Resample geom to a new axis binning.

separation(position)

Angular distance between the center of the region and the given position.

slice_by_idx(slices)

Create a new geometry by slicing the non-spatial axes.

solid_angle()

Get solid angle of the region.

squash(axis_name)

Squash geom axis.

to_bands_hdu([hdu_bands, format])

to_binsz(binsz)

Returns self

to_binsz_wcs(binsz)

Change the bin size of the underlying WCS geometry.

to_cube(axes)

Append non-spatial axes to create a higher-dimensional geometry.

to_hdulist([format, hdu_bands, hdu_region])

Convert geom to hdulist

to_image()

Remove non-spatial axes to create a 2D region.

to_wcs_geom([width_min])

Get the minimal equivalent geometry which contains the region.

union(other)

Stack a RegionGeom by making the union

upsample(factor[, axis_name])

Upsample a non-spatial dimension of the region by a given factor.

Attributes Documentation

as_energy_true

If the geom contains an energy axis rename it to energy true

axes

List of non-spatial axes.

axes_names

All axes names

binsz_wcs

Angular bin size of the underlying WcsGeom

Returns
binsz_wcs: Angle
center_coord

(astropy.coordinates.SkyCoord)

center_pix

Pixel values corresponding to the center of the region

center_skydir

Sky coordinate of the center of the region

data_shape

Shape of the Numpy data array matching this geometry.

data_shape_axes

Shape of data of the non-spatial axes and unit spatial axes.

frame

Coordinate system, either Galactic (“galactic”) or Equatorial (“icrs”).

has_energy_axis

Whether geom has an energy axis

is_allsky = False
is_flat

Whether the geom non spatial axes have length 1, i.e. if the geom is equivalent to an image.

is_hpx = False
is_image

Whether the geom is an image without extra dimensions.

is_region = True
is_regular = True
npix

Number of spatial pixels

projection = 'TAN'
region

SkyRegion object that defines the spatial component of the region geometry

wcs

WCS projection object.

width

Width of bounding box of the region.

Returns
widthQuantity

Dimensions of the region in both spatial dimensions. Units: deg

Methods Documentation

bin_volume()[source]

If the RegionGeom has a non-spatial axis, it returns the volume of the region. If not, it just returns the solid angle size.

Returns
volumeQuantity

Volume of the region.

contains(coords)[source]

Check if a given map coordinate is contained in the region. Requires the region attribute to be set.

Parameters
coordstuple, dict, MapCoord or SkyCoord

Object containing coordinate arrays we wish to check for inclusion in the region.

Returns
maskndarray

Boolean Numpy array with the same shape as the input that indicates which coordinates are inside the region.

contains_pix(pix)

Check if a given pixel coordinate is contained in the geometry.

Parameters
pixtuple

Tuple of pixel coordinates.

Returns
containmentndarray

Bool array.

contains_wcs_pix(pix)[source]

Check if a given wcs pixel coordinate is contained in the region.

Parameters
pixtuple

Tuple of pixel coordinates.

Returns
containmentndarray

Bool array.

coord_to_idx(coords, clip=False)

Convert map coordinates to pixel indices.

Parameters
coordstuple or MapCoord

Coordinate values in each dimension of the map. This can either be a tuple of numpy arrays or a MapCoord object. If passed as a tuple then the ordering should be (longitude, latitude, c_0, …, c_N) where c_i is the coordinate vector for axis i.

clipbool

Choose whether to clip indices to the valid range of the geometry. If false then indices for coordinates outside the geometry range will be set -1.

Returns
pixtuple

Tuple of pixel indices in image and band dimensions. Elements set to -1 correspond to coordinates outside the map.

coord_to_pix(coords)[source]

Convert map coordinates to pixel coordinates.

Parameters
coordstuple

Coordinate values in each dimension of the map. This can either be a tuple of numpy arrays or a MapCoord object. If passed as a tuple then the ordering should be (longitude, latitude, c_0, …, c_N) where c_i is the coordinate vector for axis i.

Returns
pixtuple

Tuple of pixel coordinates in image and band dimensions.

copy(**kwargs)

Copy and overwrite given attributes.

Parameters
**kwargsdict

Keyword arguments to overwrite in the map geometry constructor.

Returns
copyGeom

Copied map geometry.

classmethod create(region, **kwargs)[source]

Create region geometry.

The input region can be passed in the form of a ds9 string and will be parsed internally by parse. See:

Parameters
regionstr or SkyRegion

Region definition

**kwargsdict

Keyword arguments passed to RegionGeom.__init__

Returns
geomRegionGeom

Region geometry

crop()[source]

Crop the geometry at the edges.

Parameters
crop_width{sequence, array_like, int}

Number of values cropped from the edges of each axis.

Returns
geomGeom

Cropped geometry.

data_nbytes(dtype='float32')

Estimate memory usage in megabytes of the Numpy data array matching this geometry depending on the given type.

Parameters
dtypedata-type

The desired data-type for the array. Default is “float32”

Returns
memoryQuantity

Estimated memory usage in megabytes (MB)

downsample(factor, axis_name)[source]

Downsample a non-spatial dimension of the region by a given factor.

Returns
regionRegionGeom

RegionGeom with the downsampled axis.

drop(axis_name)

Drop an axis from the geom.

Parameters
axis_namestr

Name of the axis to remove.

Returns
geomGeom

New geom with the axis removed.

energy_mask(energy_min=None, energy_max=None, round_to_edge=False)

Create a mask for a given energy range.

The energy bin must be fully contained to be included in the mask.

Parameters
energy_min, energy_maxQuantity

Energy range

Returns
maskndarray

Energy mask

classmethod from_hdulist(hdulist, format='ogip', hdu=None)[source]

Read region table and convert it to region list.

Parameters
hdulistHDUList

HDU list

format{“ogip”, “ogip-arf”, “gadf”}

HDU format

Returns
geomRegionGeom

Region map geometry

classmethod from_regions(regions, **kwargs)[source]

Create region geom from list of regions

The regions are combined with union to a compound region.

Parameters
regionslist of SkyRegion or str

Regions

**kwargs: dict

Keyword arguments forwarded to RegionGeom

Returns
geomRegionGeom

Region map geometry

get_coord(mode='center', frame=None, sparse=False, axis_name=None)[source]

Get map coordinates from the geometry.

Parameters
mode{‘center’, ‘edges’}

Get center or edge coordinates for the non-spatial axes.

framestr or Frame

Coordinate frame

sparsebool

Compute sparse coordinates

axis_namestr

If mode = “edges”, the edges will be returned for this axis only.

Returns
coordMapCoord

Map coordinate object.

get_idx()[source]

Get tuple of pixel indices for this geometry.

Returns all pixels in the geometry by default. Pixel indices for a single image plane can be accessed by setting idx to the index tuple of a plane.

Parameters
idxtuple, optional

A tuple of indices with one index for each non-spatial dimension. If defined only pixels for the image plane with this index will be returned. If none then all pixels will be returned.

localbool

Flag to return local or global pixel indices. Local indices run from 0 to the number of pixels in a given image plane.

flatbool, optional

Return a flattened array containing only indices for pixels contained in the geometry.

Returns
idxtuple

Tuple of pixel index vectors with one vector for each dimension.

get_wcs_coord_and_weights(factor=10)[source]

Get the array of spatial coordinates and corresponding weights

The coordinates are the center of a pixel that intersects the region and the weights that represent which fraction of the pixel is contained in the region.

Parameters
factorint

Oversampling factor to compute the weights

Returns
region_coordMapCoord

MapCoord object with the coordinates inside the region.

weightsarray

Weights representing the fraction of each pixel contained in the region.

pad(pad_width, axis_name)

Pad the geometry at the edges.

Parameters
pad_width{sequence, array_like, int}

Number of values padded to the edges of each axis.

axis_namestr

Name of the axis to pad.

Returns
geomGeom

Padded geometry.

pix_to_coord(pix)[source]

Convert pixel coordinates to map coordinates.

Parameters
pixtuple

Tuple of pixel coordinates.

Returns
coordstuple

Tuple of map coordinates.

pix_to_idx(pix, clip=False)[source]

Convert pixel coordinates to pixel indices.

Returns -1 for pixel coordinates that lie outside of the map.

Parameters
pixtuple

Tuple of pixel coordinates.

clipbool

Choose whether to clip indices to the valid range of the geometry. If false then indices for coordinates outside the geometry range will be set -1.

Returns
idxtuple

Tuple of pixel indices.

plot_region(ax=None, **kwargs)[source]

Plot region in the sky.

Parameters
axWCSAxes

Axes to plot on. If no axes are given, the region is shown using the minimal equivalent WCS geometry.

**kwargsdict

Keyword arguments forwarded to as_artist

Returns
axWCSAxes

Axes to plot on.

replace_axis(axis)

Replace axis with a new one.

Parameters
axisMapAxis

New map axis.

Returns
mapGeom

Geom with replaced axis.

resample_axis(axis)

Resample geom to a new axis binning.

This method groups the existing bins into a new binning.

Parameters
axisMapAxis

New map axis.

Returns
mapGeom

Geom with resampled axis.

separation(position)[source]

Angular distance between the center of the region and the given position.

Parameters
positionastropy.coordinates.SkyCoord

Sky coordinate we want the angular distance to.

Returns
sepAngle

The on-sky separation between the given coordinate and the region center.

slice_by_idx(slices)

Create a new geometry by slicing the non-spatial 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
geomGeom

Sliced geometry.

solid_angle()[source]

Get solid angle of the region.

Returns
angleQuantity

Solid angle of the region. In sr. Units: sr

squash(axis_name)

Squash geom axis.

Parameters
axis_namestr

Axis to squash.

Returns
geomGeom

Geom with squashed axis.

to_bands_hdu(hdu_bands=None, format='gadf')
to_binsz(binsz)[source]

Returns self

to_binsz_wcs(binsz)[source]

Change the bin size of the underlying WCS geometry.

Parameters
binzsfloat, string or Quantity
Returns
regionRegionGeom

A RegionGeom with the same axes and region as the input, but different wcs pixelization.

to_cube(axes)[source]

Append non-spatial axes to create a higher-dimensional geometry.

Returns
regionRegionGeom

RegionGeom with the added axes.

to_hdulist(format='ogip', hdu_bands=None, hdu_region=None)[source]

Convert geom to hdulist

Parameters
format{“gadf”, “ogip”, “ogip-sherpa”}

HDU format

hdustr

Name of the HDU with the map data.

Returns
hdulistHDUList

HDU list

to_image()[source]

Remove non-spatial axes to create a 2D region.

Returns
regionRegionGeom

RegionGeom without any non-spatial axes.

to_wcs_geom(width_min=None)[source]

Get the minimal equivalent geometry which contains the region.

Parameters
width_minQuantity

Minimal width for the resulting geometry. Can be a single number or two, for different minimum widths in each spatial dimension.

Returns
wcs_geomWcsGeom

A WCS geometry object.

union(other)[source]

Stack a RegionGeom by making the union

upsample(factor, axis_name=None)[source]

Upsample a non-spatial dimension of the region by a given factor.

Returns
regionRegionGeom

RegionGeom with the upsampled axis.