World Coordinate System (astropy.wcs)

Introduction

astropy.wcs contains utilities for managing World Coordinate System (WCS) transformations in FITS files. These transformations map the pixel locations in an image to their real-world units, such as their position on the sky sphere. These transformations can work both forward (from pixel to sky) and backward (from sky to pixel).

It performs three separate classes of WCS transformations:

Each of these transformations can be used independently or together in a standard pipeline.

Getting Started

The basic workflow is as follows:

  1. from astropy import wcs

  2. Call the WCS constructor with an astropy.io.fits Header and/or HDUList object.

  3. Optionally, if the FITS file uses any deprecated or non-standard features, you may need to call one of the fix methods on the object.

  4. Use one of the following transformation methods:

    • From pixels to world coordinates:

      • all_pix2world: Perform all three transformations in series (core WCS, SIP and table lookup distortions) from pixel to world coordinates. Use this one if you’re not sure which to use.

      • wcs_pix2world: Perform just the core WCS transformation from pixel to world coordinates.

    • From world to pixel coordinates:

      • all_world2pix: Perform all three transformations (core WCS, SIP and table lookup distortions) from world to pixel coordinates, using an iterative method if necessary.

      • wcs_world2pix: Perform just the core WCS transformation from world to pixel coordinates.

    • Performing SIP transformations only:

      • sip_pix2foc: Convert from pixel to focal plane coordinates using the SIP polynomial coefficients.

      • sip_foc2pix: Convert from focal plane to pixel coordinates using the SIP polynomial coefficients.

    • Performing distortion paper transformations only:

      • p4_pix2foc: Convert from pixel to focal plane coordinates using the table lookup distortion method described in the FITS WCS distortion paper.

      • det2im: Convert from detector coordinates to image coordinates. Commonly used for narrow column correction.

For example, to convert pixel coordinates from a two dimensional image to world coordinates:

>>> from astropy.wcs import WCS
>>> w = WCS('image.fits')  
>>> lon, lat = w.all_pix2world(30, 40, 0)
>>> print(lon, lat)
31.0 41.0

Using astropy.wcs

Loading WCS information from a FITS file

This example loads a FITS file (supplied on the commandline) and uses the WCS cards in its primary header to transform.

# Load the WCS information from a fits header, and use it
# to convert pixel coordinates to world coordinates.

import numpy as np
from astropy import wcs
from astropy.io import fits
import sys


def load_wcs_from_file(filename):
    # Load the FITS hdulist using astropy.io.fits
    hdulist = fits.open(filename)

    # Parse the WCS keywords in the primary HDU
    w = wcs.WCS(hdulist[0].header)

    # Print out the "name" of the WCS, as defined in the FITS header
    print(w.wcs.name)

    # Print out all of the settings that were parsed from the header
    w.wcs.print_contents()

    # Three pixel coordinates of interest.
    # Note we've silently assumed an NAXIS=2 image here.
    # The pixel coordinates are pairs of [X, Y].
    # The "origin" argument indicates whether the input coordinates
    # are 0-based (as in Numpy arrays) or
    # 1-based (as in the FITS convention, for example coordinates
    # coming from DS9).
    pixcrd = np.array([[0, 0], [24, 38], [45, 98]], dtype=np.float64)

    # Convert pixel coordinates to world coordinates
    # The second argument is "origin" -- in this case we're declaring we
    # have 0-based (Numpy-like) coordinates.
    world = w.wcs_pix2world(pixcrd, 0)
    print(world)

    # Convert the same coordinates back to pixel coordinates.
    pixcrd2 = w.wcs_world2pix(world, 0)
    print(pixcrd2)

    # These should be the same as the original pixel coordinates, modulo
    # some floating-point error.
    assert np.max(np.abs(pixcrd - pixcrd2)) < 1e-6

    # The example below illustrates the use of "origin" to convert between
    # 0- and 1- based coordinates when executing the forward and backward
    # WCS transform.
    x = 0
    y = 0
    origin = 0
    assert (w.wcs_pix2world(x, y, origin) ==
            w.wcs_pix2world(x + 1, y + 1, origin + 1))


if __name__ == '__main__':
    load_wcs_from_file(sys.argv[-1])

Building a WCS structure programmatically

This example, rather than starting from a FITS header, sets WCS values programmatically, uses those settings to transform some points, and then saves those settings to a new FITS header.

# Set the WCS information manually by setting properties of the WCS
# object.

import numpy as np
from astropy import wcs
from astropy.io import fits

# Create a new WCS object.  The number of axes must be set
# from the start
w = wcs.WCS(naxis=2)

# Set up an "Airy's zenithal" projection
# Vector properties may be set with Python lists, or Numpy arrays
w.wcs.crpix = [-234.75, 8.3393]
w.wcs.cdelt = np.array([-0.066667, 0.066667])
w.wcs.crval = [0, -90]
w.wcs.ctype = ["RA---AIR", "DEC--AIR"]
w.wcs.set_pv([(2, 1, 45.0)])

# Three pixel coordinates of interest.
# The pixel coordinates are pairs of [X, Y].
# The "origin" argument indicates whether the input coordinates
# are 0-based (as in Numpy arrays) or
# 1-based (as in the FITS convention, for example coordinates
# coming from DS9).
pixcrd = np.array([[0, 0], [24, 38], [45, 98]], dtype=np.float64)

# Convert pixel coordinates to world coordinates.
# The second argument is "origin" -- in this case we're declaring we
# have 0-based (Numpy-like) coordinates.
world = w.wcs_pix2world(pixcrd, 0)
print(world)

# Convert the same coordinates back to pixel coordinates.
pixcrd2 = w.wcs_world2pix(world, 0)
print(pixcrd2)

# These should be the same as the original pixel coordinates, modulo
# some floating-point error.
assert np.max(np.abs(pixcrd - pixcrd2)) < 1e-6

# The example below illustrates the use of "origin" to convert between
# 0- and 1- based coordinates when executing the forward and backward
# WCS transform.
x = 0
y = 0
origin = 0
assert (w.wcs_pix2world(x, y, origin) ==
        w.wcs_pix2world(x + 1, y + 1, origin + 1))

# Now, write out the WCS object as a FITS header
header = w.to_header()

# header is an astropy.io.fits.Header object.  We can use it to create a new
# PrimaryHDU and write it to a file.
hdu = fits.PrimaryHDU(header=header)
# Save to FITS file
# hdu.writeto('test.fits')

Note

The members of the WCS object correspond roughly to the key/value pairs in the FITS header. However, they are adjusted and normalized in a number of ways that make performing the WCS transformation easier. Therefore, they can not be relied upon to get the original values in the header. To build up a FITS header directly and specifically, use astropy.io.fits.Header directly.

Validating the WCS keywords in a FITS file

Astropy includes a commandline tool, wcslint to check the WCS keywords in a FITS file:

> wcslint invalid.fits
HDU 1:
  WCS key ' ':
    - RADECSYS= 'ICRS ' / Astrometric system
      RADECSYS is non-standard, use RADESYSa.
    - The WCS transformation has more axes (2) than the image it is
      associated with (0)
    - 'celfix' made the change 'PV1_5 : Unrecognized coordinate
      transformation parameter'.

HDU 2:
  WCS key ' ':
    - The WCS transformation has more axes (3) than the image it is
      associated with (0)
    - 'celfix' made the change 'In CUNIT2 : Mismatched units type
      'length': have 'Hz', want 'm''.
    - 'unitfix' made the change 'Changed units: 'HZ      ' -> 'Hz''.

Bounds checking

Bounds checking is enabled by default, and any computed world coordinates outside of [-180°, 180°] for longitude and [-90°, 90°] in latitude are marked as invalid. To disable this behavior, use astropy.wcs.Wcsprm.bounds_check.

Supported projections

As astropy.wcs is based on wcslib, it supports the standard projections defined in the FITS WCS standard. These projection codes are specified in the second part of the CTYPEn keywords (accessible through Wcsprm.ctype), for example, RA---TAN-SIP. The supported projection codes are:

  • AZP: zenithal/azimuthal perspective

  • SZP: slant zenithal perspective

  • TAN: gnomonic

  • STG: stereographic

  • SIN: orthographic/synthesis

  • ARC: zenithal/azimuthal equidistant

  • ZPN: zenithal/azimuthal polynomial

  • ZEA: zenithal/azimuthal equal area

  • AIR: Airy’s projection

  • CYP: cylindrical perspective

  • CEA: cylindrical equal area

  • CAR: plate carrée

  • MER: Mercator’s projection

  • COP: conic perspective

  • COE: conic equal area

  • COD: conic equidistant

  • COO: conic orthomorphic

  • SFL: Sanson-Flamsteed (“global sinusoid”)

  • PAR: parabolic

  • MOL: Mollweide’s projection

  • AIT: Hammer-Aitoff

  • BON: Bonne’s projection

  • PCO: polyconic

  • TSC: tangential spherical cube

  • CSC: COBE quadrilateralized spherical cube

  • QSC: quadrilateralized spherical cube

  • HPX: HEALPix

  • XPH: HEALPix polar, aka “butterfly”

And, if built with wcslib 5.0 or later, the following polynomial distortions are supported:

  • TPV: Polynomial distortion

  • TUV: Polynomial distortion

Note

Though wcslib 5.4 and later handles SIP polynomial distortion, for backward compatibility, SIP is handled by astropy itself and methods exist to handle it specially.

Subsetting and Pixel Scales

WCS objects can be broken apart into their constituent axes using the sub function. There is also a celestial convenience function that will return a WCS object with only the celestial axes included.

The pixel scales of a celestial image or the pixel dimensions of a non-celestial image can be extracted with the utility functions proj_plane_pixel_scales and non_celestial_pixel_scales. Likewise, celestial pixel area can be extracted with the utility function proj_plane_pixel_area.

Matplotlib plots with correct WCS projection

The WCSAxes framework, previously a standalone package, allows the WCS to be used to define projections in Matplotlib. More information on using WCSAxes can be found here.

from matplotlib import pyplot as plt
from astropy.io import fits
from astropy.wcs import WCS
from astropy.utils.data import get_pkg_data_filename

filename = get_pkg_data_filename('tutorials/FITS-images/HorseHead.fits')

hdu = fits.open(filename)[0]
wcs = WCS(hdu.header)

fig = plt.figure()
fig.add_subplot(111, projection=wcs)
plt.imshow(hdu.data, origin='lower', cmap=plt.cm.viridis)
plt.xlabel('RA')
plt.ylabel('Dec')

(png, svg, pdf)

../_images/index-18.png

See Also

Reference/API

astropy.wcs Package

astropy.wcs contains utilities for managing World Coordinate System (WCS) transformations in FITS files. These transformations map the pixel locations in an image to their real-world units, such as their position on the sky sphere.

It performs three separate classes of WCS transformations:

Each of these transformations can be used independently or together in a standard pipeline.

Functions

find_all_wcs(header[, relax, keysel, fix, …])

Find all the WCS transformations in the given header.

get_include()

Get the path to astropy.wcs’s C header files.

validate(source)

Prints a WCS validation report for the given FITS file.

Classes

Auxprm

Class that contains auxiliary coordinate system information of a specialist nature.

DistortionLookupTable(*table*, *crpix*, …)

Represents a single lookup table for a distortion paper transformation.

FITSFixedWarning

The warning raised when the contents of the FITS header have been modified to be standards compliant.

InconsistentAxisTypesError()

The WCS header inconsistent or unrecognized coordinate axis type(s).

InvalidCoordinateError()

One or more of the world coordinates is invalid.

InvalidSubimageSpecificationError()

The subimage specification is invalid.

InvalidTabularParametersError()

The given tabular parameters are invalid.

InvalidTransformError()

The WCS transformation is invalid, or the transformation parameters are invalid.

NoConvergence(*args[, best_solution, …])

An error class used to report non-convergence and/or divergence of numerical methods.

NoSolutionError()

No solution can be found in the given interval.

NoWcsKeywordsFoundError()

No WCS keywords were found in the given header.

NonseparableSubimageCoordinateSystemError()

Non-separable subimage coordinate system.

SingularMatrixError()

The linear transformation matrix is singular.

Sip(*a, b, ap, bp, crpix*)

The Sip class performs polynomial distortion correction using the SIP convention in both directions.

Tabprm

A class to store the information related to tabular coordinates, i.e., coordinates that are defined via a lookup table.

WCS([header, fobj, key, minerr, relax, …])

WCS objects perform standard WCS transformations, and correct for SIP and distortion paper table-lookup transformations, based on the WCS keywords and supplementary data read from a FITS file.

WCSBase(*sip, cpdis, wcsprm, det2im*)

Wcs objects amalgamate basic WCS (as provided by wcslib), with SIP and distortion paper operations.

WcsError

Base class of all invalid WCS errors.

Wcsprm([header, key, relax, naxis, keysel, …])

Wcsprm performs the core WCS transformations.

Wtbarr

Classes to construct coordinate lookup tables from a binary table extension (BINTABLE).

Class Inheritance Diagram

Inheritance diagram of astropy.wcs.Auxprm, astropy.wcs.DistortionLookupTable, astropy.wcs.wcs.FITSFixedWarning, astropy.wcs._wcs.InconsistentAxisTypesError, astropy.wcs._wcs.InvalidCoordinateError, astropy.wcs._wcs.InvalidSubimageSpecificationError, astropy.wcs._wcs.InvalidTabularParametersError, astropy.wcs._wcs.InvalidTransformError, astropy.wcs.wcs.NoConvergence, astropy.wcs._wcs.NoSolutionError, astropy.wcs._wcs.NoWcsKeywordsFoundError, astropy.wcs._wcs.NonseparableSubimageCoordinateSystemError, astropy.wcs._wcs.SingularMatrixError, astropy.wcs.Sip, astropy.wcs.Tabprm, astropy.wcs.wcs.WCS, astropy.wcs.WCSBase, astropy.wcs._wcs.WcsError, astropy.wcs.Wcsprm, astropy.wcs.Wtbarr

astropy.wcs.utils Module

Functions

add_stokes_axis_to_wcs(wcs, add_before_ind)

Add a new Stokes axis that is uncorrelated with any other axes.

celestial_frame_to_wcs(frame[, projection])

For a given coordinate frame, return the corresponding WCS object.

wcs_to_celestial_frame(wcs)

For a given WCS, return the coordinate frame that matches the celestial component of the WCS.

proj_plane_pixel_scales(wcs)

For a WCS returns pixel scales along each axis of the image pixel at the CRPIX location once it is projected onto the “plane of intermediate world coordinates” as defined in Greisen & Calabretta 2002, A&A, 395, 1061.

proj_plane_pixel_area(wcs)

For a celestial WCS (see astropy.wcs.WCS.celestial) returns pixel area of the image pixel at the CRPIX location once it is projected onto the “plane of intermediate world coordinates” as defined in Greisen & Calabretta 2002, A&A, 395, 1061.

is_proj_plane_distorted(wcs[, maxerr])

For a WCS returns False if square image (detector) pixels stay square when projected onto the “plane of intermediate world coordinates” as defined in Greisen & Calabretta 2002, A&A, 395, 1061.It will return True if transformation from image (detector) coordinates to the focal plane coordinates is non-orthogonal or if WCS contains non-linear (e.g., SIP) distortions..

non_celestial_pixel_scales(inwcs)

Calculate the pixel scale along each axis of a non-celestial WCS, for example one with mixed spectral and spatial axes.

skycoord_to_pixel(coords, wcs[, origin, mode])

Convert a set of SkyCoord coordinates into pixels.

pixel_to_skycoord(xp, yp, wcs[, origin, …])

Convert a set of pixel coordinates into a SkyCoord coordinate.

pixel_to_pixel(wcs_in, wcs_out, *inputs)

Transform pixel coordinates in a dataset with a WCS to pixel coordinates in another dataset with a different WCS.

local_partial_pixel_derivatives(wcs, *pixel)

Return a matrix of shape (world_n_dim, pixel_n_dim) where each entry [i, j] is the partial derivative d(world_i)/d(pixel_j) at the requested pixel position.

fit_wcs_from_points(xy, world_coords[, …])

Given two matching sets of coordinates on detector and sky, compute the WCS.

Class Inheritance Diagram

Inheritance diagram of astropy.wcs.utils.custom_wcs_to_frame_mappings, astropy.wcs.utils.custom_frame_to_wcs_mappings

Acknowledgments and Licenses

wcslib is licenced under the GNU Lesser General Public License.