WCS

class mpdaf.obj.WCS(hdr=None, crpix=None, crval=(1.0, 1.0), cdelt=(1.0, 1.0), deg=False, rot=0, shape=None, cd=None, frame=None, equinox=None)[source]

Bases: object

The WCS class manages the world coordinates of the spatial axes of MPDAF images, using the pywcs package.

Note that MPDAF images are stored in python arrays that are indexed in [y,x] axis order. In general the axes of these arrays are not along celestial axes such as right-ascension and declination. They are cartesian axes of a flat map projection of the sky around the observation center, and they may be rotated away from the celestial axes. When their rotation angle is zero, the Y axis is parallel to the declination axis. However the X axis is only along the right ascension axis for observations at zero declination.

Pixels in MPDAF images are not generally square on the sky. To scale index offsets in the image to angular distances in the map projection, the Y-axis and X-axis index offsets must be scaled by different numbers. These numbers can be obtained by calling the get_axis_increments() method, which returns the angular increment per pixel increment along the Y and X axes of the array. The Y-axis increment is always positive, but the X-axis increment is negative if east is anti-clockwise of north when the X-axis pixels are plotted from left to right,

The rotation angle of the map projection, relative to the sky, can be obtained by calling the get_rot() method. This returns the angle between celestial north and the Y axis of the image, in the sense of an eastward rotation of celestial north from the Y-axis.

When the linearized coordinates of the map projection are insufficient, the celestial coordinates of one or more pixels can be queried by calling the pix2sky() method, which returns coordinates in the [dec,ra] axis order. In the other direction, the [y,x] indexes of the pixel of a given celestial coordinate can be obtained by calling the sky2pix() method.

Parameters:
hdrastropy.fits.CardList

A FITS header. If the hdr parameter is not None, the WCS object is created from the data header, and the remaining parameters are ignored.

crpixfloat or (float,float)

The FITS array indexes of the reference pixel of the image, given in the order (y,x). Note that the first pixel of the FITS image is [1,1], whereas in the python image array it is [0,0]. Thus to place the reference pixel at [ry,rx] in the python image array would require crpix=(ry+1,rx+1).

If both crpix and shape are None, then crpix is given the value (1.0,1.0) and the reference position is at index [0,0] in the python array of the image.

If crpix is None and shape is not None, then crpix is set to (shape + 1.0)/2.0, which places the reference point at the center of the image.

crvalfloat or (float,float)

The celestial coordinates of the reference pixel (ref_dec,ref_ra). If this parameter is not provided, then (0.0,0.0) is substituted.

cdeltfloat or (float,float)

If the hdr and cd parameters are both None, then this argument can be used to specify the pixel increments along the Y and X axes of the image, respectively. If this parameter is not provided, (1.0,1.0) is substituted. Note that it is conventional for cdelt[1] to be negative, such that east is plotted towards the left when the image rotation angle is zero.

degbool

If True, then cdelt and crval are in decimal degrees (CTYPE1=’RA—TAN’,CTYPE2=’DEC–TAN’,CUNIT1=CUNIT2=’deg’). If False (the default), the celestial coordinates are linear (CTYPE1=CTYPE2=’LINEAR’).

rotfloat

If the hdr and cd parameters are both None, then this argument can be used to specify a value for the rotation angle of the image. This is the angle between celestial north and the Y axis of the image, in the sense of an eastward rotation of celestial north from the Y-axis.

Along with the cdelt parameter, the rot parameter is used to construct a FITS CD rotation matrix. This is done as described in equation 189 of Calabretta, M. R., and Greisen, E. W, Astronomy & Astrophysics, 395, 1077-1122, 2002, where it serves as the value of the CROTA term.

shapeinteger or (integer,integer)

The dimensions of the image axes (optional). The dimensions are given in python order (ny,nx).

cdnumpy.ndarray

This parameter can optionally be used to specify the FITS CD rotation matrix. By default this parameter is None. However if a matrix is provided and hdr is None, then it is used instead of cdelt and rot, which are then ignored. The matrix should be ordered like

cd = numpy.array([[CD1_1, CD1_2],

[CD2_1, CD2_2]]),

where CDj_i are the names of the corresponding FITS keywords.

Attributes:
wcsastropy.wcs.WCS

The underlying object that performs most of the world coordinate conversions.

Attributes Summary

naxis1

naxis2

unit

Return the default angular unit used for sky coordinates.

Methods Summary

coord([spaxel, relative, center, polar, ...])

Return the full coordinate array of the image

copy()

Return a copy of a WCS object.

get_axis_increments([unit])

Return the displacements on the sky that result from incrementing the array indexes of the image by one along the Y and X axes, respectively.

get_cd()

Return the coordinate conversion matrix (CD).

get_center([unit])

Return the center (dec,ra) of the image array.

get_crpix1()

Return the value of the FITS CRPIX1 parameter.

get_crpix2()

Return the value of the FITS CRPIX2 parameter.

get_crval1([unit])

Return the value of the FITS CRVAL1 parameter.

get_crval2([unit])

Return the value of the FITS CRVAL2 parameter.

get_end([unit])

Return the [dec,ra] coordinates of pixel (-1,-1).

get_range([unit])

Return the minimum and maximum right-ascensions and declinations in the image array.

get_rot([unit])

Return the rotation angle of the image.

get_start([unit])

Return the [dec,ra] coordinates of pixel (0,0).

get_step([unit])

Return the angular height and width of a pixel along the Y and X axes of the image array.

info()

Print information about a WCS object.

isEqual(other[, start_atol, rot_atol])

Return True if other and self have the same attributes.

is_deg()

Return True if world coordinates are in decimal degrees.

pix2sky(x[, unit])

Convert image pixel indexes (y,x) to world coordinates (dec,ra).

rebin(factor)

Rebin to a new coordinate system.

rotate(theta)

Rotate WCS coordinates to new orientation given by theta.

sameStep(other)

Return True if other and self have the same pixel sizes.

set_axis_increments(increments[, unit])

Set the displacements on the sky that result from incrementing the array indexes of the image by one along the Y and X axes, respectively.

set_cd(cd)

Install a new coordinate transform matrix.

set_crpix1(x)

Set the value of the FITS CRPIX1 parameter.

set_crpix2(y)

Set the value of the FITS CRPIX2 parameter.

set_crval1(x[, unit])

Set the value of the CRVAL1 keyword.

set_crval2(x[, unit])

Set the value of the CRVAL2 keyword.

set_step(step[, unit])

Set the height and width of pixels on the sky.

sky2pix(x[, nearest, unit])

Convert world coordinates (dec,ra) to image pixel indexes (y,x).

to_cube_header(wave)

Generate an astropy.fits header object with WCS information and wavelength information.

to_header()

Generate an astropy.fits header object containing the WCS information.

Attributes Documentation

naxis1
naxis2
unit

Return the default angular unit used for sky coordinates.

Returns:
outastropy.units.Unit

The unit to use for coordinate angles.

Methods Documentation

coord(spaxel=False, relative=False, center=None, polar=False, unit=None, reshape=False, mask=None)[source]

Return the full coordinate array of the image

Parameters:
unitastropy.units.Unit

Unit of the spatial coordinates if None return pixels

relativebool

If True, the coordinate are given relative to the center

centertuple

center coordinate use for relative and radial mode if none, the image center is used must be given in unit

polarbool

if True, return polar coordinate (r, theta)

unitastropy.units.Unit

The units of the world coordinates

mask2D array of bool

If not None, return only non masked spaxels

reshapebool

if True, return two 1D array, if False, return two 2D array (meshgrid)

Returns:
outarray of array of float

array of [dec,ra], [p,q], [delta_dec, delta_ra], [delta_p, delta_q], [r, theta]

copy()[source]

Return a copy of a WCS object.

get_axis_increments(unit=None)[source]

Return the displacements on the sky that result from incrementing the array indexes of the image by one along the Y and X axes, respectively.

In MPDAF, images are sampled on a regular grid of square pixels that represent a flat projection of the celestial sphere. The get_axis_increments() method returns the angular width and height of these pixels on the sky, with signs that indicate whether the angle increases or decreases as one increments the array indexes. To keep plots consistent, regardless of the rotation angle of the image on the sky, the returned height is always positive, but the returned width is negative if a plot of the image with pixel 0,0 at the bottom left would place east anticlockwise of north, and positive otherwise.

Parameters:
unitastropy.units.Unit

The angular units of the returned values.

Returns:
outnumpy.ndarray

(dy,dx). These are the angular increments of pixels along the Y and X axes of the image. The returned values are either in the unit specified by the ‘unit’ input parameter, or in the unit specified by the self.unit property.

get_cd()[source]

Return the coordinate conversion matrix (CD).

This is a 2x2 matrix that can be used to convert from the column and row indexes of a pixel in the image array to a coordinate within a flat map-projection of the celestial sphere. For example, if the celestial coordinates of the observation are right-ascension and declination, and r and d denote their gnonomic TAN projection onto a flat plane, then a pixel at row and column [col,row] has [r,d] coordinates given by:

(r,d) = np.dot(get_cd(), (col - get_crpix1(), row - get_crpix2())
Returns:
outnump.ndarray

A 2D array containing the coordinate transformation matrix, arranged such that the elements described in the FITS standard are arranged as follows:

[[CD_11, CD_12]
 [CD_21, CD_22]]
get_center(unit=None)[source]

Return the center (dec,ra) of the image array.

Parameters:
unitastropy.units.Unit

The units of the returned angles.

Returns:
centernumpy.ndarray

The center, arranged as [dec, ra]. The returned values are either in the units specified in the ‘unit’ input parameter, or in the units stored in the self.unit property.

get_crpix1()[source]

Return the value of the FITS CRPIX1 parameter.

CRPIX1 contains the index of the reference position of the image along the X-axis of the image. Beware that this is a FITS array index, which is 1 greater than the corresponding python array index. For example, a crpix value of 1 denotes a python array index of 0. The reference pixel index is a floating point value that can indicate a position between two pixels. It can also indicate an index that is outside the bounds of the array.

Returns:
outfloat

The value of the FITS CRPIX1 parameter.

get_crpix2()[source]

Return the value of the FITS CRPIX2 parameter.

CRPIX2 contains the index of the reference position of the image along the Y-axis of the image. Beware that this is a FITS array index, which is 1 greater than the corresponding python array index. For example, a crpix value of 1 denotes a python array index of 0. The reference pixel index is a floating point value that can indicate a position between two pixels. It can also indicate an index that is outside the bounds of the array.

Returns:
outfloat

The value of the FITS CRPIX2 parameter.

get_crval1(unit=None)[source]

Return the value of the FITS CRVAL1 parameter.

CRVAL1 contains the coordinate reference value of the first image axis (eg. right-ascension).

Parameters:
unitastropy.units.Unit

The angular units to give the return value.

Returns:
outfloat

The value of CRVAL1 in the specified angular units. If the units are not given, then the unit in the self.unit property is used.

get_crval2(unit=None)[source]

Return the value of the FITS CRVAL2 parameter.

CRVAL2 contains the coordinate reference value of the second image axis (eg. declination).

Parameters:
unitastropy.units.Unit

The angular units to give the return value.

Returns:
outfloat

The value of CRVAL2 in the specified angular units. If the units are not given, then the unit in the self.unit property is used.

get_end(unit=None)[source]

Return the [dec,ra] coordinates of pixel (-1,-1).

Parameters:
unitastropy.units.Unit

The angular units of the returned coordinates.

Returns:
outnumpy.ndarray

The equatorial coordinate of pixel [-1,-1], ordered as, [dec,ra]. If a value was given to the optional ‘unit’ argument, the angular unit specified there will be used for the return value. Otherwise the unit stored in the self.unit property will be used.

get_range(unit=None)[source]

Return the minimum and maximum right-ascensions and declinations in the image array.

Specifically a list is returned with the following contents:

[dec_min, ra_min, dec_max, ra_max]

Note that if the Y axis of the image is not parallel to the declination axis, then the 4 returned values will all come from different corners of the image. In particular, note that this means that the coordinates [dec_min,ra_min] and [dec_max,ra_max] will only coincide with pixels in the image if the Y axis is aligned with the declination axis. Otherwise they will be outside the bounds of the image.

Parameters:
unitastropy.units.Unit

The units of the returned angles.

Returns:
outnumpy.ndarray

The range of right ascensions and declinations, arranged as [dec_min, ra_min, dec_max, ra_max]. The returned values are either in the units specified in the ‘unit’ input parameter, or in the units stored in the self.unit property.

get_rot(unit=Unit('deg'))[source]

Return the rotation angle of the image.

The angle is defined such that a rotation angle of zero aligns north along the positive Y axis, and a positive rotation angle rotates north away from the Y axis, in the sense of a rotation from north to east.

Note that the rotation angle is defined in a flat map-projection of the sky. It is what would be seen if the pixels of the image were drawn with their pixel widths scaled by the angular pixel increments returned by the get_axis_increments() method.

Parameters:
unitastropy.units.Unit

The unit to give the returned angle (degrees by default).

Returns:
outfloat

The angle between celestial north and the Y axis of the image, in the sense of an eastward rotation of celestial north from the Y-axis.

get_start(unit=None)[source]

Return the [dec,ra] coordinates of pixel (0,0).

Parameters:
unitastropy.units.Unit

The angular units of the returned coordinates.

Returns:
outnumpy.ndarray

The equatorial coordinate of pixel [0,0], ordered as: [dec,ra]. If a value was given to the optional ‘unit’ argument, the angular unit specified there will be used for the return value. Otherwise the unit stored in the self.unit property will be used.

get_step(unit=None)[source]

Return the angular height and width of a pixel along the Y and X axes of the image array.

In MPDAF, images are sampled on a regular grid of square pixels that represent a flat projection of the celestial sphere. The get_step() method returns the angular width and height of these pixels on the sky.

See also get_axis_increments().

Parameters:
unitastropy.units.Unit

The angular units of the returned values.

Returns:
outnumpy.ndarray

(dy,dx). These are the angular height and width of pixels along the Y and X axes of the image. The returned values are either in the unit specified by the ‘unit’ input parameter, or in the unit specified by the self.unit property.

info()[source]

Print information about a WCS object.

isEqual(other, start_atol=1e-06, rot_atol=1e-06)[source]

Return True if other and self have the same attributes.

Beware that if the two wcs objects have the same world coordinate characteristics, but come from images of different dimensions, the objects will be considered different.

Parameters:
otherWCS

The wcs object to be compared to self.

start_atolfloat

Absolute tolerance for the test with WCS.get_start, in degrees.

start_rotfloat

Absolute tolerance for the test with WCS.get_rot, in degrees.

Returns:
outbool

True if the two WCS objects have the same attributes.

is_deg()[source]

Return True if world coordinates are in decimal degrees.

(CTYPE1=’RA—TAN’,CTYPE2=’DEC–TAN’,CUNIT1=CUNIT2=’deg)

pix2sky(x, unit=None)[source]

Convert image pixel indexes (y,x) to world coordinates (dec,ra).

Parameters:
xarray

An (n,2) array of image pixel indexes. These should be python array indexes, ordered like (y,x) and with 0,0 denoting the lower left pixel of the image.

unitastropy.units.Unit

The units of the world coordinates.

Returns:
out(n,2) array of dec- and ra- world coordinates.
rebin(factor)[source]

Rebin to a new coordinate system.

This is a helper function for the Image and Cube rebin() methods.

Parameters:
factor(integer,integer)

Factor in y and x.

Returns:
outWCS
rotate(theta)[source]

Rotate WCS coordinates to new orientation given by theta.

Analog to astropy.wcs.WCS.rotateCD, which is deprecated since version 1.3 (see https://github.com/astropy/astropy/issues/5175).

Parameters:
thetafloat

Rotation in degree.

sameStep(other)[source]

Return True if other and self have the same pixel sizes.

Parameters:
otherWCS

The wcs object to compare to self.

Returns:
outbool

True if the two arrays of axis step increments are equal.

set_axis_increments(increments, unit=None)[source]

Set the displacements on the sky that result from incrementing the array indexes of the image by one along the Y and X axes, respectively.

In MPDAF, images are sampled on a regular grid of square pixels that represent a flat projection of the celestial sphere. The set_axis_increments() method changes the angular width and height of these pixels on the sky, with signs that indicate whether the angle increases or decreases as one increments the array indexes. To keep plots consistent, regardless of the rotation angle of the image on the sky, the height should always be positive, and the width should be negative if a plot of the image with pixel 0,0 at the bottom left would place east anticlockwise of north, and positive otherwise.

Parameters:
incrementsnumpy.ndarray

(dy,dx). These are the desired angular increments of pixels along the Y and X axes of the image. These should either be in the unit specified by the ‘unit’ input parameter, or, if unit=None, in the unit specified by the self.unit property.

unitastropy.units.Unit

The angular units of the specified increments.

set_cd(cd)[source]

Install a new coordinate transform matrix.

This is a 2x2 matrix that is used to convert from the row and column indexes of a pixel in the image array to a coordinate within a flat map-projection of the celestial sphere. It is formerly described in the FITS standard. The matrix should be ordered like:

cd = numpy.array([[CD1_1, CD1_2],
                  [CD2_1, CD2_2]]),

where CDj_i are the names of the corresponding FITS keywords.

Parameters:
cdnumpy.ndarray

The 2x2 coordinate conversion matrix, with its elements ordered for multiplying a column vector in FITS (x,y) axis order.

set_crpix1(x)[source]

Set the value of the FITS CRPIX1 parameter.

This sets the reference pixel index along the X-axis of the image.

This is a floating point value which can denote a position between pixels. It is specified with the FITS indexing convention, where FITS pixel 1 is equivalent to pixel 0 in python arrays. In general subtract 1 from x to get the corresponding floating-point pixel index along axis 1 of the image array. In cases where x is an integer, the corresponding row in the python data array that contains the image is data[:, x-1].

Parameters:
xfloat

The index of the reference pixel along the X axis.

set_crpix2(y)[source]

Set the value of the FITS CRPIX2 parameter.

This sets the reference pixel index along the Y-axis of the image.

This is a floating point value which can denote a position between pixels. It is specified with the FITS indexing convention, where FITS pixel 1 is equivalent to pixel 0 in python arrays. In general subtract 1 from y to get the corresponding floating-point pixel index along axis 0 of the image array. In cases where y is an integer, the corresponding column in the python data array that contains the image is data[y-1, :].

Parameters:
yfloat

The index of the reference pixel along the Y axis.

set_crval1(x, unit=None)[source]

Set the value of the CRVAL1 keyword.

It indicates the coordinate reference value along the first image axis (eg. right ascension).

Parameters:
xfloat

The value of the reference pixel on the first axis.

unitastropy.units.Unit

The angular units of the world coordinates.

set_crval2(x, unit=None)[source]

Set the value of the CRVAL2 keyword.

It indicates the coordinate reference value along the second image axis (eg. declination).

Parameters:
xfloat

The value of the reference pixel on the second axis.

unitastropy.units.Unit

The angular units of the world coordinates.

set_step(step, unit=None)[source]

Set the height and width of pixels on the sky.

In MPDAF, images are sampled on a regular grid of square pixels that represent a flat projection of the celestial sphere. The set_step() method changes the angular width and height of these pixels on the sky.

Parameters:
steparray_like

(h,w). These are the desired angular height and width of pixels along the Y and X axes of the image. These should either be in the unit specified by the ‘unit’ input parameter, or, if unit=None, in the unit specified by the self.unit property.

unitastropy.units.Unit

The angular units of the specified increments.

sky2pix(x, nearest=False, unit=None)[source]

Convert world coordinates (dec,ra) to image pixel indexes (y,x).

If nearest=True; returns the nearest integer pixel.

Parameters:
xarray

An (n,2) array of dec- and ra- world coordinates.

nearestbool

If nearest is True returns the nearest integer pixel in place of the decimal pixel.

unitastropy.units.Unit

The units of the world coordinates

Returns:
out(n,2) array of image pixel indexes. These are

python array indexes, ordered like (y,x) and with 0,0 denoting the lower left pixel of the image.

to_cube_header(wave)[source]

Generate an astropy.fits header object with WCS information and wavelength information.

to_header()[source]

Generate an astropy.fits header object containing the WCS information.