NDCubeABC#
- class ndcube.ndcube.NDCubeABC[source]#
Bases:
NDDataBase
Attributes Summary
Returns the WCS physical types that vary along each array axis.
The WCS transform for the NDCube, including the coordinates specified in
.extra_coords
.Coordinates not described by
NDCubeABC.wcs
which vary along one or more axes.Coordinate metadata which applies to the whole cube.
Methods Summary
axis_world_coords
(*axes[, pixel_corners, wcs])Returns objects representing the world coordinates of pixel centers for a desired axes.
axis_world_coords_values
(*axes[, ...])Returns the world coordinate values of all pixels for desired axes.
crop
(*points[, wcs, keepdims])Crop using real world coordinates.
crop_by_values
(*points[, units, wcs, keepdims])Crop using real world coordinates.
Attributes Documentation
- array_axis_physical_types#
Returns the WCS physical types that vary along each array axis.
Returns an iterable of tuples where each tuple corresponds to an array axis and holds strings denoting the WCS physical types associated with that array axis. Since multiple physical types can be associated with one array axis, tuples can be of different lengths. Likewise, as a single physical type can correspond to multiple array axes, the same physical type string can appear in multiple tuples.
The physical types returned by this property are drawn from the
combined_wcs
property so they include the coordinates contained inextra_coords
.
- combined_wcs#
The WCS transform for the NDCube, including the coordinates specified in
.extra_coords
.This transform should implement the high level wcsapi, and have
pixel_n_dim
equal to the number of array dimensions in theNDCube
. The number of world dimensions should be equal to the number of world dimensions inself.wcs
and inself.extra_coords
combined.
- extra_coords#
Coordinates not described by
NDCubeABC.wcs
which vary along one or more axes.
- global_coords#
Coordinate metadata which applies to the whole cube.
Methods Documentation
- abstract axis_world_coords(*axes: int | str, pixel_corners: bool = False, wcs: BaseHighLevelWCS | ExtraCoordsABC | None = None) Iterable[Any] [source]#
Returns objects representing the world coordinates of pixel centers for a desired axes.
- Parameters:
axes (
int
orstr
, or multipleint
orstr
, optional) – Axis number in numpy ordering or unique substring ofndcube.NDCube.wcs.world_axis_physical_types
of axes for which real world coordinates are desired. Not specifying axes inputs causes results for all axes to be returned.pixel_corners (
bool
, optional) – IfTrue
then instead of returning the coordinates at the centers of the pixels, the coordinates at the pixel corners will be returned. This increases the size of the output by 1 in all dimensions as all corners are returned.wcs (
astropy.wcs.wcsapi.BaseHighLevelWCS
, optional) – The WCS object to used to calculate the world coordinates. Although technically this can be any valid WCS, it will typically beself.wcs
,self.extra_coords
, orself.combined_wcs
combining both the WCS and extra coords. Default=self.wcs
- Returns:
axes_coords (iterable) – An iterable of “high level” objects giving the real world coords for the axes requested by user. For example, a tuple of
SkyCoord
objects. The types returned are determined by the WCS object. The dimensionality of these objects should match that of their corresponding array dimensions, unlesspixel_corners=True
in which case the length along each axis will be 1 greater than the number of pixels.
Examples
>>> NDCube.axis_world_coords('lat', 'lon') >>> NDCube.axis_world_coords(2)
- abstract axis_world_coords_values(*axes: int | str, pixel_corners: bool = False, wcs: BaseHighLevelWCS | ExtraCoordsABC | None = None) Iterable[Quantity] [source]#
Returns the world coordinate values of all pixels for desired axes. In contrast to
ndcube.NDCube.axis_world_coords()
, this method returnsQuantity
objects. Which only provide units rather than full coordinate metadata provided by high-level coordinate objects.- Parameters:
axes (
int
orstr
, or multipleint
orstr
, optional) – Axis number in numpy ordering or unique substring ofndcube.NDCube.wcs.world_axis_physical_types
of axes for which real world coordinates are desired. axes=None implies all axes will be returned.pixel_corners (
bool
, optional) – IfTrue
then coordinates at pixel corners will be returned rather than at pixel centers. This increases the size of the output along each dimension by 1 as all corners are returned.wcs (
BaseHighLevelWCS
orExtraCoordsABC
, optional) – The WCS object to be used to calculate the world coordinates. Although technically this can be any valid WCS, it will typically beself.wcs
,self.extra_coords
, orself.combined_wcs
, combing both the WCS and extra coords. Defaults to the.wcs
property.
- Returns:
axes_coords (
tuple
ofQuantity
) – An iterable of raw coordinate values for all pixels for the requested axes. The returned units are determined by the WCS object. The dimensionality of these objects should match that of their corresponding array dimensions, unlesspixel_corners=True
in which case the length along each axis will be 1 greater than the number of pixels.
Examples
>>> NDCube.axis_world_coords_values('lat', 'lon') >>> NDCube.axis_world_coords_values(2)
- abstract crop(*points: Iterable[Any], wcs: BaseHighLevelWCS | ExtraCoordsABC | None = None, keepdims: bool = False) NDCubeABC [source]#
Crop using real world coordinates. This method crops the NDCube to the smallest bounding box in pixel space that contains all the provided world coordinate points.
This function takes the points defined as high-level astropy coordinate objects such as
SkyCoord
,SpectralCoord
, etc.- Parameters:
points (iterable of iterables) – Tuples of high level coordinate objects e.g.
SkyCoord
. Each iterable of coordinate objects represents a single location in the data array in real world coordinates.The coordinates of the points as they are passed to
world_to_array_index
. Therefore their number and order must be compatible with the API of that method, i.e. they must be passed in world order.wcs (
BaseHighLevelWCS
orExtraCoordsABC
) – The WCS to use to calculate the pixel coordinates based on the input. Will default to the.wcs
property if not given. While any valid WCS could be used it is expected that either the.wcs
or.extra_coords
properties will be used.keepdims (
bool
, optional) – IfFalse
and if cropping results in length-1 dimensions, these are sliced away in output cube. IfTrue
, length-1 dimensions are kept. Default=False
- Returns:
Examples
>>> # An example of cropping a region of interest on the Sun from a 3-D image-time cube >>> point1 = [SkyCoord(-50*u.deg, -40*u.deg, frame=frames.HeliographicStonyhurst), None] >>> point2 = [SkyCoord(0*u.deg, -6*u.deg, frame=frames.HeliographicStonyhurst), None] >>> NDCube.crop(point1, point2)
- abstract crop_by_values(*points: Iterable[Quantity | float], units: Iterable[str | Unit] | None = None, wcs: BaseHighLevelWCS | ExtraCoordsABC | None = None, keepdims: bool = False) NDCubeABC [source]#
Crop using real world coordinates. This method crops the NDCube to the smallest bounding box in pixel space that contains all the provided world coordinate points.
This function takes points as iterables of low-level coordinate objects, i.e.
Quantity
objects. This differs fromcrop()
which takes high-level coordinate objects requiring all the relevant coordinate information such as coordinate frame etc. Hence this method’s API is more basic but less explicit.- Parameters:
points (iterable) – Tuples of coordinate values, the length of the tuples must be equal to the number of world dimensions. These points are passed to
wcs.world_to_array_index_values
so their units and order must be compatible with that method.units (
str
orUnit
) – If the inputs are set without units, the user must set the units inside this argument asstr
orUnit
objects. The length of the iterable must equal the number of world dimensions and must have the same order as the coordinate points.wcs (
BaseHighLevelWCS
orExtraCoordsABC
) – The WCS to use to calculate the pixel coordinates based on the input. Will default to the.wcs
property if not given. While any valid WCS could be used it is expected that either the.wcs
or.extra_coords
properties will be used.keepdims (
bool
, optional) – IfFalse
and if cropping results in length-1 dimensions, these are sliced away in output cube. IfTrue
, length-1 dimensions are kept. Default=False
- Returns:
Examples
>>> # An example of cropping a region of interest on the Sun from a 3-D image-time cube >>> NDCube.crop_by_values((-600, -600, 0), (0, 0, 0), units=(u.arcsec, u.arcsec, u.s))