WcsNDMap#

class gammapy.maps.WcsNDMap[source]#

Bases: WcsMap

WCS map with any number of non-spatial dimensions.

This class uses an ND numpy array to store map values. For maps with non-spatial dimensions and variable pixel size it will allocate an array with dimensions commensurate with the largest image plane.

Parameters:
geomWcsGeom

WCS geometry object.

datandarray

Data array. If none then an empty array will be allocated.

dtypestr, optional

Data type, default is float32

metadict

Dictionary to store meta data.

unitstr or Unit

The map unit

Methods Summary

binary_dilate(width[, kernel, use_fft])

Binary dilation of boolean mask adding a given margin.

binary_erode(width[, kernel, use_fft])

Binary erosion of boolean mask removing a given margin.

convolve(kernel[, method, mode])

Convolve map with a kernel.

crop(crop_width)

Crop the spatial dimensions of the map.

cutout(position, width[, mode, odd_npix, ...])

Create a cutout around a given position.

cutout_and_mask_region([region])

Compute cutout and mask for a given region of the map.

downsample(factor[, preserve_counts, ...])

Downsample the spatial dimension by a given factor.

fill_by_idx(idx[, weights])

Fill pixels at idx with given weights.

from_hdu(hdu[, hdu_bands, format])

Make a WcsNDMap object from a FITS HDU.

get_by_idx(idx)

Return map values at the given pixel indices.

interp_by_coord(coords[, method, ...])

Interpolate map values at the given map coordinates.

interp_by_pix(pix[, method, fill_value, ...])

Interpolate map values at the given pixel coordinates.

mask_contains_region(region)

Check if input region is contained in a boolean mask map.

plot([ax, fig, add_cbar, stretch, axes_loc, ...])

Plot image on matplotlib WCS axes.

plot_mask([ax])

Plot the mask as a shaded area.

set_by_idx(idx, vals)

Set pixels at idx with given vals.

smooth(width[, kernel])

Smooth the map.

stack(other[, weights, nan_to_num])

Stack cutout into map.

to_region_nd_map([region, func, weights, method])

Get region ND map in a given region.

to_region_nd_map_histogram([region, ...])

Convert map into region map by histogramming.

upsample(factor[, order, preserve_counts, ...])

Upsample the spatial dimension by a given factor.

Methods Documentation

binary_dilate(width, kernel='disk', use_fft=True)[source]#

Binary dilation of boolean mask adding a given margin.

Parameters:
widthtuple of Quantity

Angular sizes of the margin in (lon, lat) in that specific order. If only one value is passed, the same margin is applied in (lon, lat).

kernel{‘disk’, ‘box’}, optional

Kernel shape. Default is “disk”.

use_fftbool, optional

Use scipy.signal.fftconvolve if True. Otherwise, use scipy.ndimage.binary_dilation. Default is True.

Returns:
mapWcsNDMap

Dilated mask map.

binary_erode(width, kernel='disk', use_fft=True)[source]#

Binary erosion of boolean mask removing a given margin.

Parameters:
widthQuantity, str or float

If a float is given it interpreted as width in pixels. If an (angular) quantity is given it converted to pixels using geom.wcs.wcs.cdelt. The width corresponds to radius in case of a disk kernel, and the side length in case of a box kernel.

kernel{‘disk’, ‘box’}, optional

Kernel shape. Default is “disk”.

use_fftbool, optional

Use scipy.signal.fftconvolve if True. Otherwise, use scipy.ndimage.binary_erosion. Default is True.

Returns:
mapWcsNDMap

Eroded mask map.

convolve(kernel, method='fft', mode='same')[source]#

Convolve map with a kernel.

If the kernel is two-dimensional, it is applied to all image planes likewise. If the kernel is higher dimensional, it should either match the map in the number of dimensions or the map must be an image (no non-spatial axes). In that case, the corresponding kernel is selected and applied to every image plane or to the single input image respectively.

Parameters:
kernelPSFKernel or numpy.ndarray

Convolution kernel.

methodstr, optional

The method used by convolve. Default is ‘fft’.

modestr, optional

The convolution mode used by convolve. Default is ‘same’.

Returns:
mapWcsNDMap

Convolved map.

crop(crop_width)[source]#

Crop the spatial dimensions of the map.

Parameters:
crop_width{sequence, array_like, int}

Number of pixels cropped from the edges of each axis. Defined analogously to pad_with from numpy.pad.

Returns:
mapMap

Cropped map.

cutout(position, width, mode='trim', odd_npix=False, min_npix=1)[source]#

Create a cutout around a given position.

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”.

odd_npixbool, optional

Force width to odd number of pixels. Default is False.

min_npixbool, optional

Force width to a minimum number of pixels. Default is 1.

Returns:
cutoutWcsNDMap

Cutout map.

cutout_and_mask_region(region=None)[source]#

Compute cutout and mask for a given region of the map.

The function will estimate the minimal size of the cutout, which encloses the region.

Parameters:
region: `~regions.Region`, optional

Extended region. Default is None.

Returns:
cutout, masktuple of WcsNDMap

Cutout and mask map.

downsample(factor, preserve_counts=True, axis_name=None, weights=None)[source]#

Downsample the spatial dimension by a given factor.

Parameters:
factorint

Downsampling factor.

preserve_countsbool, optional

Preserve the integral over each bin. This should be set to True if the map is an integral quantity (e.g. counts) and False if the map is a differential quantity (e.g. intensity). Default is True.

axis_namestr, optional

Which axis to downsample. By default, spatial axes are downsampled. Default is None.

Returns:
mapMap

Downsampled map.

fill_by_idx(idx, weights=None)[source]#

Fill pixels at idx with given weights.

Parameters:
idxtuple

Tuple of pixel index arrays for each dimension of the map. Tuple should be ordered as (I_lon, I_lat, I_0, …, I_n) for WCS maps and (I_hpx, I_0, …, I_n) for HEALPix maps.

weightsndarray, optional

Weights vector. If None, weights are set to 1. Default is None.

classmethod from_hdu(hdu, hdu_bands=None, format=None)[source]#

Make a WcsNDMap object from a FITS HDU.

Parameters:
hduBinTableHDU or ImageHDU

The map FITS HDU.

hdu_bandsBinTableHDU

The BANDS table HDU.

format{‘gadf’, ‘fgst-ccube’,’fgst-template’}

FITS format convention.

Returns:
mapWcsNDMap

WCS map.

get_by_idx(idx)[source]#

Return map values at the given pixel indices.

Parameters:
idxtuple

Tuple of pixel index arrays for each dimension of the map. Tuple should be ordered as (I_lon, I_lat, I_0, …, I_n) for WCS maps and (I_hpx, I_0, …, I_n) for HEALPix maps.

Returns:
valsndarray

Array of pixel values. numpy.nan is used to flag coordinates outside the map.

interp_by_coord(coords, method='linear', fill_value=None, values_scale='lin')[source]#

Interpolate map values at the given map coordinates.

Parameters:
coordstuple, dict or MapCoord

Coordinate arrays for each dimension of the map. Tuple should be ordered as (lon, lat, x_0, …, x_n) where x_i are coordinates for non-spatial dimensions of the map. “lon” and “lat” are optional and will be taken at the center of the region by default.

method{“linear”, “nearest”}

Method to interpolate data values. By default linear interpolation is performed.

fill_valueNone or float value

The value to use for points outside of the interpolation domain. If None, values outside the domain are extrapolated.

values_scale{“lin”, “log”, “sqrt”}

Optional value scaling. Default is “lin”.

Returns:
valsndarray

Interpolated pixel values.

interp_by_pix(pix, method='linear', fill_value=None, values_scale='lin')[source]#

Interpolate map values at the given pixel coordinates.

Parameters:
pixtuple

Tuple of pixel coordinate arrays for each dimension of the map. Tuple should be ordered as (p_lon, p_lat, p_0, …, p_n) where p_i are pixel coordinates for non-spatial dimensions of the map.

method{“linear”, “nearest”}

Method to interpolate data values. Default is “linear”.

fill_valuefloat, optional

The value to use for points outside the interpolation domain. If None, values outside the domain are extrapolated. Default is None.

Returns:
valsndarray

Interpolated pixel values.

mask_contains_region(region)[source]#

Check if input region is contained in a boolean mask map.

Parameters:
region: `~regions.SkyRegion` or `~regions.PixRegion`

Region or list of Regions (pixel or sky regions accepted).

Returns:
containedbool

Whether region is contained in the mask.

plot(ax=None, fig=None, add_cbar=False, stretch='linear', axes_loc=None, kwargs_colorbar=None, **kwargs)[source]#

Plot image on matplotlib WCS axes.

Parameters:
axWCSAxes, optional

WCS axis object to plot on. Default is None.

figFigure, optional

Figure object. Default is None.

add_cbarbool, optional

Add color bar. Default is False.

stretchstr, optional
Passed to astropy.visualization.simple_norm.

Default is “linear”.

axes_locdict, optional

Keyword arguments passed to append_axes.

kwargs_colorbardict, optional

Keyword arguments passed to colorbar.

**kwargsdict

Keyword arguments passed to imshow.

Returns:
axWCSAxes

WCS axes object.

plot_mask(ax=None, **kwargs)[source]#

Plot the mask as a shaded area.

Parameters:
axWCSAxes, optional

WCS axis object to plot on. Default is None.

**kwargsdict

Keyword arguments passed to contourf.

Returns:
axWCSAxes, optional

WCS axis object to plot on.

set_by_idx(idx, vals)[source]#

Set pixels at idx with given vals.

Parameters:
idxtuple

Tuple of pixel index arrays for each dimension of the map. Tuple should be ordered as (I_lon, I_lat, I_0, …, I_n) for WCS maps and (I_hpx, I_0, …, I_n) for HEALPix maps.

valsndarray

Values vector.

smooth(width, kernel='gauss', **kwargs)[source]#

Smooth the map.

Iterates over 2D image planes, processing one at a time.

Parameters:
widthQuantity, str or float

Smoothing width given as quantity or float. If a float is given it interpreted as smoothing width in pixels. If an (angular) quantity is given it converted to pixels using geom.wcs.wcs.cdelt. It corresponds to the standard deviation in case of a Gaussian kernel, the radius in case of a disk kernel, and the side length in case of a box kernel.

kernel{‘gauss’, ‘disk’, ‘box’}, optional

Kernel shape. Default is “gauss”.

kwargsdict

Keyword arguments passed to uniform_filter (‘box’), gaussian_filter (‘gauss’) or convolve (‘disk’).

Returns:
imageWcsNDMap

Smoothed image (a copy, the original object is unchanged).

stack(other, weights=None, nan_to_num=True)[source]#

Stack cutout into map.

Parameters:
otherWcsNDMap

Other map to stack.

weightsWcsNDMap, optional

Array to be used as weights. The spatial geometry must be equivalent to other and additional axes must be broadcastable. Default is None.

nan_to_num: bool, optional

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

to_region_nd_map(region=None, func=<function nansum>, weights=None, method='nearest')[source]#

Get region ND map in a given region.

By default, the whole map region is considered.

Parameters:
region: `~regions.Region` or `~astropy.coordinates.SkyCoord`, optional

Region. Default is None.

funcnumpy.func, optional

Function to reduce the data. Default is np.nansum. For boolean Map, use np.any or np.all.

weightsWcsNDMap, optional

Array to be used as weights. Weights are multiplied with the input data and should be defined accordingly, in case a weighted average is desired. If a mask is given, the mask is applied to the input map. The geometry must be equivalent. Default is None.

method{“nearest”, “linear”}, optional

How to interpolate if a position is given. Default is “nearest”.

Returns:
spectrumRegionNDMap

Spectrum in the given region.

to_region_nd_map_histogram(region=None, bins_axis=None, nbin=100, density=False)[source]#

Convert map into region map by histogramming.

By default, it creates a linearly spaced axis with 100 bins between (-max(abs(data)), max(abs(data))) within the given region.

Parameters:
region: `~regions.Region`, optional

Region to histogram over. Default is None.

bins_axisMapAxis, optional

Binning of the histogram. Default is None.

nbinint, optional

Number of bins to use if no bins_axis is given. Default is 100.

densitybool, optional

Normalize integral of the histogram to 1. Default is False.

Returns:
region_mapRegionNDMap

Region map with histogram.

Examples

This is how to use the method to create energy dependent histograms:

from gammapy.maps import MapAxis, Map
import numpy as np

random_state = np.random.RandomState(seed=0)

energy_axis = MapAxis.from_energy_bounds("1 TeV", "10 TeV", nbin=3)

data = Map.create(axes=[energy_axis], width=10, unit="cm2 s-1", binsz=0.02)
data.data = random_state.normal(
    size=data.data.shape, loc=0, scale=np.array([1.0, 2.0, 3.0]).reshape((-1, 1, 1))
)

hist = data.to_region_nd_map_histogram()
hist.plot(axis_name="bins")
upsample(factor, order=0, preserve_counts=True, axis_name=None)[source]#

Upsample the spatial dimension by a given factor.

Parameters:
factorint

Upsampling factor.

orderint, optional

Order of the interpolation used for upsampling. Default is 0.

preserve_countsbool, optional

Preserve the integral over each bin. This should be true if the map is an integral quantity (e.g. counts) and false if the map is a differential quantity (e.g. intensity). Default is True.

axis_namestr, optional

Which axis to upsample. By default, spatial axes are upsampled. Default is None.

Returns:
mapMap

Upsampled map.

__init__(geom, data=None, dtype='float32', meta=None, unit='')[source]#
classmethod __new__(*args, **kwargs)#