RegionGeom#
- class gammapy.maps.RegionGeom(region, axes=None, wcs=None, binsz_wcs='0.1 deg')[source]#
Bases:
gammapy.maps.geom.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
- region
SkyRegion
Region object.
- axeslist of
MapAxis
Non-spatial data axes.
- wcs
WCS
Optional wcs object to project the region if needed.
- binsz_wcs
float
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.
- region
Attributes Summary
If the geom contains an energy axis rename it to energy true
List of non-spatial axes.
All axes names
Angular bin size of the underlying
WcsGeom
Pixel values corresponding to the center of the region
Sky coordinate of the center of the region
Shape of the Numpy data array matching this geometry.
Shape of data of the non-spatial axes and unit spatial axes.
Coordinate system, either Galactic ("galactic") or Equatorial ("icrs").
Whether geom has an energy axis
Whether regions are all point regions
Whether the geom non spatial axes have length 1, equivalent to an image.
Whether the geom is an image without extra dimensions.
Number of spatial pixels
SkyRegion
object that defines the spatial component of the region geometryWCS projection object.
Width of bounding box of the region.
Methods Summary
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
is_allclose
(other[, rtol_axes, atol_axes])Compare two data IRFs for equivalency
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, kwargs_point, path_effect])Plot region in the sky.
rename_axes
(names, new_names)Rename axes contained in the geometry
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.
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
- center_coord#
- 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_all_point_sky_regions#
Whether regions are all point regions
- is_allsky = False#
- is_flat#
Whether the geom non spatial axes have length 1, 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'#
- wcs#
WCS projection object.
- width#
Width of bounding box of the region.
- Returns
- width
Quantity
Dimensions of the region in both spatial dimensions. Units:
deg
- width
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
- volume
Quantity
Volume of the region.
- volume
- contains(coords)[source]#
Check if a given map coordinate is contained in the region. Requires the
region
attribute to be set.For
PointSkyRegion
the method always returns true.
- contains_pix(pix)#
Check if a given pixel coordinate is contained in the geometry.
- Parameters
- pixtuple
Tuple of pixel coordinates.
- Returns
- containment
ndarray
Bool array.
- containment
- contains_wcs_pix(pix)[source]#
Check if a given wcs pixel coordinate is contained in the region.
For
PointSkyRegion
the method always returns true.- Parameters
- pixtuple
Tuple of pixel coordinates.
- Returns
- containment
ndarray
Bool array.
- containment
- 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.
- coordstuple or
- 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
- copy
Geom
Copied map geometry.
- copy
- 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__
- regionstr or
- Returns
- geom
RegionGeom
Region geometry
- geom
- 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
- geom
Geom
Cropped geometry.
- geom
- 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
- memory
Quantity
Estimated memory usage in megabytes (MB)
- memory
- downsample(factor, axis_name)[source]#
Downsample a non-spatial dimension of the region by a given factor.
- Returns
- region
RegionGeom
RegionGeom with the downsampled axis.
- region
- drop(axis_name)#
Drop an axis from the geom.
- Parameters
- axis_namestr
Name of the axis to remove.
- Returns
- geom
Geom
New geom with the axis removed.
- geom
- 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.
- classmethod from_hdulist(hdulist, format='ogip', hdu=None)[source]#
Read region table and convert it to region list.
- Parameters
- hdulist
HDUList
HDU list
- format{“ogip”, “ogip-arf”, “gadf”}
HDU format
- hdulist
- Returns
- geom
RegionGeom
Region map geometry
- geom
- 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
- regionslist of
- Returns
- geom
RegionGeom
Region map geometry
- geom
- 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
- coord
MapCoord
Map coordinate object.
- coord
- 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_coord
MapCoord
MapCoord object with the coordinates inside the region.
- weights
array
Weights representing the fraction of each pixel contained in the region.
- region_coord
- is_allclose(other, rtol_axes=1e-06, atol_axes=1e-06)[source]#
Compare two data IRFs for equivalency
- Parameters
- other
RegionGeom
Geom to compare against.
- rtol_axesfloat
Relative tolerance for the axes comparison.
- atol_axesfloat
Relative tolerance for the axes comparison.
- other
- Returns
- is_allclosebool
Whether the geometry is all close.
- 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
- geom
Geom
Padded geometry.
- geom
- 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_point=None, path_effect=None, **kwargs)[source]#
Plot region in the sky.
- Parameters
- ax
WCSAxes
Axes to plot on. If no axes are given, the region is shown using the minimal equivalent WCS geometry.
- kwargs_pointdict
Keyword arguments passed to
Line2D
for plotting of point sources- path_effect
PathEffect
Path effect applied to artists and lines.
- **kwargsdict
Keyword arguments forwarded to
as_artist
- ax
- Returns
- ax
WCSAxes
Axes to plot on.
- ax
- rename_axes(names, new_names)#
Rename axes contained in the geometry
- Parameters
- nameslist or str
Names of the axes.
- new_nameslist or str
New names of the axes (list must be of same length than
names
).
- Returns
- geom
Geom
Renamed geometry.
- geom
- replace_axis(axis)#
Replace axis with a new one.
- resample_axis(axis)#
Resample geom to a new axis binning.
This method groups the existing bins into a new binning.
- separation(position)[source]#
Angular distance between the center of the region and the given position.
- Parameters
- position
astropy.coordinates.SkyCoord
Sky coordinate we want the angular distance to.
- position
- Returns
- sep
Angle
The on-sky separation between the given coordinate and the region center.
- sep
- slice_by_idx(slices)#
Create a new geometry by slicing the non-spatial axes.
- solid_angle()[source]#
Get solid angle of the region.
- Returns
- angle
Quantity
Solid angle of the region. In sr. Units:
sr
- angle
- squash(axis_name)#
Squash geom axis.
- Parameters
- axis_namestr
Axis to squash.
- Returns
- geom
Geom
Geom with squashed axis.
- geom
- to_bands_hdu(hdu_bands=None, format='gadf')#
- to_binsz_wcs(binsz)[source]#
Change the bin size of the underlying WCS geometry.
- Parameters
- binzsfloat, string or
Quantity
- binzsfloat, string or
- Returns
- region
RegionGeom
A RegionGeom with the same axes and region as the input, but different wcs pixelization.
- region
- to_cube(axes)[source]#
Append non-spatial axes to create a higher-dimensional geometry.
- Returns
- region
RegionGeom
RegionGeom with the added axes.
- region
- 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
- hdulist
HDUList
HDU list
- hdulist
- to_image()[source]#
Remove non-spatial axes to create a 2D region.
- Returns
- region
RegionGeom
RegionGeom without any non-spatial axes.
- region
- to_wcs_geom(width_min=None)[source]#
Get the minimal equivalent geometry which contains the region.
- Parameters
- width_min
Quantity
Minimal width for the resulting geometry. Can be a single number or two, for different minimum widths in each spatial dimension.
- width_min
- Returns
- wcs_geom
WcsGeom
A WCS geometry object.
- wcs_geom
- upsample(factor, axis_name=None)[source]#
Upsample a non-spatial dimension of the region by a given factor.
- Returns
- region
RegionGeom
RegionGeom with the upsampled axis.
- region