RegionNDMap#
- class gammapy.maps.RegionNDMap(geom, data=None, dtype='float32', meta=None, unit='')[source]#
Bases:
gammapy.maps.core.Map
N-dimensional region map. A
RegionNDMap
owns aRegionGeom
instance as well as a data array containing the values associated to that region in the sky along the non-spatial axis, usually an energy axis. The spatial dimensions of aRegionNDMap
are reduced to a single spatial bin with an arbitrary shape, and any extra dimensions are described by an arbitrary number of non-spatial axes.- Parameters
- geom
RegionGeom
Region geometry object.
- data
ndarray
Data array. If none then an empty array will be allocated.
- dtypestr, optional
Data type, default is float32
- meta
dict
Dictionary to store meta data.
- unitstr or
Unit
The map unit
- geom
Attributes Summary
Data array (
ndarray
)Map geometry (
Geom
)Whether map is mask with bool dtype
Map meta (
dict
)Map data times unit (
Quantity
)Map unit (
Unit
)Methods Summary
apply_edisp
(edisp)Apply energy dispersion to map.
coadd
(map_in[, weights])Add the contents of
map_in
to this map.copy
(**kwargs)Copy map instance and overwrite given attributes, except for geometry.
create
(region[, axes, dtype, meta, unit, ...])Create an empty region map object.
crop
()Crop the spatial dimensions of the map.
cumsum
(axis_name)Compute cumulative sum along a given axis
cutout
(*args, **kwargs)Return self
downsample
(factor[, preserve_counts, ...])Downsample the non-spatial dimension by a given factor.
fill_by_coord
(coords[, weights])Fill pixels at
coords
with givenweights
.fill_by_idx
(idx[, weights])Fill pixels at
idx
with givenweights
.fill_by_pix
(pix[, weights])Fill pixels at
pix
with givenweights
.fill_events
(events)Fill event coordinates (
EventList
).from_geom
(geom[, meta, data, unit, dtype])Generate an empty map from a
Geom
instance.from_hdulist
(hdulist[, format, ogip_column, hdu])Create from
HDUList
.from_stack
(maps[, axis, axis_name])Create Map from list of images and a non-spatial axis.
from_table
(table[, format, colname])Create region map from table
get_by_coord
(coords[, fill_value])Return map values at the given map coordinates.
get_by_idx
(idxs)Return map values at the given pixel indices.
get_by_pix
(pix[, fill_value])Return map values at the given pixel coordinates.
get_image_by_coord
(coords)Return spatial map at the given axis coordinates.
get_image_by_idx
(idx)Return spatial map at the given axis pixel indices.
get_image_by_pix
(pix)Return spatial map at the given axis pixel coordinates
get_spectrum
(*args, **kwargs)Return self
integral
(axis_name, coords, **kwargs)Compute integral along a given axis
interp_by_coord
(coords, **kwargs)Interpolate map values at the given map coordinates.
interp_by_pix
(pix, **kwargs)Interpolate map values at the given pixel coordinates.
interp_to_geom
(geom[, preserve_counts, ...])Interpolate map to input geometry.
is_allclose
(other[, rtol_axes, atol_axes])Compare two Maps for close equivalency
iter_by_axis
(axis_name[, keepdims])"Iterate over a given axis
iter_by_axis_data
(axis_name)Iterate data by axis
iter_by_image
([keepdims])Iterate over image planes of a map.
Iterate over image planes of the map.
mask_nearest_position
(position)Given a sky coordinate return nearest valid position in the mask
normalize
([axis_name])Normalise data in place along a given axis.
pad
(pad_width[, axis_name, mode, cval, method])Pad the spatial dimensions of the map.
plot
([ax, axis_name])Plot the data contained in region map along the non-spatial axis.
plot_grid
([figsize, ncols])Plot map as a grid of subplots for non-spatial axes
plot_hist
([ax])Plot as histogram.
Plot map with interactive widgets to explore the non spatial axes.
plot_mask
([ax])Plot the mask as a shaded area in a xmin-xmax range
plot_region
([ax])Plot region
read
(filename[, format, ogip_column, hdu])Read from file.
reduce
(axis_name[, func, keepdims, weights])Reduce map over a single non-spatial axis
reduce_over_axes
([func, keepdims, ...])Reduce map over non-spatial axes
rename_axes
(names, new_names)Rename the Map axes.
reproject_to_geom
(geom[, preserve_counts, ...])Reproject map to input geometry.
resample
(geom[, weights, preserve_counts])Resample pixels to
geom
with givenweights
.resample_axis
(axis[, weights, ufunc])Resample map to a new axis by grouping and reducing smaller bins by a given ufunc
sample_coord
(n_events[, random_state])Sample position and energy of events.
set_by_coord
(coords, vals)Set pixels at
coords
with givenvals
.set_by_idx
(idx, value)Set pixels at
idx
with givenvals
.set_by_pix
(pix, vals)Set pixels at
pix
with givenvals
.slice_by_idx
(slices)Slice sub map from map object.
split_by_axis
(axis_name)Split a Map along an axis into multiple maps.
stack
(other[, weights, nan_to_num])Stack other region map into map.
sum_over_axes
([axes_names, keepdims, weights])To sum map values over all non-spatial axes.
to_cube
(axes)Append non-spatial axes to create a higher-dimensional Map.
to_hdulist
([format, hdu, hdu_bands, hdu_region])Convert to
HDUList
.to_region_nd_map
(*args, **kwargs)to_table
([format])Convert to
Table
.to_unit
(unit)Convert map to different unit
upsample
(factor[, order, preserve_counts, ...])Upsample the non-spatial dimension by a given factor.
write
(filename[, overwrite, format, hdu])Write map to file
Attributes Documentation
- is_mask#
Whether map is mask with bool dtype
- tag = 'map'#
Methods Documentation
- apply_edisp(edisp)#
Apply energy dispersion to map. Requires energy axis.
- Parameters
- edisp
gammapy.irf.EDispKernel
Energy dispersion matrix
- edisp
- Returns
- map
WcsNDMap
Map with energy dispersion applied.
- map
- coadd(map_in, weights=None)#
Add the contents of
map_in
to this map.This method can be used to combine maps containing integral quantities (e.g. counts) or differential quantities if the maps have the same binning.
- Parameters
- map_in
Map
Input map.
- weights: `Map` or `~numpy.ndarray`
The weight factors while adding
- map_in
- copy(**kwargs)#
Copy map instance and overwrite given attributes, except for geometry.
- Parameters
- **kwargsdict
Keyword arguments to overwrite in the map constructor.
- Returns
- copy
Map
Copied Map.
- copy
- classmethod create(region, axes=None, dtype='float32', meta=None, unit='', wcs=None, binsz_wcs='0.1deg', data=None)[source]#
Create an empty region map object.
- Parameters
- regionstr or
SkyRegion
Region specification
- axeslist of
MapAxis
Non spatial axes.
- dtypestr
Data type, default is ‘float32’
- unitstr or
Unit
Data unit.
- meta
dict
Dictionary to store meta data.
- wcs
WCS
WCS projection to use for local projections of the region
- binsz_wcs: `~astropy.units.Quantity` or str
Bin size used for the default WCS, if wcs=None.
- data
ndarray
Data array
- regionstr or
- Returns
- map
RegionNDMap
Region map
- map
- cumsum(axis_name)#
Compute cumulative sum along a given axis
- Parameters
- axis_namestr
Along which axis to integrate.
- Returns
- cumsum
Map
Map with cumulative sum
- cumsum
- downsample(factor, preserve_counts=True, axis_name=None, weights=None)[source]#
Downsample the non-spatial dimension by a given factor.
By default the first axes is downsampled.
- Parameters
- factorint
Downsampling factor.
- preserve_countsbool
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).
- axis_namestr
Which axis to downsample. Default is “energy”.
- weights
RegionNDMap
Contains the weights to apply to the axis to reduce. Default is just weighs of one.
- Returns
- map
RegionNDMap
Downsampled region map.
- map
- fill_by_coord(coords, weights=None)#
Fill pixels at
coords
with givenweights
.
- fill_by_idx(idx, weights=None)[source]#
Fill pixels at
idx
with givenweights
.- 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.
- weights
ndarray
Weights vector. Default is weight of one.
- fill_by_pix(pix, weights=None)#
Fill pixels at
pix
with givenweights
.- Parameters
- pixtuple
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. Pixel indices can be either float or integer type. Float indices will be rounded to the nearest integer.
- weights
ndarray
Weights vector. Default is weight of one.
- static from_geom(geom, meta=None, data=None, unit='', dtype='float32')#
Generate an empty map from a
Geom
instance.- Parameters
- geom
Geom
Map geometry.
- data
numpy.ndarray
data array
- meta
dict
Dictionary to store meta data.
- unitstr or
Unit
Data unit.
- geom
- Returns
- map_out
Map
Map object
- map_out
- classmethod from_hdulist(hdulist, format='gadf', ogip_column=None, hdu=None, **kwargs)[source]#
Create from
HDUList
.- Parameters
- hdulist
HDUList
HDU list.
- format{“gadf”, “ogip”, “ogip-arf”}
Format specification
- ogip_column{“COUNTS”, “QUALITY”, “BACKSCAL”}
If format ‘ogip’ is chosen which table hdu column to read.
- hdustr
Name or index of the HDU with the map data.
- hdulist
- Returns
- region_nd_map
RegionNDMap
Region map.
- region_nd_map
- classmethod from_stack(maps, axis=None, axis_name=None)#
Create Map from list of images and a non-spatial axis.
The image geometries must be aligned, except for the axis that is stacked.
- Parameters
- Returns
- map
Map
Map with additional non-spatial axis.
- map
- classmethod from_table(table, format='', colname=None)[source]#
Create region map from table
- Parameters
- table
Table
Table with input data
- format{“gadf-sed”, “lightcurve”, “profile”}
Format to use
- colnamestr
Column name to take the data from.
- table
- Returns
- region_map
RegionNDMap
Region map
- region_map
- get_by_coord(coords, fill_value=nan)#
Return map values at the given map coordinates.
- Parameters
- coordstuple 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.
- fill_valuefloat
Value which is returned if the position is outside of the projection footprint
- coordstuple or
- Returns
- vals
ndarray
Values of pixels in the map. np.nan used to flag coords outside of map.
- vals
- get_by_idx(idxs)[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
- vals
ndarray
Array of pixel values. np.nan used to flag coordinate outside of map
- vals
- get_by_pix(pix, fill_value=nan)#
Return map values at the given pixel coordinates.
- Parameters
- pixtuple
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. Pixel indices can be either float or integer type.
- fill_valuefloat
Value which is returned if the position is outside of the projection footprint
- Returns
- vals
ndarray
Array of pixel values. np.nan used to flag coordinates outside of map
- vals
- get_image_by_coord(coords)#
Return spatial map at the given axis coordinates.
- Parameters
- coordstuple or dict
Tuple should be ordered as (x_0, …, x_n) where x_i are coordinates for non-spatial dimensions of the map. Dict should specify the axis names of the non-spatial axes such as {‘axes0’: x_0, …, ‘axesn’: x_n}.
- Returns
- map_out
Map
Map with spatial dimensions only.
- map_out
See also
Examples
import numpy as np from gammapy.maps import Map, MapAxis from astropy.coordinates import SkyCoord from astropy import units as u # Define map axes energy_axis = MapAxis.from_edges( np.logspace(-1., 1., 4), unit='TeV', name='energy', ) time_axis = MapAxis.from_edges( np.linspace(0., 10, 20), unit='h', name='time', ) # Define map center skydir = SkyCoord(0, 0, frame='galactic', unit='deg') # Create map m_wcs = Map.create( map_type='wcs', binsz=0.02, skydir=skydir, width=10.0, axes=[energy_axis, time_axis], ) # Get image by coord tuple image = m_wcs.get_image_by_coord(('500 GeV', '1 h')) # Get image by coord dict with strings image = m_wcs.get_image_by_coord({'energy': '500 GeV', 'time': '1 h'}) # Get image by coord dict with quantities image = m_wcs.get_image_by_coord({'energy': 0.5 * u.TeV, 'time': 1 * u.h})
- get_image_by_idx(idx)#
Return spatial map at the given axis pixel indices.
- Parameters
- idxtuple
Tuple of scalar indices for each non spatial dimension of the map. Tuple should be ordered as (I_0, …, I_n).
- Returns
- map_out
Map
Map with spatial dimensions only.
- map_out
See also
- get_image_by_pix(pix)#
Return spatial map at the given axis pixel coordinates
- Parameters
- pixtuple
Tuple of scalar pixel coordinates for each non-spatial dimension of the map. Tuple should be ordered as (I_0, …, I_n). Pixel coordinates can be either float or integer type.
- Returns
- map_out
Map
Map with spatial dimensions only.
- map_out
See also
- integral(axis_name, coords, **kwargs)#
Compute integral along a given axis
This method uses interpolation of the cumulative sum.
- interp_by_coord(coords, **kwargs)[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.
- coordstuple, dict or
- Returns
- vals
ndarray
Interpolated pixel values.
- vals
- interp_by_pix(pix, **kwargs)[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. 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.
- Returns
- vals
ndarray
Interpolated pixel values.
- vals
- interp_to_geom(geom, preserve_counts=False, fill_value=0, **kwargs)#
Interpolate map to input geometry.
- Parameters
- geom
Geom
Target Map geometry
- preserve_countsbool
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)
- **kwargsdict
Keyword arguments passed to
Map.interp_by_coord
- geom
- Returns
- interp_map
Map
Interpolated Map
- interp_map
- is_allclose(other, rtol_axes=0.001, atol_axes=1e-06, **kwargs)#
Compare two Maps for close equivalency
- Parameters
- other
gammapy.maps.Map
The Map to compare against
- rtol_axesfloat
Relative tolerance for the axes comparison.
- atol_axesfloat
Relative tolerance for the axes comparison.
- **kwargsdict
keywords passed to
numpy.allclose
- other
- Returns
- is_allclosebool
Whether the Map is all close.
- iter_by_axis(axis_name, keepdims=False)#
“Iterate over a given axis
- Yields
- map
Map
Map iteration.
- map
See also
iter_by_image
iterate by image returning a map
- iter_by_axis_data(axis_name)[source]#
Iterate data by axis
- Parameters
- axis_namestr
Axis name
- Returns
- idx, datatuple,
Quantity
Data and index
- idx, datatuple,
- iter_by_image(keepdims=False)#
Iterate over image planes of a map.
- Parameters
- keepdimsbool
Keep dimensions.
- Yields
- map
Map
Map iteration.
- map
See also
iter_by_image_data
iterate by image returning data and index
- iter_by_image_data()#
Iterate over image planes of the map.
The image plane index is in data order, so that the data array can be indexed directly.
- Yields
- (data, idx)tuple
Where
data
is anumpy.ndarray
view of the image plane data, andidx
is a tuple of int, the index of the image plane.
See also
iter_by_image
iterate by image returning a map
- mask_nearest_position(position)#
Given a sky coordinate return nearest valid position in the mask
If the mask contains additional axes, the mask is reduced over those.
- normalize(axis_name=None)#
Normalise data in place along a given axis.
- Parameters
- axis_namestr
Along which axis to normalize.
- pad(pad_width, axis_name=None, mode='constant', cval=0, method='linear')#
Pad the spatial dimensions of the map.
- Parameters
- pad_width{sequence, array_like, int}
Number of pixels padded to the edges of each axis.
- axis_namestr
Which axis to downsample. By default spatial axes are padded.
- mode{‘edge’, ‘constant’, ‘interp’}
Padding mode. ‘edge’ pads with the closest edge value. ‘constant’ pads with a constant value. ‘interp’ pads with an extrapolated value.
- cvalfloat
Padding value when mode=’consant’.
- Returns
- map
Map
Padded map.
- map
- plot(ax=None, axis_name=None, **kwargs)[source]#
Plot the data contained in region map along the non-spatial axis.
- Parameters
- ax
Axis
Axis used for plotting
- axis_namestr
Which axis to plot on the x axis. Extra axes will be plotted as additional lines.
- **kwargsdict
Keyword arguments passed to
errorbar
- ax
- Returns
- ax
Axis
Axis used for plotting
- ax
- plot_grid(figsize=None, ncols=3, **kwargs)#
Plot map as a grid of subplots for non-spatial axes
- Parameters
- figsizetuple of int
Figsize to plot on
- ncolsint
Number of columns to plot
- **kwargsdict
Keyword arguments passed to
Map.plot
.
- Returns
- axes
ndarray
ofAxes
Axes grid
- axes
- plot_interactive()[source]#
Plot map with interactive widgets to explore the non spatial axes.
- Parameters
- rc_paramsdict
Passed to
matplotlib.rc_context(rc=rc_params)
to style the plot.- **kwargsdict
Keyword arguments passed to
WcsNDMap.plot
.
Examples
You can try this out e.g. using a Fermi-LAT diffuse model cube with an energy axis:
from gammapy.maps import Map m = Map.read("$GAMMAPY_DATA/fermi_3fhl/gll_iem_v06_cutout.fits") m.plot_interactive(add_cbar=True, stretch="sqrt")
If you would like to adjust the figure size you can use the
rc_params
argument:rc_params = {'figure.figsize': (12, 6), 'font.size': 12} m.plot_interactive(rc_params=rc_params)
- plot_region(ax=None, **kwargs)[source]#
Plot region
- Parameters
- ax
WCSAxes
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
- ax
- classmethod read(filename, format='gadf', ogip_column=None, hdu=None)[source]#
Read from file.
- Parameters
- filename
pathlib.Path
or str Filename.
- format{“gadf”, “ogip”, “ogip-arf”}
Which format to use.
- ogip_column{None, “COUNTS”, “QUALITY”, “BACKSCAL”}
If format ‘ogip’ is chosen which table hdu column to read.
- hdustr
Name or index of the HDU with the map data.
- filename
- Returns
- region_map
RegionNDMap
Region nd map
- region_map
- reduce(axis_name, func=<ufunc 'add'>, keepdims=False, weights=None)#
Reduce map over a single non-spatial axis
- Parameters
- Returns
- map_out
Map
Map with the given non-spatial axes reduced
- map_out
- reduce_over_axes(func=<ufunc 'add'>, keepdims=False, axes_names=None, weights=None)#
Reduce map over non-spatial axes
- Parameters
- Returns
- map_out
Map
Map with non-spatial axes reduced
- map_out
- rename_axes(names, new_names)#
Rename the Map axes.
- 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
Map
Renamed Map.
- geom
- reproject_to_geom(geom, preserve_counts=False, precision_factor=10)#
Reproject map to input geometry.
- Parameters
- geom
Geom
Target Map geometry
- preserve_countsbool
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)
- precision_factorint
Minimal factor between the bin size of the output map and the oversampled base map. Used only for the oversampling method.
- geom
- Returns
- output_map
Map
Reprojected Map
- output_map
- resample(geom, weights=None, preserve_counts=True)#
Resample pixels to
geom
with givenweights
.- Parameters
- geom
Geom
Target Map geometry
- weights
ndarray
Weights vector. Default is weight of one. Must have same shape as the data of the map.
- preserve_countsbool
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)
- geom
- Returns
- resampled_map
Map
Resampled map
- resampled_map
- resample_axis(axis, weights=None, ufunc=<ufunc 'add'>)#
Resample map to a new axis by grouping and reducing smaller bins by a given ufunc
By default, the map content are summed over the smaller bins. Other numpy ufunc can be used, e.g.
numpy.logical_and
ornumpy.logical_or
.
- sample_coord(n_events, random_state=0)#
Sample position and energy of events.
- Parameters
- n_eventsint
Number of events to sample.
- random_state{int, ‘random-seed’, ‘global-rng’,
RandomState
} Defines random number generator initialisation. Passed to
get_random_state
.
- Returns
- coords
MapCoord
object. Sequence of coordinates and energies of the sampled events.
- coords
- set_by_coord(coords, vals)#
Set pixels at
coords
with givenvals
.
- set_by_idx(idx, value)[source]#
Set pixels at
idx
with givenvals
.- 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.
- vals
ndarray
Values vector.
- set_by_pix(pix, vals)#
Set pixels at
pix
with givenvals
.- Parameters
- pixtuple
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. Pixel indices can be either float or integer type. Float indices will be rounded to the nearest integer.
- vals
ndarray
Values vector.
- slice_by_idx(slices)#
Slice sub map from map object.
- split_by_axis(axis_name)#
Split a Map along an axis into multiple maps.
- Parameters
- axis_namestr
Name of the axis to split
- Returns
- mapslist
A list of
Map
- stack(other, weights=None, nan_to_num=True)[source]#
Stack other region map into map.
- Parameters
- other
RegionNDMap
Other map to stack
- weights
RegionNDMap
Array to be used as weights. The spatial geometry must be equivalent to
other
and additional axes must be broadcastable.- nan_to_num: bool
Non-finite values are replaced by zero if True (default).
- other
- sum_over_axes(axes_names=None, keepdims=True, weights=None)#
To sum map values over all non-spatial axes.
- Parameters
- keepdimsbool, optional
If this is set to true, the axes which are summed over are left in the map with a single bin
- axes_names: list of str
Names of MapAxis to reduce over. If None, all will summed over
- weights
Map
Weights to be applied. The Map should have the same geometry.
- Returns
- map_out
Map
Map with non-spatial axes summed over
- map_out
- to_cube(axes)#
Append non-spatial axes to create a higher-dimensional Map.
This will result in a Map with a new geometry with N+M dimensions where N is the number of current dimensions and M is the number of axes in the list. The data is reshaped onto the new geometry
- Parameters
- axeslist
Axes that will be appended to this Map. The axes should have only one bin
- Returns
- map
WcsNDMap
new map
- map
- to_hdulist(format='gadf', hdu='SKYMAP', hdu_bands=None, hdu_region=None)[source]#
Convert to
HDUList
.- Parameters
- format{“gadf”, “ogip”, “ogip-sherpa”, “ogip-arf”, “ogip-arf-sherpa”}
Format specification
- hdustr
Name of the HDU with the map data, used for “gadf” format.
- hdu_bandsstr
Name or index of the HDU with the BANDS table, used for “gadf” format.
- hdu_regionstr
Name or index of the HDU with the region table.
- Returns
- hdulist
HDUList
HDU list
- hdulist
- to_table(format='gadf')[source]#
Convert to
Table
.Data format specification: PHA
- Parameters
- format{“gadf”, “ogip”, “ogip-arf”, “ogip-arf-sherpa”}
Format specification
- Returns
- table
Table
Table
- table
- to_unit(unit)#
Convert map to different unit
- Parameters
- unit
Unit
or str New unit
- unit
- Returns
- map
Map
Map with new unit and converted data
- map
- upsample(factor, order=0, preserve_counts=True, axis_name=None)[source]#
Upsample the non-spatial dimension by a given factor.
By default the first axes is upsampled.
- Parameters
- factorint
Upsampling factor.
- orderint
Order of the interpolation used for upsampling.
- preserve_countsbool
Preserve the integral over each bin. This should be true if the RegionNDMap is an integral quantity (e.g. counts) and false if the RegionNDMap is a differential quantity (e.g. intensity).
- axis_namestr
Which axis to upsample. Default is “energy”.
- Returns
- map
RegionNDMap
Upsampled region map.
- map
- write(filename, overwrite=False, format='gadf', hdu='SKYMAP')[source]#
Write map to file
- Parameters
- filename
pathlib.Path
or str Filename.
- format{“gadf”, “ogip”, “ogip-sherpa”, “ogip-arf”, “ogip-arf-sherpa”}
Which format to use.
- overwritebool
Overwrite existing files?
- filename