Under the hood, there are 3 separate classes that perform different parts of the transformation:
Additionally, the class WCS aggregates all of these transformations together in a pipeline:
- Detector to image plane correction (by a pair of DistortionLookupTable objects).
- SIP distortion correction (by an underlying Sip object)
- Paper IV table-lookup distortion correction (by a pair of DistortionLookupTable objects).
- wcslib WCS transformation (by a Wcsprm object)
WCS objects perform standard WCS transformations, and correct for SIP and Paper IV table-lookup distortions, based on the WCS keywords and supplementary data read from a FITS file.
header: A PyFITS header object. If header is not provided, the object will be initialized to default values.
fobj: A PyFITS file (hdulist) object. It is needed when header keywords point to a Paper IV Lookup table distortion stored in a different extension.
key: A string. The name of a particular WCS transform to use. This may be either ' ' or 'A'-'Z' and corresponds to the "a" part of the CTYPEia cards. key may only be provided if header is also provided.
minerr: A floating-point value. The minimum value a distortion correction must have in order to be applied. If the value of CQERRja is smaller than minerr, the corresponding distortion is not applied.
relax: Degree of permissiveness:
- False: Recognize only FITS keywords defined by the published WCS standard.
- True: Admit all recognized informal extensions of the WCS standard.
naxis: int or sequence. Extracts specific coordinate axes using sub. If a header is provided, and naxis is not None, naxis will be passed to sub in order to select specific axes from the header. See sub for more details about this parameter.
Exceptions:
Transforms pixel coordinates to sky coordinates by doing all of the following in order:
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the pixel coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
Calculates the footprint of the image on the sky.
A footprint is defined as the positions of the corners of the image on the sky after all available distortions have been applied.
Returns a (4, 2) array of (x, y) coordinates.
Return a shallow copy of the object.
Convenience method so user doesn’t have to import the copy stdlib module.
The pre-linear transformation distortion lookup table, CPDIS1.
The pre-linear transformation distortion lookup table, CPDIS2.
Return a deep copy of the object.
Convenience method so user doesn’t have to import the copy stdlib module.
Convert detector coordinates to image plane coordinates using Paper IV table-lookup distortion correction.
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the pixel coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
Writes out a ds9 style regions file. It can be loaded directly by ds9.
Convert pixel coordinates to focal plane coordinates using Paper IV table-lookup distortion correction.
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the pixel coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
Convert pixel coordinates to focal plane coordinates using the SIP polynomial distortion convention and Paper IV table-lookup distortion correction.
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the pixel coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
Convert focal plane coordinates to pixel coordinates using the SIP polynomial distortion convention.
Paper IV table lookup distortion correction is not applied, even if that information existed in the FITS file that initialized this WCS object.
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the focal plane coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
Convert pixel coordinates to focal plane coordinates using the SIP polynomial distortion convention.
Paper IV table lookup distortion correction is not applied, even if that information existed in the FITS file that initialized this WCS object. To correct for that, use pix2foc or p4_pix2foc.
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the pixel coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
sub(axes) -> WCS object
Extracts the coordinate description for a subimage from a WCS object.
The world coordinate system of the subimage must be separable in the sense that the world coordinates at any point in the subimage must depend only on the pixel coordinates of the axes extracted. In practice, this means that the PCi_ja matrix of the original image must not contain non-zero off-diagonal terms that associate any of the subimage axes with any of the non-subimage axes.
Coordinate axes types may be specified using either strings or special integer constants. The available types are:
- 'longitude' / WCSSUB_LONGITUDE: Celestial longitude
- 'latitude' / WCSSUB_LATITUDE: Celestial latitude
- 'cubeface' / WCSSUB_CUBEFACE: Quadcube CUBEFACE axis
- 'spectral' / WCSSUB_SPECTRAL: Spectral axis
- 'stokes' / WCSSUB_STOKES: Stokes axis
Returns a WCS object.
Exceptions:
Generate a pyfits header object with the WCS information stored in this object.
Warning
This function does not write out SIP or Paper IV distortion keywords, yet, only the core WCS support by wcslib.
The output header will almost certainly differ from the input in a number of respects:
- The output header only contains WCS-related keywords. In particular, it does not contain syntactically-required keywords such as SIMPLE, NAXIS, BITPIX, or END.
- Deprecated (e.g. CROTAn) or non-standard usage will be translated to standard (this is partially dependent on whether fix was applied).
- Quantities will be converted to the units used internally, basically SI with the addition of degrees.
- Floating-point quantities may be given to a different decimal precision.
- Elements of the PCi_j matrix will be written if and only if they differ from the unit matrix. Thus, if the matrix is unity then no elements will be written.
- Additional keywords such as WCSAXES, CUNITia, LONPOLEa and LATPOLEa may appear.
- The original keycomments will be lost, although to_header tries hard to write meaningful comments.
- Keyword order may be changed.
Returns a pyfits Header object.
Transforms pixel coordinates to sky coordinates by doing only the basic wcslib transformation. No SIP or Paper IV table lookup distortion correction is applied. To perform distortion correction, see all_pix2sky, sip_pix2foc, p4_pix2foc, or pix2foc.
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the sky coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
Transforms sky coordinates to pixel coordinates, using only the basic wcslib WCS transformation. No SIP or Paper IV table lookup distortion is applied.
Either two or three arguments may be provided.
- 2 arguments: An N x 2 array of x- and y-coordinates, and an origin.
- 3 arguments: 2 one-dimensional arrays of x and y coordinates, and an origin.
Here, origin is the coordinate in the upper left corner of the image. In FITS and Fortran standards, this is 1. In Numpy and C standards this is 0.
Returns the pixel coordinates. If the input was a single array and origin, a single array is returned, otherwise a tuple of arrays is returned.
Exceptions:
str
Character code for alternate coordinate descriptions. For example, the "a" in keyword names such as CTYPEia). This is a space character for the primary coordinate description, or one of the 26 upper-case letters, A-Z.
double array[2][2]
The CDi_ja linear transformation matrix.
For historical compatibility, two alternate specifications of the CDi_ja and CROTAia keywords are supported. Although these may not formally co-exist with PCi_ja, the approach here is simply to ignore them if given in conjunction with PCi_ja.
has_pci_ja, has_cdi_ja and has_crotaia can be used to determine which of these alternatives are present in the header.
CDi_ja and CROTAia keywords, if found, are to be stored in the cd and crota arrays which are dimensioned similarly to pc and cdelt.
These alternate specifications of the linear transformation matrix are translated immediately to PCi_ja by set and are nowhere visible to the lower-level routines. In particular, set resets cdelt to unity if CDi_ja is present (and no PCi_ja). If no CROTAia is associated with the latitude axis, set reverts to a unity PCi_ja matrix.
double array[naxis]
Coordinate increments (CDELTia) for each coord axis.
An undefined value is represented by NaN.
boolean
If True, an offset will be applied to (x, y) to force (x,y) = (0,0) at the fiducial point.
Translates AIPS-convention celestial projection types, -NCP and -GLS.
Returns 0 for success; -1 if no change required.
list of strings
A list of the coordinate axis names, from CNAMEia.
int array[naxis]
An array recording the column numbers for each axis in a pixel list.
int
Where the coordinate representation is associated with an image-array column in a FITS binary table, this property may be used to record the relevant column number.
It should be set to zero for an image header or pixel list.
double array[naxis]
The random error in each coordinate axis, CRDERia.
An undefined value is represented by NaN.
double array[naxis]
CROTAia keyvalues for each coordinate axis.
CROTAia is an alternate specification of the linear transformation matrix, maintained for historical compatibility. See cd for the current method.
double array[naxis]
Coordinate reference pixels (CRPIXja) for each pixel axis.
double array[naxis]
Coordinate reference values (CRVALia) for each coordinate axis.
double array[naxis]
The systematic error in the coordinate value axes, CSYERia.
An undefined value is represented by NaN.
list of strings[naxis]
List of CTYPEia keyvalues.
The ctype keyword values must be in upper case and there must be zero or one pair of matched celestial axis types, and zero or one spectral axis.
int
Index into the pixcrd (pixel coordinate) array for the CUBEFACE axis. This is used for quadcube projections where the cube faces are stored on a separate axis.
The quadcube projections (TSC, CSC, QSC) may be represented in FITS in either of two ways:
The six faces may be laid out in one plane and numbered as follows:
0 4 3 2 1 4 3 2 5Faces 2, 3 and 4 may appear on one side or the other (or both). The sky-to-pixel routines map faces 2, 3 and 4 to the left but the pixel-to-sky routines accept them on either side.
The COBE convention in which the six faces are stored in a three-dimensional structure using a CUBEFACE axis indexed from 0 to 5 as above.
These routines support both methods; set determines which is being used by the presence or absence of a CUBEFACE axis in ctype. p2s and s2p translate the CUBEFACE axis representation to the single plane representation understood by the lower-level projection routines.
list of strings[naxis]
List of CUNITia keyvalues which define the units of measurement of the CRVALia, CDELTia and CDi_ja keywords.
As CUNITia is an optional header keyword, cunit may be left blank but otherwise is expected to contain a standard units specification as defined by WCS Paper I. Utility function wcsutrn, (not currently wrapped for Python) is available to translate commonly used non-standard units specifications but this must be done as a separate step before invoking set.
For celestial axes, if cunit is not blank, set uses wcsunits to parse it and scale cdelt, crval, and cd to decimal degrees. It then resets cunit to "deg".
For spectral axes, if cunit is not blank, set uses wcsunits to parse it and scale cdelt, crval, and cd to SI units. It then resets cunit accordingly.
set ignores cunit for other coordinate types; cunit may be used to label coordinate values.
Fixes WCS keyvalues for malformed cylindrical projections.
Returns 0 for success; -1 if no change required.
string
Representative mid-point of the date of observation in ISO format, yyyy-mm-ddThh:mm:ss.
See also
dateobs
string
Start of the date of observation in ISO format, yyyy-mm-ddThh:mm:ss.
See also
dateavg
Translates the old DATE-OBS date format to year-2000 standard form (yyyy-mm-ddThh:mm:ss) and derives MJD-OBS from it if not already set. Alternatively, if mjdobs is set and dateobs isn’t, then datfix derives dateobs from it. If both are set but disagree by more than half a day then ValueError is raised.
Returns 0 for success; -1 if no change required.
double
The equinox associated with dynamical equatorial or ecliptic coordinate systems, EQUINOXa (or EPOCH in older headers). Not applicable to ICRS equatorial or ecliptic coordinates.
An undefined value is represented by NaN.
fix(translate_units=’‘, naxis=0)
Applies all of the corrections handled separately by datfix, unitfix, celfix, spcfix and cylfix.
translate_units: string. Do potentially unsafe translations of non-standard unit strings.
Although "S" is commonly used to represent seconds, its translation to "s" is potentially unsafe since the standard recognizes "S" formally as Siemens, however rarely that may be used. The same applies to "H" for hours (Henry), and "D" for days (Debye).
This string controls what to do in such cases, and is case-insensitive.
If the string contains "s", translate "S" to "s".
If the string contains "h", translate "H" to "h".
If the string contains "d", translate "D" to "d".
Thus '' doesn’t do any unsafe translations, whereas 'shd' does all of them.
naxis: int array[naxis]. Image axis lengths. If this array is set to zero or None, then cylfix will not be invoked.
Returns a dictionary containing the following keys, each referring to a status string for each of the sub-fix functions that were called:
get_ps() -> list of tuples
Returns PSi_ma keywords for each i and m. Returned as a list of tuples of the form (i, m, value):
- i: int. Axis number, as in PSi_ma, (i.e. 1-relative)
- m: int. Parameter number, as in PSi_ma, (i.e. 0-relative)
- value: string. Parameter value.
See also
set_ps
get_pv() -> list of tuples
Returns PVi_ma keywords for each i and m. Returned as a list of tuples of the form (i, m, value):
- i: int. Axis number, as in PVi_ma, (i.e. 1-relative)
- m: int. Parameter number, as in PVi_ma, (i.e. 0-relative)
- value: string. Parameter value.
See also
set_pv
has_cdi_ja() -> bool
Returns True if CDi_ja is present. CDi_ja is an alternate specification of the linear transformation matrix, maintained for historical compatibility.
Matrix elements in the IRAF convention are equivalent to the product CDi_ja = CDELTia * PCi_ja, but the defaults differ from that of the PCi_ja matrix. If one or more CDi_ja keywords are present then all unspecified CDi_ja default to zero. If no CDi_ja (or CROTAia) keywords are present, then the header is assumed to be in PCi_ja form whether or not any PCi_ja keywords are present since this results in an interpretation of CDELTia consistent with the original FITS specification.
While CDi_ja may not formally co-exist with PCi_ja, it may co-exist with CDELTia and CROTAia which are to be ignored.
See also
cd
has_crotaia() -> bool
Returns True if CROTAia is present. CROTAia is an alternate specification of the linear transformation matrix, maintained for historical compatibility.
In the AIPS convention, CROTAia may only be associated with the latitude axis of a celestial axis pair. It specifies a rotation in the image plane that is applied after the CDELTia; any other CROTAia keywords are ignored.
CROTAia may not formally co-exist with PCi_ja. CROTAia and CDELTia may formally co-exist with CDi_ja but if so are to be ignored.
See also
cd
has_pci_ja() -> bool
Returns True if PCi_ja is present. PCi_ja is the recommended way to specify the linear transformation matrix.
See also
cd
double array[2][2] (read-only)
Inverse of the matrix containing the product of the CDELTia diagonal matrix and the PCi_ja matrix.
Warning
This value may not be correct until after set is called.
is_unity() -> bool
Returns True if the linear transformation matrix (cd) is unity.
Warning
This value may not be correct until after set is called.
int (read-only)
The index into the sky coordinate array containing latitude values.
double
The native latitude of the celestial pole, LATPOLEa (deg).
string (read-only)
Celestial axis type for latitude, e.g. RA.
Warning
This value may not be correct until after set is called.
int (read-only)
The index into the sky coordinate array containing longitude values.
string (read-only)
Celestial axis type for longitude, e.g. DEC.
Warning
This value may not be correct until after set is called.
double
The native longitude of the celestial pole, LONPOLEa (deg).
mix(mixpix, mixcel, vspan, vstep, viter, world, pixcrd, origin) -> dict
Given either the celestial longitude or latitude plus an element of the pixel coordinate, solves for the remaining elements by iterating on the unknown celestial coordinate element using s2p.
Returns dictionary with the following keys:
Exceptions:
See also
lat, lng
double
Modified Julian Date (MJD = JD - 2400000.5), MJD-AVG, corresponding to DATE-AVG.
An undefined value is represented by NaN.
See also
mjdobs
double
Modified Julian Date (MJD = JD - 2400000.5), MJD-OBS, corresponding to DATE-OBS.
An undefined value is represented by NaN.
See also
mjdavg
string
The name given to the coordinate representation WCSNAMEa.
int (read-only)
The number of axes (pixel and coordinate), given by the NAXIS or WCSAXESa keyvalues.
The number of coordinate axes is determined at parsing time, and can not be subsequently changed.
It is determined from the highest of the following:
- NAXIS
- WCSAXESa
- The highest axis number in any parameterized WCS keyword. The keyvalue, as well as the keyword, must be syntactically valid otherwise it will not be considered.
If none of these keyword types is present, i.e. if the header only contains auxiliary WCS keywords for a particular coordinate representation, then no coordinate description is constructed for it.
This value may differ for different coordinate representations of the same image.
double array[3]
Location of the observer in a standard terrestrial reference frame, OBSGEO-X, OBSGEO-Y, OBSGEO-Z (in meters).
An undefined value is represented by NaN.
p2s(pixcrd, origin) -> dict
Converts pixel to sky coordinates.
Returns a dictionary with the following keys:
Exceptions:
See also
lat, lng
double array[2][2]
The PCi_ja (pixel coordinate) transformation matrix. The order is:
[[PC1_1, PC1_2],
[PC2_1, PC2_2]]
double
The native latitude of the fiducial point, i.e. the point whose celestial coordinates are given in ref[1:2]. If undefined (NaN) the initialization routine, set, will set this to a projection-specific default.
See also
theta0
double array[2][2] (read-only)
Matrix containing the product of the CDELTia diagonal matrix and the PCi_ja matrix.
Warning
This value may not be correct until after set is called.
print_contents()
Print the contents of the WCS object to stdout. Probably only useful for debugging purposes, and may be removed in the future.
string
The equatorial or ecliptic coordinate system type, RADESYSa.
double
Rest frequency (Hz) from RESTFRQa.
An undefined value is represented by NaN.
double
Rest wavelength (m) from RESTWAVa.
An undefined value is represented by NaN.
s2p(sky, origin) -> dict
Transforms sky coordinates to pixel coordinates.
Returns a dictionary with the following keys:
Exceptions:
See also
lat, lng
set()
Sets up a WCS object for use according to information supplied within it.
Note that this routine need not be called directly; it will be invoked by p2s and s2p if necessary.
Some attributes that are based on other attributes (such as lattyp on ctype) may not be correct until after set is called.
set strips off trailing blanks in all string members.
Among all the obvious details, set recognizes the NCP projection and converts it to the equivalent SIN projection and it also recognizes GLS as a synonym for SFL. It does alias translation for the AIPS spectral types (FREQ-LSR, FELO-HEL, etc.) but without changing the input header keywords.
Exceptions:
set_ps(list)
Sets PSi_ma keywords for each i and m. The input must be a sequence of tuples of the form (i, m, value):
- i: int. Axis number, as in PSi_ma, (i.e. 1-relative)
- m: int. Parameter number, as in PSi_ma, (i.e. 0-relative)
- value: string. Parameter value.
See also
get_ps
set_pv(list)
Sets PVi_ma keywords for each i and m. The input must be a sequence of tuples of the form (i, m, value):
- i: int. Axis number, as in PVi_ma, (i.e. 1-relative)
- m: int. Parameter number, as in PVi_ma, (i.e. 0-relative)
- value: string. Parameter value.
spcfix() -> int
Translates AIPS-convention spectral coordinate types. {FREQ, VELO, FELO}-{OBS, HEL, LSR} (e.g. FREQ-LSR, VELO-OBS, FELO-HEL)
Returns 0 for success; -1 if no change required.
int (read-only)
The index containing the spectral axis values.
string
Spectral reference frame (standard of rest), SPECSYSa.
See also
ssysobs, velosys.
sptr(ctype, i=-1)
Translates the spectral axis in a WCS object. For example, a FREQ axis may be translated into ZOPT-F2W and vice versa.
Exceptions:
string
The actual spectral reference frame in which there is no differential variation in the spectral coordinate across the field-of-view, SSYSOBSa.
See also
specsys, velosys
string
The spectral reference frame (standard of rest) in which the redshift was measured, SSYSSRCa.
sub(axes) -> WCS object
Extracts the coordinate description for a subimage from a WCS object.
The world coordinate system of the subimage must be separable in the sense that the world coordinates at any point in the subimage must depend only on the pixel coordinates of the axes extracted. In practice, this means that the PCi_ja matrix of the original image must not contain non-zero off-diagonal terms that associate any of the subimage axes with any of the non-subimage axes.
Coordinate axes types may be specified using either strings or special integer constants. The available types are:
- 'longitude' / WCSSUB_LONGITUDE: Celestial longitude
- 'latitude' / WCSSUB_LATITUDE: Celestial latitude
- 'cubeface' / WCSSUB_CUBEFACE: Quadcube CUBEFACE axis
- 'spectral' / WCSSUB_SPECTRAL: Spectral axis
- 'stokes' / WCSSUB_STOKES: Stokes axis
Returns a WCS object.
Exceptions:
double
The native longitude of the fiducial point, i.e. the point whose celestial coordinates are given in ref[1:2]. If undefined (NaN) the initialization routine, set, will set this to a projection-specific default.
See also
phi0
to_header(relax=False) -> string
to_header translates a WCS object into a FITS header.
- If the colnum member is non-zero then a binary table image array header will be produced.
- Otherwise, if the colax member is set non-zero then a pixel list header will be produced.
- Otherwise, a primary image or image extension header will be produced.
The output header will almost certainly differ from the input in a number of respects:
- The output header only contains WCS-related keywords. In particular, it does not contain syntactically-required keywords such as SIMPLE, NAXIS, BITPIX, or END.
- Deprecated (e.g. CROTAn) or non-standard usage will be translated to standard (this is partially dependent on whether fix was applied).
- Quantities will be converted to the units used internally, basically SI with the addition of degrees.
- Floating-point quantities may be given to a different decimal precision.
- Elements of the PCi_j matrix will be written if and only if they differ from the unit matrix. Thus, if the matrix is unity then no elements will be written.
- Additional keywords such as WCSAXES, CUNITia, LONPOLEa and LATPOLEa may appear.
- The original keycomments will be lost, although to_header tries hard to write meaningful comments.
- Keyword order may be changed.
Keywords can be translated between the image array, binary table, and pixel lists forms by manipulating the colnum or colax members of the WCS object.
relax: Degree of permissiveness:
- False: Recognize only FITS keywords defined by the published WCS standard.
- True: Admit all recognized informal extensions of the WCS standard.
Returns a raw FITS header as a string.
unitfix(translate_units=’‘) -> int
Translates non-standard CUNITia keyvalues. For example, DEG -> deg, also stripping off unnecessary whitespace.
translate_units: string. Do potentially unsafe translations of non-standard unit strings.
Although "S" is commonly used to represent seconds, its recognizes "S" formally as Siemens, however rarely that may be translation to "s" is potentially unsafe since the standard used. The same applies to "H" for hours (Henry), and "D" for days (Debye).
This string controls what to do in such cases, and is case-insensitive.
Thus '' doesn’t do any unsafe translations, whereas 'shd' does all of them.
Returns 0 for success; -1 if no change required.
double
The angle in degrees that should be used to decompose an observed velocity into radial and transverse components.
An undefined value is represented by NaN.
double
The relative radial velocity (m/s) between the observer and the selected standard of rest in the direction of the celestial reference coordinate, VELOSYSa.
An undefined value is represented by NaN.
See also
specsys, ssysobs
double
The redshift, ZSOURCEa, of the source.
An undefined value is represented by NaN.
DistortionLookupTable(table, crpix, crval, cdelt)
Represents a single lookup table for a Paper IV distortion transformation.
double array[naxis]
Coordinate increments (CDELTia) for each coord axis.
An undefined value is represented by NaN.
double array[naxis]
Coordinate reference pixels (CRPIXja) for each pixel axis.
double array[naxis]
Coordinate reference values (CRVALia) for each coordinate axis.
float array
The array data for the DistortionLookupTable.
get_offset(x, y) -> (x, y)
Returns the offset from the distortion table for pixel point (x, y).
Sip(a, b, ap, bp, crpix)
The Sip class performs polynomial distortion correction using the SIP convention in both directions.
Shupe, D. L., M. Moshir, J. Li, D. Makovoz and R. Narron. 2005. “The SIP Convention for Representing Distortion in FITS Image Headers.” ADASS XIV.
double array[a_order+1][a_order+1]
The SIP A_i_j matrix used for pixel to focal plane transformation.
Its values may be changed in place, but it may not be resized, without creating a new Sip object.
double array[ap_order+1][ap_order+1]
The SIP AP_i_j matrix used for focal plane to pixel transformation. Its values may be changed in place, but it may not be resized, without creating a new Sip object.
double array[b_order+1][b_order+1]
The SIP B_i_j matrix used for pixel to focal plane transformation. Its values may be changed in place, but it may not be resized, without creating a new Sip object.
double array[bp_order+1][bp_order+1]
The SIP BP_i_j matrix used for focal plane to pixel transformation. Its values may be changed in place, but it may not be resized, without creating a new Sip object.
double array[naxis]
Coordinate reference pixels (CRPIXja) for each pixel axis.
sip_foc2pix(foccrd, origin) -> double array[ncoord][nelem]
Convert focal plane coordinates to pixel coordinates using the SIP polynomial distortion convention.
Returns an array of pixel coordinates.
Exceptions:
sip_pix2foc(pixcrd, origin) -> double array[ncoord][nelem]
Convert pixel coordinates to focal plane coordinates using the SIP polynomial distortion convention.
Returns an array of focal plane coordinates.
Exceptions: