This example loads a FITS file (supplied on the command line) and uses the FITS keywords in its primary header to create a WCS and 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

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

# Parse the WCS keywords in the primary HDU

# 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__':