GWCS Documentation

GWCS is a package for managing the World Coordinate System (WCS) of astronomical data.

Introduction & Motivation for GWCS

The mapping from ‘pixel’ coordinates to corresponding ‘real-world’ coordinates (e.g. celestial coordinates, spectral wavelength) is crucial to relating astronomical data to the phenomena they describe. Images and other types of data often come encoded with information that describes this mapping – this is referred to as the ‘World Coordinate System’ or WCS. The term WCS is often used to refer specifically to the most widely used ‘FITS implementation of WCS’, but here unless specified WCS refers to the broader concept of relating pixel ⟷ world. (See the discussion in APE14 for more on this topic).

The FITS WCS standard, currently the most widely used method of encoding WCS in data, describes a set of required FITS header keywords and allowed values that describe how pixel ⟷ world transformations should be done. This current paradigm of encoding data with only instructions on how to relate pixel to world, separate from the transformation machinery itself, has several limitations:

  • Limited flexibility. WCS keywords and their values are rigidly defined so that the instructions are unambiguous. This places limitations on, for example, describing geometric distortion in images since only a handful of distortion models are defined in the FITS standard (and therefore can be encoded in FITS headers as WCS information).

  • Separation of data from transformation pipelines. The machinery that transforms pixel ⟷ world does not exist along side the data – there is merely a roadmap for how one would do the transformation. External packages and libraries (e.g wcslib, or its Python interface astropy.wcs) must be written to interpret the instructions and execute the transformation. These libraries don’t allow easy access to coordinate frames along the course of the full pixel to world transformation pipeline. Additionally, since these libraries can only interpret FITS WCS information, any custom ‘WCS’ definitions outside of FITS require the user to write their own transformation pipelines.

  • Incompatibility with varying file formats. New file formats that are becoming more widely used in place of FITS to store astronomical data, like the ASDF format, also require a method of encoding WCS information. FITS WCS and the accompanying libraries are adapted for FITS only. A more flexible interface would be agnostic to file type, as long as the necessary information is present.

The GWCS package and GWCS object is a generalized WCS implementation that mitigates these limitations. The goal of the GWCS package is to provide a flexible toolkit for expressing and evaluating transformations between pixel and world coordinates, as well as intermediate frames along the course of this transformation.The GWCS object supports a data model which includes the entire transformation pipeline from input pixel coordinates to world coordinates (and vice versa). The basis of the GWCS object is astropy modeling. Models that describe the pixel ⟷ world transformations can be chained, joined or combined with arithmetic operators using the flexible framework of compound models in modeling. This approach allows for easy access to intermediate frames. In the case of a celestial output frame coordinates provides further transformations between standard celestial coordinate frames. Spectral output coordinates are instances of Quantity and can be transformed to other units with the tools in that package. Time coordinates are instances of Time. GWCS supports transforms initialized with Quantity objects ensuring automatic unit conversion.

Pixel Conventions and Definitions

This API assumes that integer pixel values fall at the center of pixels (as assumed in the FITS-WCS standard, see Section 2.1.4 of Greisen et al., 2002, A&A 446, 747), while at the same time matching the Python 0-index philosophy. That is, the first pixel is considered pixel 0, but pixel coordinates (0, 0) are the center of that pixel. Hence the first pixel spans pixel values -0.5 to 0.5.

There are two main conventions for ordering pixel coordinates. In the context of 2-dimensional imaging data/arrays, one can either think of the pixel coordinates as traditional Cartesian coordinates (which we call x and y here), which are usually given with the horizontal coordinate (x) first, and the vertical coordinate (y) second, meaning that pixel coordinates would be given as (x, y). Alternatively, one can give the coordinates by first giving the row in the data, then the column, i.e. (row, column). While the former is a more common convention when e.g. plotting (think for example of the Matplotlib scatter(x, y) method), the latter is the convention used when accessing values from e.g. Numpy arrays that represent images (image[row, column]).

The GWCS object assumes Cartesian order (x, y), however the Common Interface for World Coordinate System - APE 14 accepts both conventions. The order of the pixel coordinates ((x, y) vs (row, column)) in the Common API depends on the method or property used, and this can normally be determined from the property or method name. Properties and methods containing pixel assume (x, y) ordering, while properties and methods containing array assume (row, column) ordering.

Installation

gwcs requires:

To install from source:

git clone https://github.com/spacetelescope/gwcs.git
cd gwcs
python setup.py install

To install the latest release:

pip install gwcs

The latest release of GWCS is also available as part of astroconda.

Getting Started

The WCS data model represents a pipeline of transformations between two coordinate frames, the final one usually a physical coordinate system. It is represented as a list of steps executed in order. Each step defines a starting coordinate frame and the transform to the next frame in the pipeline. The last step has no transform, only a frame which is the output frame of the total transform. As a minimum a WCS object has an input_frame (defaults to “detector”), an output_frame and the transform between them.

The WCS is validated using the ASDF Standard and serialized to file using the asdf package. There are two ways to save the WCS to a file:

A step by step example of constructing an imaging GWCS object.

The following example shows how to construct a GWCS object equivalent to a FITS imaging WCS without distortion, defined in this FITS imaging header:

WCSAXES =                    2 / Number of coordinate axes
WCSNAME = '47 Tuc     '        / Coordinate system title
CRPIX1  =               2048.0 / Pixel coordinate of reference point
CRPIX2  =               1024.0 / Pixel coordinate of reference point
PC1_1   =   1.290551569736E-05 / Coordinate transformation matrix element
PC1_2   =  5.9525007864732E-06 / Coordinate transformation matrix element
PC2_1   =  5.0226382102765E-06 / Coordinate transformation matrix element
PC2_2   = -1.2644844123757E-05 / Coordinate transformation matrix element
CDELT1  =                  1.0 / [deg] Coordinate increment at reference point
CDELT2  =                  1.0 / [deg] Coordinate increment at reference point
CUNIT1  = 'deg'                / Units of coordinate increment and value
CUNIT2  = 'deg'                / Units of coordinate increment and value
CTYPE1  = 'RA---TAN'           / TAN (gnomonic) projection + SIP distortions
CTYPE2  = 'DEC--TAN'           / TAN (gnomonic) projection + SIP distortions
CRVAL1  =        5.63056810618 / [deg] Coordinate value at reference point
CRVAL2  =      -72.05457184279 / [deg] Coordinate value at reference point
LONPOLE =                180.0 / [deg] Native longitude of celestial pole
LATPOLE =      -72.05457184279 / [deg] Native latitude of celestial pole
RADESYS = 'ICRS'                / Equatorial coordinate system

The following imports are generally useful:

>>> import numpy as np
>>> from astropy.modeling import models
>>> from astropy import coordinates as coord
>>> from astropy import units as u
>>> from gwcs import wcs
>>> from gwcs import coordinate_frames as cf

The forward_transform is constructed as a combined model using astropy.modeling. The frames are subclasses of CoordinateFrame. Although strings are acceptable as coordinate_frames it is recommended this is used only in testing/debugging.

Using the modeling package create a combined model to transform detector coordinates to ICRS following the FITS WCS standard convention.

First, create a transform which shifts the input x and y coordinates by CRPIX. We subtract 1 from the CRPIX values because the first pixel is considered pixel 1 in FITS WCS:

>>> shift_by_crpix = models.Shift(-(2048 - 1)*u.pix) & models.Shift(-(1024 - 1)*u.pix)

Create a transform which rotates the inputs using the PC matrix.

>>> matrix = np.array([[1.290551569736E-05, 5.9525007864732E-06],
...                    [5.0226382102765E-06 , -1.2644844123757E-05]])
>>> rotation = models.AffineTransformation2D(matrix * u.deg,
...                                          translation=[0, 0] * u.deg)
>>> rotation.input_units_equivalencies = {"x": u.pixel_scale(1*u.deg/u.pix),
...                                       "y": u.pixel_scale(1*u.deg/u.pix)}
>>> rotation.inverse = models.AffineTransformation2D(np.linalg.inv(matrix) * u.pix,
...                                                  translation=[0, 0] * u.pix)
>>> rotation.inverse.input_units_equivalencies = {"x": u.pixel_scale(1*u.pix/u.deg),
...                                               "y": u.pixel_scale(1*u.pix/u.deg)}

Create a tangent projection and a rotation on the sky using CRVAL.

>>> tan = models.Pix2Sky_TAN()
>>> celestial_rotation =  models.RotateNative2Celestial(5.63056810618*u.deg, -72.05457184279*u.deg, 180*u.deg)

>>> det2sky = shift_by_crpix | rotation | tan | celestial_rotation
>>> det2sky.name = "linear_transform"

Create a detector coordinate frame and a celestial ICRS frame.

>>> detector_frame = cf.Frame2D(name="detector", axes_names=("x", "y"),
...                             unit=(u.pix, u.pix))
>>> sky_frame = cf.CelestialFrame(reference_frame=coord.ICRS(), name='icrs',
...                               unit=(u.deg, u.deg))

This WCS pipeline has only one step - from detector to sky:

>>> pipeline = [(detector_frame, det2sky),
...             (sky_frame, None)
...            ]
>>> wcsobj = wcs.WCS(pipeline)
>>> print(wcsobj)
  From      Transform
-------- ----------------
detector linear_transform
    icrs             None

To convert a pixel (x, y) = (1, 2) to sky coordinates, call the WCS object as a function:

>>> sky = wcsobj(1*u.pix, 2*u.pix, with_units=True)
>>> print(sky)
<SkyCoord (ICRS): (ra, dec) in deg
  (5.52515954, -72.05190935)>

The invert() method evaluates the backward_transform() if available, otherwise applies an iterative method to calculate the reverse coordinates.

>>> wcsobj.invert(sky)
(<Quantity 1. pix>, <Quantity 2. pix>)

Save a WCS object as a pure ASDF file

>>> from asdf import AsdfFile
>>> tree = {"wcs": wcsobj}
>>> wcs_file = AsdfFile(tree)
>>> wcs_file.write_to("imaging_wcs.asdf")

Listing of imaging_wcs.asdf

Save a WCS object as an ASDF extension in a FITS file

>>> from astropy.io import fits
>>> from asdf import fits_embed
>>> hdul = fits.open("example_imaging.fits")
>>> hdul.info()
Filename: example_imaging.fits
No.    Name      Ver    Type      Cards   Dimensions   Format
0  PRIMARY       1 PrimaryHDU     775   ()
1  SCI           1 ImageHDU        71   (600, 550)   float32
>>> tree = {"sci": hdul[1].data,
...         "wcs": wcsobj}
>>> fa = fits.embed.AsdfInFits(hdul, tree)
>>> fa.write_to("imaging_with_wcs_in_asdf.fits")
>>> fits.info("imaging_with_wcs_in_asdf.fits")
Filename: example_with_wcs.asdf
No.    Name      Ver    Type      Cards   Dimensions   Format
0  PRIMARY       1 PrimaryHDU     775   ()
1  SCI           1 ImageHDU        71   (600, 550)   float32
2  ASDF          1 BinTableHDU     11   1R x 1C   [5086B]

Reading a WCS object from a file

ASDF is used to read a WCS object from a pure ASDF file or from an ASDF extension in a FITS file.

>>> import asdf
>>> asdf_file = asdf.open("imaging_wcs.asdf")
>>> wcsobj = asdf_file.tree['wcs']

>>> import asdf
>>> fa = asdf.open("imaging_with_wcs_in_asdf.fits")
>>> wcsobj = fa.tree["wcs"]

Other Examples

Adding distortion to the imaging example

Let’s expand the WCS created in Getting Started by adding a polynomial distortion correction.

Because the polynomial models in modeling do not support units yet, this example will use transforms without units. At the end the units associated with the output frame are used to create a SkyCoord object.

The imaging example without units:

>>> import numpy as np
>>> from astropy.modeling import models
>>> from astropy import coordinates as coord
>>> from astropy import units as u
>>> from gwcs import wcs
>>> from gwcs import coordinate_frames as cf

>>> crpix = (2048, 1024)
>>> shift_by_crpix = models.Shift(-crpix[0]) & models.Shift(-crpix[1])
>>> matrix = np.array([[1.290551569736E-05, 5.9525007864732E-06],
...                    [5.0226382102765E-06 , -1.2644844123757E-05]])
>>> rotation = models.AffineTransformation2D(matrix)
>>> rotation.inverse = models.AffineTransformation2D(np.linalg.inv(matrix))
>>> tan = models.Pix2Sky_TAN()
>>> celestial_rotation =  models.RotateNative2Celestial(5.63056810618, -72.05457184279, 180)
>>> det2sky = shift_by_crpix | rotation | tan | celestial_rotation
>>> det2sky.name = "linear_transform"
>>> detector_frame = cf.Frame2D(name="detector", axes_names=("x", "y"),
...                             unit=(u.pix, u.pix))
>>> sky_frame = cf.CelestialFrame(reference_frame=coord.ICRS(), name='icrs',
...                               unit=(u.deg, u.deg))
>>> pipeline = [(detector_frame, det2sky),
...             (sky_frame, None)
...            ]
>>> wcsobj = wcs.WCS(pipeline)
>>> print(wcsobj)
  From      Transform
-------- ----------------
detector linear_transform
    icrs             None

First create distortion corrections represented by a polynomial model of fourth degree. The example uses the astropy Polynomial2D and Mapping models.

>>> poly_x = models.Polynomial2D(4)
>>> poly_x.parameters = [0, 1, 8.55e-06, -4.73e-10, 2.37e-14, 0, -5.20e-06,
...                      -3.98e-11, 1.97e-15, 2.17e-06, -5.23e-10, 3.47e-14,
...                      1.08e-11, -2.46e-14, 1.49e-14]
>>> poly_y = models.Polynomial2D(4)
>>> poly_y.parameters = [0, 0, -1.75e-06, 8.57e-11, -1.77e-14, 1, 6.18e-06,
...                      -5.09e-10, -3.78e-15, -7.22e-06, -6.17e-11,
...                      -3.66e-14, -4.18e-10, 1.22e-14, -9.96e-15]
>>> distortion = ((models.Shift(-crpix[0]) & models.Shift(-crpix[1])) |
...               models.Mapping((0, 1, 0, 1)) | (poly_x & poly_y) |
...               (models.Shift(crpix[0]) & models.Shift(crpix[1])))
>>> distortion.name = "distortion"

Create an intermediate frame for distortion free coordinates.

>>> undistorted_frame = cf.Frame2D(name="undistorted_frame", unit=(u.pix, u.pix),
...                                axes_names=("undist_x", "undist_y"))

Using the example in Getting Started, add the distortion correction to the WCS pipeline and initialize the WCS.

>>> pipeline = [(detector_frame, distortion),
...             (undistorted_frame, det2sky),
...             (sky_frame, None)
...             ]
>>> wcsobj = wcs.WCS(pipeline)
>>> print(wcsobj)
       From          Transform
----------------- ----------------
         detector       distortion
undistorted_frame linear_transform
             icrs             None

Finally, save this WCS to an ASDF file:

>>> from asdf import AsdfFile
>>> tree = {"wcs": wcsobj}
>>> wcs_file = AsdfFile(tree)
>>> wcs_file.write_to("imaging_wcs_wdist.asdf")

An IFU Example - managing a discontiguous WCS

An IFU image represents the projection of several slices on a detector. Between the slices there are pixels which don’t belong to any slice. In general each slice has a unique WCS transform. There are two ways to represent this kind of transforms in GWCS depending on the way the instrument is calibrated and the available information.

Using a pixel to slice mapping

In this case a pixel map associating each pixel with a slice label (number or string) is available. The image below represents the projection of the slits of an IFU on a detector with a size (500, 1000). Slices are labeled from 1 to 6, while label 0 is reserved for pixels between the slices.

_images/ifu-regions.png

There are several models in GWCS which are useful in creating a WCS.

Given (x, y) pixel indices, LabelMapperArray returns labels (int or str) associated with these indices. RegionsSelector maps labels with transforms. It uses the LabelMapperArray to map these transforms to pixel indices.

A step by step example of constructing the WCS for an IFU with 6 slits follows.

First, import the usual packages.

>>> import numpy as np
>>> from astropy.modeling import models
>>> from astropy import coordinates as coord
>>> from astropy import units as u
>>> from gwcs import wcs, selector
>>> from gwcs import coordinate_frames as cf

The output frame is common for all slits and is a composite frame with two subframes, CelestialFrame and SpectralFrame.

>>> sky_frame = cf.CelestialFrame(name='icrs', reference_frame=coord.ICRS(), axes_order=(0, 2))
>>> spec_frame = cf.SpectralFrame(name='wave', unit=(u.micron,), axes_order=(1,), axes_names=('lambda',))
>>> cframe = cf.CompositeFrame([sky_frame, spec_frame], name='world')
>>> det = cf.Frame2D(name='detector')

All slices have the same input and output frames, however each slices has a different model transforming from pixels to world coordinates (RA, lambda, dec). For the sake of brevity this example uses a simple shift transform for each slice. Detailed examples of how to create more realistic transforms are available in Adding distortion to the imaging example.

>>> transforms = {}
>>> for i in range(1, 7):
...     transforms[i] = models.Mapping([0, 0, 1]) | models.Shift(i * 0.1) & models.Shift(i * 0.2) & models.Scale(i * 0.1)

One way to initialize LabelMapperArray is to pass it the shape of the array and the vertices of each slit on the detector {label: vertices} see :meth: from_vertices. In this example the mask is an array with the size of the detector where each item in the array corresponds to a pixel on the detector and its value is the slice number (label) this pixel belongs to.

Assuming the array is stored in ASDF format, create the mask:

Create the pixel to world transform for the entire IFU:

>>> regions_transform = selector.RegionsSelector(inputs=['x','y'],
...                                              outputs=['ra', 'dec', 'lam'],
...                                              selector=transforms,
...                                              label_mapper=mask,
...                                              undefined_transform_value=np.nan)

The WCS object now can evaluate simultaneously the transforms of all slices.

>>> wifu = wcs.WCS(forward_transform=regions_transform, output_frame=cframe, input_frame=det)
>>> y, x = mask.mapper.shape
>>> y, x = np.mgrid[:y, :x]
>>> r, d, l = wifu(x, y)

or of single slices.

The set_input() method returns the forward_transform for a specific label.

>>> wifu.forward_transform.set_input(4)(1, 2)
    (1.4, 1.8, 0.8)

Custom model storing transforms in a dictionary

In case a pixel to slice mapping is not available, one can write a custom mdoel storing transforms in a dictionary. The model would look like this:

from astropy.modeling.core import Model
from astropy.modeling.parameters import Parameter


class CustomModel(Model):

    inputs = ('label', 'x', 'y')
    outputs = ('xout', 'yout')


    def __init__(self, labels, transforms):
        super().__init__()
        self.labels = labels
        self.models = models

    def evaluate(self, label, x, y):
        index = self.labels.index(label)
        return self.models[index](x, y)

Using gwcs

Common Interface for World Coordinate System - APE 14

To improve interoperability between packages, the Astropy Project and other interested parties have collaboratively defined a standardized application programming interface (API) for world coordinate system objects to be used in Python. This API is described in the Astropy Proposal for Enhancements (APE) 14: A shared Python interface for World Coordinate Systems.

The base classes that define the low- (BaseLowLevelWCS) and high- (BaseHighLevelWCS) level APIs are in astropy. GWCS implements both APIs. Once a gWCS object is created the API methods will be available. It is recommended that applications use the Common API to ensure transparent use of GWCS and FITSWCS objects.

Using the WCS object

This section uses the imaging_wcs_wdist.asdf created in Adding distortion to the imaging example to read in a WCS object and demo its methods.

>>> import asdf
>>> asdf_file = asdf.open("imaging_wcs_wdist.asdf")
>>> wcsobj = asdf_file.tree["wcs"]
>>> print(wcsobj)    
         From          Transform
----------------- ----------------
         detector       distortion
undistorted_frame linear_transform
             icrs             None

Inspecting Available Coordinate Frames

To see what frames are defined:

>>> print(wcsobj.available_frames)
['detector', 'undistorted_frame', 'icrs']
>>> wcsobj.input_frame
<Frame2D(name="detector", unit=(Unit("pix"), Unit("pix")), axes_names=('x', 'y'), axes_order=(0, 1))>
>>> wcsobj.output_frame
<CelestialFrame(name="icrs", unit=(Unit("deg"), Unit("deg")), axes_names=('lon', 'lat'), axes_order=(0, 1), reference_frame=<ICRS Frame>)>

Because the output_frame is a CoordinateFrame object we can get the result of the WCS transform as an SkyCoord object and transform them to other standard coordinate frames supported by astropy.coordinates.

>>> skycoord = wcsobj(1, 2, with_units=True)
>>> print(skycoord)
<SkyCoord (ICRS): (ra, dec) in deg
    (5.50090023, -72.04553535)>
>>> print(skycoord.transform_to("galactic"))
<SkyCoord (Galactic): (l, b) in deg
    (306.12713109, -44.8996588)>

Using Bounding Box

The WCS object has an attribute bounding_box (default value of None) which describes the range of acceptable values for each input axis.

>>> wcsobj.bounding_box = ((0, 2048), (0, 1000))
>>> wcsobj((2,3), (1020, 980))
[array([       nan, 5.54527989]), array([         nan, -72.06454341])]

The WCS object accepts a boolean flag called with_bounding_box with default value of True. Output values which are outside the bounding_box are set to NaN. There are cases when this is not desirable and with_bounding_box=False should be passes.

Calling the footprint() returns the footprint on the sky.

>>> wcsobj.footprint()

Manipulating Transforms

Some methods allow managing the transforms in a more detailed manner.

Transforms between frames can be retrieved and evaluated separately.

>>> dist = wcsobj.get_transform('detector', 'undistorted_frame')
>>> dist(1, 2)    
(-292.4150238489997, -616.8680129899999)

Transforms in the pipeline can be replaced by new transforms.

>>> new_transform = models.Shift(1) & models.Shift(1.5) | distortion
>>> wcsobj.set_transform('detector', 'undistorted_frame', new_transform)
>>> wcsobj(1, 2)         
(5.501064280097802, -72.04557376712566)

A transform can be inserted before or after a frame in the pipeline.

>>> scale = models.Scale(2) & models.Scale(1)
>>> wcsobj.insert_transform('icrs', scale, after=False)
>>> wcsobj(1, 2)          
(11.002128560195604, -72.04557376712566)

Inverse Transformations

Often, it is useful to be able to compute inverse transformation that converts coordinates from the output frame back to the coordinates in the input frame.

In this section, for illustration purpose, we will be using the same 2D imaging WCS from imaging_wcs_wdist.asdf created in Adding distortion to the imaging example whose forward transformation converts image coordinates to world coordinates and inverse transformation converts world coordinates back to image coordinates.

>>> wcsobj = asdf.open(get_pkg_data_filename('imaging_wcs_wdist.asdf')).tree['wcs']

The most general method available for computing inverse coordinate transformation is the WCS.invert() method. This method uses automatic or user-supplied analytical inverses whenever available to convert coordinates from the output frame to the input frame. When analytical inverse is not available as is the case for the wcsobj above, a numerical solution will be attempted using WCS.numerical_inverse().

Default parameters used by WCS.numerical_inverse() or WCS.invert() methods should be acceptable in most situations:

>>> world = wcsobj(350, 200)
>>> print(wcsobj.invert(*world))  # convert a single point
(349.9999994163172, 200.00000017679295)
>>> world = wcsobj([2, 350, -5000], [2, 200, 6000])
>>> print(wcsobj.invert(*world))  # convert multiple points at once
(array([ 2.00000000e+00,  3.49999999e+02, -5.00000000e+03]), array([1.99999972e+00, 2.00000002e+02, 6.00000000e+03])

By default, parameter quiet is set to True in WCS.numerical_inverse() and so it will return results “as is” without warning us about possible loss of accuracy or about divergence of the iterative process.

In order to catch these kind of errors that can occur during numerical inversion, we need to turn off quiet mode and be prepared to catch gwcs.wcs.WCS.NoConvergence exceptions. In the next example, let’s also add a point far away from the image for which numerical inverse fails.

>>> from gwcs import NoConvergence
>>> world = wcsobj([-85000, 2, 350, 3333, -5000], [-55000, 2, 200, 1111, 6000],
...                with_bounding_box=False)
>>> try:
...     x, y = wcsobj.invert(*world, quiet=False, maxiter=40,
...                          detect_divergence=True, with_bounding_box=False)
... except NoConvergence as e:
...     print(f"Indices of diverging points: {e.divergent}")
...     print(f"Indices of poorly converging points: {e.slow_conv}")
...     print(f"Best solution:\n{e.best_solution}")
...     print(f"Achieved accuracy:\n{e.accuracy}")
Indices of diverging points: [0]
Indices of poorly converging points: [4]
Best solution:
[[ 1.38600585e+11  6.77595594e+11]
 [ 2.00000000e+00  1.99999972e+00]
 [ 3.49999999e+02  2.00000002e+02]
 [ 3.33300000e+03  1.11100000e+03]
 [-4.99999985e+03  5.99999985e+03]]
Achieved accuracy:
[[8.56497375e+02 5.09216089e+03]
 [6.57962988e-06 3.70445289e-07]
 [5.31656943e-06 2.72052603e-10]
 [6.81557583e-06 1.06560533e-06]
 [3.96365344e-04 6.41822468e-05]]

WCS User Tools

grid_from_bounding_box is a function which returns a grid of input points based on the bounding_box of the WCS.

>>> from gwcs.wcstools import grid_from_bounding_box
>>> bounding_box = ((0, 4096), (0, 2048))
>>> x, y = grid_from_bounding_box(bounding_box)
>>> ra, dec = w(x, y)

The wcstools module contains functions of general usability.

wcs_from_fiducial is a function which given a fiducial in some coordinate system, returns a WCS object.

>>> from gwcs.wcstools import wcs_from_fiducial
>>> from astropy import coordinates as coord
>>> from astropy import units as u
>>> from astropy.modeling import models
To create a WCS from a pointing on the sky, as a minimum pass a sky coordinate and a projection to the function.
>>> fiducial = coord.SkyCoord(5.46 * u.deg, -72.2 * u.deg, frame='icrs')
>>> tan = models.Pix2Sky_TAN()

Any additional transforms are prepended to the projection and sky rotation.

>>> trans = models.Shift(-2048) & models.Shift(-1024) | models.Scale(1.38*10**-5) & models.Scale(1.38*10**-5)
>>> w = wcs_from_fiducial(fiducial, projection=tan, transform=trans)
>>> w(2048, 1024)  
    (5.46, -72.2)

Listing of imaging_wcs.asdf

Listing of imaging_wcs.asdf:

#ASDF 1.0.0
#ASDF_STANDARD 1.2.0
%YAML 1.1
%TAG ! tag:stsci.edu:asdf/
--- !core/asdf-1.1.0
asdf_library: !core/software-1.0.0 {author: Space Telescope Science Institute,
homepage: 'http://github.com/spacetelescope/asdf', name: asdf, version: 2.2.0.dev1526}
history:
  extensions:
  - !core/extension_metadata-1.0.0
    extension_class: asdf.extension.BuiltinExtension
    software: {name: asdf, version: 2.2.0.dev1526}
  - !core/extension_metadata-1.0.0
    extension_class: astropy.io.misc.asdf.extension.AstropyExtension
    software: {name: astropy, version: 3.2.dev23222}
  - !core/extension_metadata-1.0.0
    extension_class: astropy.io.misc.asdf.extension.AstropyAsdfExtension
    software: {name: astropy, version: 3.2.dev23222}
  - !core/extension_metadata-1.0.0
    extension_class: gwcs.extension.GWCSExtension
    software: {name: gwcs, version: 0.10.dev417}
wcs: !<tag:stsci.edu:gwcs/wcs-1.0.0>
  name: ''
  steps:
  - !<tag:stsci.edu:gwcs/step-1.0.0>
    frame: !<tag:stsci.edu:gwcs/frame2d-1.0.0>
      axes_names: [x, y]
      name: detector
      unit: [!unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel]
    transform: !transform/compose-1.1.0
      forward:
      - !transform/remap_axes-1.1.0
        mapping: [0, 1, 0, 1]
      - !transform/concatenate-1.1.0
        forward:
        - !transform/polynomial-1.1.0
          coefficients: !core/ndarray-1.0.0
            data:
            - [0.0, 0.5, 0.6000000000000001, 0.7000000000000001, 0.8]
            - [0.1, 0.9, 1.0, 1.1, 0.0]
            - [0.2, 1.2000000000000002, 1.3, 0.0, 0.0]
            - [0.30000000000000004, 1.4000000000000001, 0.0, 0.0, 0.0]
            - [0.4, 0.0, 0.0, 0.0, 0.0]
            datatype: float64
            shape: [5, 5]
        - !transform/polynomial-1.1.0
          coefficients: !core/ndarray-1.0.0
            data:
            - [0.0, 1.0, 1.2000000000000002, 1.4000000000000001, 1.6]
            - [0.2, 1.8, 2.0, 2.2, 0.0]
            - [0.4, 2.4000000000000004, 2.6, 0.0, 0.0]
            - [0.6000000000000001, 2.8000000000000003, 0.0, 0.0, 0.0]
            - [0.8, 0.0, 0.0, 0.0, 0.0]
            datatype: float64
            shape: [5, 5]
  - !<tag:stsci.edu:gwcs/step-1.0.0>
    frame: !<tag:stsci.edu:gwcs/frame2d-1.0.0>
      axes_names: [undist_x, undist_y]
      name: undistorted_frame
      unit: [!unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel]
    transform: !transform/compose-1.1.0
      forward:
      - !transform/compose-1.1.0
        forward:
        - !transform/compose-1.1.0
          forward:
          - !transform/concatenate-1.1.0
            forward:
            - !transform/shift-1.2.0 {offset: -2048.0}
            - !transform/shift-1.2.0 {offset: -1024.0}
          - !transform/affine-1.2.0
            inverse: !transform/affine-1.2.0
              matrix: !core/ndarray-1.0.0
                data:
                - [65488.318039522, 30828.31712434267]
                - [26012.509548778366, -66838.34993781192]
                datatype: float64
                shape: [2, 2]
              translation: !core/ndarray-1.0.0
                data: [0.0, 0.0]
                datatype: float64
                shape: [2]
            matrix: !core/ndarray-1.0.0
              data:
              - [1.290551569736e-05, 5.9525007864732e-06]
              - [5.0226382102765e-06, -1.2644844123757e-05]
              datatype: float64
              shape: [2, 2]
            translation: !core/ndarray-1.0.0
              data: [0.0, 0.0]
              datatype: float64
              shape: [2]
        - !transform/gnomonic-1.1.0 {direction: pix2sky}
      - !transform/rotate3d-1.2.0 {phi: 5.63056810618, psi: 180.0, theta: -72.05457184279}
      inverse: !transform/compose-1.1.0
        forward:
        - !transform/rotate3d-1.2.0 {direction: celestial2native, phi: 5.63056810618,
          psi: 180.0, theta: -72.05457184279}
        - !transform/compose-1.1.0
          forward:
          - !transform/gnomonic-1.1.0 {direction: sky2pix}
          - !transform/compose-1.1.0
            forward:
            - !transform/affine-1.2.0
              matrix: !core/ndarray-1.0.0
                data:
                - [65488.318039522, 30828.31712434267]
                - [26012.509548778366, -66838.34993781192]
                datatype: float64
                shape: [2, 2]
              translation: !core/ndarray-1.0.0
                data: [0.0, 0.0]
                datatype: float64
                shape: [2]
            - !transform/concatenate-1.1.0
              forward:
              - !transform/shift-1.2.0 {offset: 2048.0}
              - !transform/shift-1.2.0 {offset: 1024.0}
      name: linear_transform
  - !<tag:stsci.edu:gwcs/step-1.0.0>
    frame: !<tag:stsci.edu:gwcs/celestial_frame-1.0.0>
      axes_names: [lon, lat]
      name: icrs
      reference_frame: !<tag:astropy.org:astropy/coordinates/frames/icrs-1.1.0>
        frame_attributes: {}
      unit: [!unit/unit-1.0.0 deg, !unit/unit-1.0.0 deg]
...

WCS validation

The WCS is validated when an object is read in or written to a file. However, this happens transparently to the end user and knowing the details of the validation machinery is not necessary to use or construct a WCS object.

GWCS uses the Advanced Scientific Data Format (ASDF) to validate the transforms, coordinate frames and the overall WCS object structure. ASDF makes use of abstract data type definitions called schemas. The serialization and deserialization happens in classes, referred to as tags. Most of the transform schemas live in the asdf-standard package while most of the transform tags live in astropy. GWCS Schema Definitions are available for the WCS object, coordinate frames and some WCS specific transforms.

Packages using GWCS may create their own transforms and schemas and register them as an Asdf Extension. If those are of general use, it is recommended they be included in astropy.

GWCS Schema Definitions

WCS object

wcs-1.0.0

A system for describing generalized world coordinate transformations.

Description

ASDF WCS is a way of specifying transformations (usually from detector space to world coordinate space and back) by using the transformations in the transform-schema module.

Outline

Schema DefinitionsOriginal Schema

Schema Definitions

This type is an object with the following properties:
  • name
    stringRequired
    A descriptive name for this WCS.

    No length restriction
  • steps
    arrayRequired
    A list of steps in the forward transformation from detector to world coordinates. The inverse transformation is determined automatically by reversing this list, and inverting each of the individual transforms according to the rules described in inverse.

    No length restriction
    Items in the array are restricted to the following types:

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/wcs-1.0.0"
tag: "tag:stsci.edu:gwcs/wcs-1.0.0"

title: >
  A system for describing generalized world coordinate transformations.
description: >
  ASDF WCS is a way of specifying transformations (usually from
  detector space to world coordinate space and back) by using the
  transformations in the `transform-schema` module.
type: object
properties:
  name:
    description: |
      A descriptive name for this WCS.
    type: string

  steps:
    description: |
      A list of steps in the forward transformation from detector to
      world coordinates.
      The inverse transformation is determined automatically by
      reversing this list, and inverting each of the individual
      transforms according to the rules described in
      [inverse](https://asdf-standard.readthedocs.io/en/latest/generated/stsci.edu/asdf/transform/transform-1.1.0.html#inverse).
    type: array
    items:
      $ref: step-1.0.0

required: [name, steps]
additionalProperties: true
step-1.0.0

Describes a single step of a WCS transform pipeline.

Outline

Schema DefinitionsExamplesOriginal Schema

Schema Definitions

This type is an object with the following properties:
  • frame
    objectRequired
    The frame of the inputs to the transform.

    This node must validate against any of the following:

  • transform
    object
    The transform from this step to the next one. The last step in a WCS should not have a transform, but exists only to describe the frames and units of the final output axes.

    This node must validate against any of the following:

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/step-1.0.0"
tag: "tag:stsci.edu:gwcs/step-1.0.0"

title: >
  Describes a single step of a WCS transform pipeline.
description: >
examples: []

type: object
properties:
  frame:
    description: |
      The frame of the inputs to the transform.
    anyOf:
      - type: string
      - $ref: frame-1.0.0

  transform:
    description: |
      The transform from this step to the next one.  The
      last step in a WCS should not have a transform, but
      exists only to describe the frames and units of the
      final output axes.
    anyOf:
      - $ref: "tag:stsci.edu:asdf/transform/transform-1.1.0"
      - type: 'null'
    default: null

required: [frame]

Coordinate Frames

celestial_frame-1.0.0

Represents a celestial frame.

Outline

Schema DefinitionsOriginal Schema

Schema Definitions

This node must validate against all of the following:

  • This type is an object with the following properties:
    • axes_names
      object
      This node has no type definition (unrestricted)
    • axes_order
      object
      This node has no type definition (unrestricted)
    • unit
      object
      This node has no type definition (unrestricted)

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/celestial_frame-1.0.0"
tag: "tag:stsci.edu:gwcs/celestial_frame-1.0.0"

title: >
  Represents a celestial frame.

allOf:
  - type: object
    properties:
      axes_names:
        minItems: 2
        maxItems: 3

      axes_order:
        minItems: 2
        maxItems: 3

      unit:
        minItems: 2
        maxItems: 3

  - $ref: frame-1.0.0
frame-1.0.0

The base class of all coordinate frames.

Description

These objects are designed to be nested in arbitrary ways to build up transformation pipelines out of a number of low-level pieces.

Outline

Schema DefinitionsExamplesOriginal Schema

Schema Definitions

This type is an object with the following properties:
  • name
    stringRequired
    A user-friendly name for the frame.

    No length restriction
  • axes_order
    array
    The order of the axes.

    No length restriction
    Items in the array are restricted to the following types:
    integer
  • axes_names
    array
    The name of each axis in this frame.

    No length restriction

    Items in the array must be any of the following types:

    • string
      No length restriction
    • null
  • reference_frame
    tag:astropy.org:astropy/coordinates/frames/baseframe-1.0.0
    The reference frame.

  • unit
    array
    Units for each axis.

    No length restriction
    Items in the array are restricted to the following types:
  • axis_physical_types
    array
    An iterable of strings describing the physical type for each world axis. These should be names from the VO UCD1+ controlled Vocabulary (http://www.ivoa.net/documents/latest/UCDlist.html).

    No length restriction
    Items in the array are restricted to the following types:
    string
    No length restriction

Examples

A celestial frame in the ICRS reference frame. :

!<tag:stsci.edu:gwcs/celestial_frame-1.0.0>
  axes_names: [lon, lat]
  name: CelestialFrame
  reference_frame: !<tag:astropy.org:astropy/coordinates/frames/icrs-1.1.0>
    frame_attributes: {}
  unit: [!unit/unit-1.0.0 deg, !unit/unit-1.0.0 deg]

A pixel frame in three dimensions :

!<tag:stsci.edu:gwcs/frame-1.0.0>
  axes_names: [raster position, slit position, wavelength]
  axes_order: [0, 1, 2]
  axes_type: [SPATIAL, SPATIAL, SPECTRAL]
  name: pixel
  naxes: 3
  unit: [!unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel]

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/frame-1.0.0"
title: |
  The base class of all coordinate frames.

description: |
  These objects are designed to be nested in arbitrary ways to build up
  transformation pipelines out of a number of low-level pieces.

examples:
  -
    - |
        A celestial frame in the ICRS reference frame.
    - |
        !<tag:stsci.edu:gwcs/celestial_frame-1.0.0>
          axes_names: [lon, lat]
          name: CelestialFrame
          reference_frame: !<tag:astropy.org:astropy/coordinates/frames/icrs-1.1.0>
            frame_attributes: {}
          unit: [!unit/unit-1.0.0 deg, !unit/unit-1.0.0 deg]

  -
    - |
        A pixel frame in three dimensions
    - |
        !<tag:stsci.edu:gwcs/frame-1.0.0>
          axes_names: [raster position, slit position, wavelength]
          axes_order: [0, 1, 2]
          axes_type: [SPATIAL, SPATIAL, SPECTRAL]
          name: pixel
          naxes: 3
          unit: [!unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel]

type: object
properties:
  name:
    description: |
      A user-friendly name for the frame.
    type: string

  axes_order:
    description: |
      The order of the axes.
    type: array
    items:
      type: integer

  axes_names:
    description: |
      The name of each axis in this frame.
    type: array
    items:
      anyOf:
        - type: string
        - type: 'null'

  reference_frame:
    description: |
      The reference frame.
    $ref: "tag:astropy.org:astropy/coordinates/frames/baseframe-1.0.0"

  unit:
    description: |
      Units for each axis.
    type: array
    items:
      $ref: "tag:stsci.edu:asdf/unit/unit-1.0.0"

  axis_physical_types:
    description: |
      An iterable of strings describing the physical type for each world axis.
      These should be names from the VO UCD1+ controlled Vocabulary
      (http://www.ivoa.net/documents/latest/UCDlist.html).
    type: array
    items:
      type:
        string
    
required: [name]
additionalProperties: true
spectral_frame-1.0.0

Represents a spectral frame.

Outline

Schema DefinitionsOriginal Schema

Schema Definitions

This node must validate against all of the following:

  • This type is an object with the following properties:
    • reference_position
      object
      The position of the reference frame.

      This node has no type definition (unrestricted)
    • axes_names
      object
      This node has no type definition (unrestricted)
    • axes_order
      object
      This node has no type definition (unrestricted)
    • unit
      object
      This node has no type definition (unrestricted)

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/spectral_frame-1.0.0"
tag: "tag:stsci.edu:gwcs/spectral_frame-1.0.0"

title: >
  Represents a spectral frame.
allOf:
  - type: object
    properties:
      reference_position:
        description: |
          The position of the reference frame.
        enum: [geocenter, barycenter, heliocenter]
        default: geocenter

      axes_names:
        minItems: 1
        maxItems: 1

      axes_order:
        minItems: 1
        maxItems: 1

      unit:
        minItems: 1
        maxItems: 1

  - $ref: frame-1.0.0
frame2d-1.0.0

Represents a 2D frame.

Outline

Schema DefinitionsExamplesOriginal Schema

Schema Definitions

This node must validate against all of the following:

  • This type is an object with the following properties:
    • axes_names
      object
      This node has no type definition (unrestricted)
    • axes_order
      object
      This node has no type definition (unrestricted)
    • unit
      object
      This node has no type definition (unrestricted)

Examples

A two dimensional spatial frame :

!<tag:stsci.edu:gwcs/frame2d-1.0.0>
  axes_names: [lon, lat]
  name: Frame2D
  unit: [!unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel]

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/frame2d-1.0.0"
tag: "tag:stsci.edu:gwcs/frame2d-1.0.0"

title: >
  Represents a 2D frame.

examples:
  -
    - |
        A two dimensional spatial frame
    - |
        !<tag:stsci.edu:gwcs/frame2d-1.0.0>
          axes_names: [lon, lat]
          name: Frame2D
          unit: [!unit/unit-1.0.0 pixel, !unit/unit-1.0.0 pixel]


allOf:
  - type: object
    properties:
      axes_names:
        minItems: 2
        maxItems: 2

      axes_order:
        minItems: 2
        maxItems: 2

      unit:
        minItems: 2
        maxItems: 2

  - $ref: frame-1.0.0
temporal_frame-1.0.0

Represents a temporal frame.

Outline

Schema DefinitionsOriginal Schema

Schema Definitions

This type is an object with the following properties:
  • name
    stringRequired
    A user-friendly name for the frame.

    No length restriction
  • axes_order
    array
    The order of the axes.

    No length restriction
    Items in the array are restricted to the following types:
    integer
  • axes_names
    array
    The name of each axis in this frame.

    No length restriction

    Items in the array must be any of the following types:

    • string
      No length restriction
    • null
  • reference_frame
    tag:stsci.edu:asdf/time/time-1.1.0
    The reference frame.

  • unit
    array
    Units for each axis.

    No length restriction
    Items in the array are restricted to the following types:
  • axis_physical_types
    array
    An iterable of strings describing the physical type for each world axis. These should be names from the VO UCD1+ controlled Vocabulary (http://www.ivoa.net/documents/latest/UCDlist.html).

    No length restriction
    Items in the array are restricted to the following types:
    string
    No length restriction

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/temporal_frame-1.0.0"
tag: "tag:stsci.edu:gwcs/temporal_frame-1.0.0"

title: >
  Represents a temporal frame.


type: object
properties:
  name:
    description: |
      A user-friendly name for the frame.
    type: string

  axes_order:
    description: |
      The order of the axes.
    type: array
    items:
      type: integer

  axes_names:
    description: |
      The name of each axis in this frame.
    type: array
    items:
      anyOf:
        - type: string
        - type: 'null'

  reference_frame:
    description: |
      The reference frame.
    $ref: "tag:stsci.edu:asdf/time/time-1.1.0"

  unit:
    description: |
      Units for each axis.
    type: array
    items:
      $ref: "tag:stsci.edu:asdf/unit/unit-1.0.0"

  axis_physical_types:
    description: |
      An iterable of strings describing the physical type for each world axis.
      These should be names from the VO UCD1+ controlled Vocabulary
      (http://www.ivoa.net/documents/latest/UCDlist.html).
    type: array
    items:
      type:
        string

required: [name]
composite_frame-1.0.0

Represents a set of frames.

Outline

Schema DefinitionsOriginal Schema

Schema Definitions

This node must validate against all of the following:

  • This type is an object with the following properties:
    • name
      string
      Name of composite frame.

      No length restriction
    • frames
      array
      List of frames in the composite frame.

      No length restriction

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/composite_frame-1.0.0"
tag: "tag:stsci.edu:gwcs/composite_frame-1.0.0"

title: >
  Represents a set of frames.
allOf:
  - type: object
    properties:
      name:
        description:
          Name of composite frame.
        type: string

      frames:
        description:
          List of frames in the composite frame.
        type: array

Transforms

label_mapper-1.0.0

Represents a mapping from a coordinate value to a label.

Description

A label mapper instance maps inputs to a label. It is used together with regions_selector. The label_mapper returns the label corresponding to given inputs. The regions_selector returns the transform corresponding to this label. This maps inputs (e.g. pixels on a detector) to transforms uniquely.

Outline

Schema DefinitionsExamplesOriginal Schema

Schema Definitions

This node must validate against all of the following:

  • This type is an object with the following properties:
    • mapper
      objectRequired
      A mapping of inputs to labels. In the general case this is a astropy.modeling.core.Model.

      It could be a numpy array with the shape of the detector/observation. Pixel values are of type integer or string and represent region labels. Pixels which are not within any region have value no_label.

      It could be a dictionary which maps tuples to labels or floating point numbers to labels.

      This node must validate against any of the following:

    • inputs
      array
      Names of inputs.

      No length restriction
      Items in the array are restricted to the following types:
      string
      No length restriction
    • inputs_mapping
      tag:stsci.edu:asdf/transform/transform-1.1.0
      mapping

    • atol
      number
      absolute tolerance to compare keys in mapper.

    • no_label
      object
      Fill in value for missing output.

      This node must validate against any of the following:

      • number
      • string
        No length restriction

Examples

Map array indices are to labels.:

!<tag:stsci.edu:gwcs/label_mapper-1.0.0>
  mapper: !core/ndarray-1.0.0
    data:
    - [1, 0, 2]
    - [1, 0, 2]
    - [1, 0, 2]
    datatype: int64
    shape: [3, 3]
    no_label: 0

Map numbers dictionary to transforms which return labels.:

!<tag:stsci.edu:gwcs/label_mapper-1.0.0>
  atol: 1.0e-08
  inputs: [x, y]
  inputs_mapping: !transform/remap_axes-1.1.0
      mapping: [0]
      n_inputs: 2
  mapper: !!omap
    - !!omap
      labels: [-1.67833272, -1.9580548, -1.118888]
    - !!omap
      models:
      - !transform/shift-1.1.0 {offset: 6.0}
      - !transform/shift-1.1.0 {offset: 2.0}
      - !transform/shift-1.1.0 {offset: 4.0}
  no_label: 0

Map a number within a range of numbers to transforms which return labels.:

!<tag:stsci.edu:gwcs/label_mapper-1.0.0>
  mapper: !!omap
  - !!omap
    labels:
    - [3.2, 4.1]
    - [2.67, 2.98]
    - [1.95, 2.3]
  - !!omap
    models:
    - !transform/shift-1.1.0 {offset: 6.0}
    - !transform/shift-1.1.0 {offset: 2.0}
    - !transform/shift-1.1.0 {offset: 4.0}
  inputs: [x, y]
  inputs_mapping: !transform/remap_axes-1.1.0
    mapping: [0]
    n_inputs: 2

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/label_mapper-1.0.0"
tag: "tag:stsci.edu:gwcs/label_mapper-1.0.0"
title: >
  Represents a mapping from a coordinate value to a label.
description: |
  A label mapper instance maps inputs to a label.  It is used together
  with
  [regions_selector](ref:regions_selector-1.0.0). The
  [label_mapper](ref:label_mapper-1.0.0)
  returns the label corresponding to given inputs. The
  [regions_selector](ref:regions_selector-1.0.0)
  returns the transform corresponding to this label. This maps inputs
  (e.g. pixels on a detector) to transforms uniquely.

examples:
  -
    - Map array indices are to labels.

    - |
        !<tag:stsci.edu:gwcs/label_mapper-1.0.0>
          mapper: !core/ndarray-1.0.0
            data:
            - [1, 0, 2]
            - [1, 0, 2]
            - [1, 0, 2]
            datatype: int64
            shape: [3, 3]
            no_label: 0

  -
    - Map numbers dictionary to transforms which return labels.

    - |
        !<tag:stsci.edu:gwcs/label_mapper-1.0.0>
          atol: 1.0e-08
          inputs: [x, y]
          inputs_mapping: !transform/remap_axes-1.1.0
              mapping: [0]
              n_inputs: 2
          mapper: !!omap
            - !!omap
              labels: [-1.67833272, -1.9580548, -1.118888]
            - !!omap
              models:
              - !transform/shift-1.1.0 {offset: 6.0}
              - !transform/shift-1.1.0 {offset: 2.0}
              - !transform/shift-1.1.0 {offset: 4.0}
          no_label: 0

  -
    - Map a number within a range of numbers to transforms which return labels.

    - |
        !<tag:stsci.edu:gwcs/label_mapper-1.0.0>
          mapper: !!omap
          - !!omap
            labels:
            - [3.2, 4.1]
            - [2.67, 2.98]
            - [1.95, 2.3]
          - !!omap
            models:
            - !transform/shift-1.1.0 {offset: 6.0}
            - !transform/shift-1.1.0 {offset: 2.0}
            - !transform/shift-1.1.0 {offset: 4.0}
          inputs: [x, y]
          inputs_mapping: !transform/remap_axes-1.1.0
            mapping: [0]
            n_inputs: 2

allOf:
  - $ref: "tag:stsci.edu:asdf/transform/transform-1.1.0"
  - type: object
    properties:
      mapper:
        description: |
          A mapping of inputs to labels.
          In the general case this is a `astropy.modeling.core.Model`.

          It could be a numpy array with the shape of the detector/observation.
          Pixel values are of type integer or string and represent
          region labels. Pixels which are not within any region have value ``no_label``.

          It could be a dictionary which maps tuples to labels or floating point numbers to labels.

        anyOf:
          - $ref: "tag:stsci.edu:asdf/core/ndarray-1.0.0"
          - $ref: "tag:stsci.edu:asdf/transform/transform-1.1.0"
          - type: object
            properties:
              labels:
                type: array
                items:
                  anyOf:
                    - type: number
                    - type: array
                      items:
                        type: number
                      minLength: 2
                      maxLength: 2
              models:
                type: array
                items:
                  $ref: "tag:stsci.edu:asdf/transform/transform-1.1.0"

      inputs:
        type: array
        items:
          type: string
        description: |
          Names of inputs.
      inputs_mapping:
        $ref: "tag:stsci.edu:asdf/transform/transform-1.1.0"
        description: |
          [mapping](https://asdf-standard.readthedocs.io/en/latest/generated/stsci.edu/asdf/transform/remap_axes-1.1.0.html)
      atol:
        type: number
        description: |
          absolute tolerance to compare keys in mapper.
      no_label:
        description: |
          Fill in value for missing output.
        anyOf:
          - type: number
          - type: string

    required: [mapper]
regions_selector-1.0.0

Represents a discontinuous transform.

Description

Maps regions to transgorms and evaluates the transforms with the corresponding inputs.

Outline

Schema DefinitionsExamplesOriginal Schema

Schema Definitions

This node must validate against all of the following:

  • This type is an object with the following properties:
    • label_mapper
      ./label_mapper-1.0.0Required
      An instance of label_mapper-1.1.0

    • inputs
      arrayRequired
      Names of inputs.

      No length restriction
      Items in the array are restricted to the following types:
      string
      No length restriction
    • outputs
      arrayRequired
      Names of outputs.

      No length restriction
      Items in the array are restricted to the following types:
      string
      No length restriction
    • selector
      objectRequired
      A mapping of regions to trransforms.

      This type is an object with the following properties:
      • labels
        array
        An array of unique region labels.

        No length restriction
        Items in the array are restricted to the following types:
        [‘integer’, ‘string’]
      • transforms
        array
        A transform for each region. The order should match the order of labels.

        No length restriction
        Items in the array are restricted to the following types:
    • undefined_transform_value
      number
      Value to be returned if there’s no transform defined for the inputs.

Examples

Create a regions_selector schema for 2 regions, labeled “1” and “2”.:

!<tag:stsci.edu:gwcs/regions_selector-1.0.0>
  inputs: [x, y]
  label_mapper: !<tag:stsci.edu:gwcs/label_mapper-1.0.0>
    mapper: !core/ndarray-1.0.0
      datatype: int8
      data:
      - [0, 1, 1, 0, 2, 0]
      - [0, 1, 1, 0, 2, 0]
      - [0, 1, 1, 0, 2, 0]
      - [0, 1, 1, 0, 2, 0]
      - [0, 1, 1, 0, 2, 0]
      datatype: int64
      shape: [5, 6]
  no_label: 0
  outputs: [ra, dec, lam]
  selector: !!omap
  - !!omap
    labels: [1, 2]
  - !!omap
    transforms:
    - !transform/compose-1.1.0
      forward:
      - !transform/remap_axes-1.1.0
        mapping: [0, 1, 1]
      - !transform/concatenate-1.1.0
        forward:
        - !transform/concatenate-1.1.0
          forward:
          - !transform/shift-1.1.0 {offset: 1.0}
          - !transform/shift-1.1.0 {offset: 2.0}
        - !transform/shift-1.1.0 {offset: 3.0}
    - !transform/compose-1.1.0
      forward:
      - !transform/remap_axes-1.1.0
        mapping: [0, 1, 1]
      - !transform/concatenate-1.1.0
        forward:
        - !transform/concatenate-1.1.0
          forward:
          - !transform/scale-1.1.0 {factor: 2.0}
          - !transform/scale-1.1.0 {factor: 3.0}
        - !transform/scale-1.1.0 {factor: 3.0}
  undefined_transform_value: .nan

Original Schema 

%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://stsci.edu/schemas/gwcs/regions_selector-1.0.0"
tag: "tag:stsci.edu:gwcs/regions_selector-1.0.0"
title: >
  Represents a discontinuous transform.
description: |
  Maps regions to transgorms and evaluates the transforms with the corresponding inputs.

examples:
  -
    - Create a regions_selector schema for 2 regions, labeled "1" and "2".
    - |
        !<tag:stsci.edu:gwcs/regions_selector-1.0.0>
          inputs: [x, y]
          label_mapper: !<tag:stsci.edu:gwcs/label_mapper-1.0.0>
            mapper: !core/ndarray-1.0.0
              datatype: int8
              data:
              - [0, 1, 1, 0, 2, 0]
              - [0, 1, 1, 0, 2, 0]
              - [0, 1, 1, 0, 2, 0]
              - [0, 1, 1, 0, 2, 0]
              - [0, 1, 1, 0, 2, 0]
              datatype: int64
              shape: [5, 6]
          no_label: 0
          outputs: [ra, dec, lam]
          selector: !!omap
          - !!omap
            labels: [1, 2]
          - !!omap
            transforms:
            - !transform/compose-1.1.0
              forward:
              - !transform/remap_axes-1.1.0
                mapping: [0, 1, 1]
              - !transform/concatenate-1.1.0
                forward:
                - !transform/concatenate-1.1.0
                  forward:
                  - !transform/shift-1.1.0 {offset: 1.0}
                  - !transform/shift-1.1.0 {offset: 2.0}
                - !transform/shift-1.1.0 {offset: 3.0}
            - !transform/compose-1.1.0
              forward:
              - !transform/remap_axes-1.1.0
                mapping: [0, 1, 1]
              - !transform/concatenate-1.1.0
                forward:
                - !transform/concatenate-1.1.0
                  forward:
                  - !transform/scale-1.1.0 {factor: 2.0}
                  - !transform/scale-1.1.0 {factor: 3.0}
                - !transform/scale-1.1.0 {factor: 3.0}
          undefined_transform_value: .nan


allOf:
  - $ref: "tag:stsci.edu:asdf/transform/transform-1.1.0"
  - type: object
    properties:
      label_mapper:
        description: |
          An instance of
          [label_mapper-1.1.0](ref:label_mapper-1.0.0)
        $ref: "./label_mapper-1.0.0"
      inputs:
        description: |
          Names of inputs.
        type: array
        items:
          type: string
      outputs:
        description: |
          Names of outputs.
        type: array
        items:
          type: string
      selector:
        description: |
          A mapping of regions to trransforms.
        type: object
        properties:
          labels:
            description: |
              An array of unique region labels.
            type: array
            items:
              type:
                - integer
                - string
          transforms:
            description: |
              A transform for each region. The order should match the order of labels.
            type: array
            items:
              $ref: "tag:stsci.edu:asdf/transform/transform-1.1.0"
      undefined_transform_value:
        description: |
          Value to be returned if there's no transform defined for the inputs.
        type: number
    required: [label_mapper, inputs, outputs, selector]

Fitting a WCS to input pixels & sky positions

Suppose we have an image where we have centroid positions for a number of sources, and we have matched these positions to an external catalog to obtain (RA, Dec). If this data is missing or has inaccurate WCS information, it is useful to fit or re-fit a GWCS object with this matched list of coordinate pairs to be able to transform between pixel and sky.

This example shows how to use the wcs_from_points tool to fit a WCS to a matched set of pixel and sky positions. Along with arrays of the (x,y) pixel position in the image and the matched sky coordinates, the fiducial point for the projection must be supplied as a SkyCoord object. Additionally, the projection type must be specified from the available projections in projcode.

Geometric distortion can also be fit to the input coordinates - the distortion type (2D polynomial, chebyshev, legendre) and the degree can be supplied to fit this component of the model.

The following example will show how to fit a WCS, including a 4th degree 2D polynomial, to a set of input pixel positions of sources in an image and their corresponding positions on the sky obtained from a catalog.

Import the wcs_from_points function,

>>> from gwcs.wcstools import wcs_from_points

along with some useful general imports.

>>> from astropy.coordinates import SkyCoord
>>> from astropy.io import ascii
>>> import astropy.units as u
>>> import numpy as np

A collection of 20 matched coordinate pairs in x, y, RA, and Dec stored in two arrays, will be used to fit the WCS information. The function requires tuples of arrays.

>>> xy = (np.array([2810.156, 2810.156,  650.236, 1820.927, 3425.779, 2750.369,
...                 212.422, 1146.91 ,   27.055, 2100.888,  648.149,   22.212,
...                 2003.314,  727.098,  248.91 ,  409.998, 1986.931,  128.925,
...                 1106.654, 1502.67 ]),
...       np.array([1670.347, 1670.347,  360.325,  165.663,  900.922,  700.148,
...                 1416.235, 1372.364,  398.823,  580.316,  317.952,  733.984,
...                 339.024,  234.29 , 1241.608,  293.545, 1794.522, 1365.706,
...                 583.135,   25.306]))
>>> radec = (np.array([246.75001315, 246.75001315, 246.72033646, 246.72303144,
...                    246.74164072, 246.73540614, 246.73379121, 246.73761455,
...                    246.7179495 , 246.73051123, 246.71970072, 246.7228646 ,
...                    246.72647213, 246.7188386 , 246.7314031 , 246.71821002,
...                    246.74785534, 246.73265223, 246.72579817, 246.71943263]),
...          np.array([43.48690547,  43.48690547,  43.46792989,  43.48075238,
...                    43.49560501,  43.48903538,  43.46045875,  43.47030776,
...                    43.46132376,  43.48252763,  43.46802566,  43.46035331,
...                    43.48218262,  43.46908299,  43.46131665,  43.46560591,
...                    43.47791234,  43.45973025,  43.47208325,  43.47779988]))

We can now choose the reference point on the sky for the projection. This is passed in as a SkyCoord object so that information about the celestial frame and units is given as well. The input world coordinates are passed in as unitless arrays, and so are assumed to be of the same unit and frame as the fiducial point.

>>> proj_point = SkyCoord(246.7368408, 43.480712949, frame = 'icrs', unit = (u.deg,u.deg))

We can now call the function that returns a GWCS object corresponding to the best fit parameters that relate the input pixels and sky coordinates with a TAN projection centered at the reference point we specified, with a distortion model (degree 4 polynomial). This function will return a GWCS object that can be used to transform between coordinate frames.

>>> gwcs_obj = wcs_from_points(xy, radec, proj_point)

This GWCS object contains parameters for a TAN projection, rotation, scale, skew and a polynomial fit to x and y that represent the best-fit to the input coordinates. With WCS information associated with the data now, we can easily work in both pixel and sky space, and transform between frames.

The GWCS object, which by default when called executes for forward transformation, can be used to convert coordinates from pixel to world.

>>> gwcs_obj(36.235,642.215)    
(246.72158004206716, 43.46075091731673)
Or equivalently
>>> gwcs_obj.forward_transform(36.235,642.215)  
(246.72158004206716, 43.46075091731673)

See also

Reference/API

gwcs.wcs Module

Classes

WCS([forward_transform, input_frame, ...])

Basic WCS class.

NoConvergence(*args[, best_solution, ...])

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

Class Inheritance Diagram

Inheritance diagram of gwcs.wcs.WCS, gwcs.wcs.NoConvergence

gwcs.coordinate_frames Module

Defines coordinate frames and ties them to data axes.

Classes

Frame2D([axes_order, unit, axes_names, ...])

A 2D coordinate frame.

CelestialFrame([axes_order, ...])

Celestial Frame Representation

SpectralFrame([axes_order, reference_frame, ...])

Represents Spectral Frame

CompositeFrame(frames[, name])

Represents one or more frames.

CoordinateFrame(naxes, axes_type, axes_order)

Base class for Coordinate Frames.

TemporalFrame(reference_frame[, unit, ...])

A coordinate frame for time axes.

Class Inheritance Diagram

Inheritance diagram of gwcs.coordinate_frames.Frame2D, gwcs.coordinate_frames.CelestialFrame, gwcs.coordinate_frames.SpectralFrame, gwcs.coordinate_frames.CompositeFrame, gwcs.coordinate_frames.CoordinateFrame, gwcs.coordinate_frames.TemporalFrame

gwcs.wcstools Module

Functions

wcs_from_fiducial(fiducial[, ...])

Create a WCS object from a fiducial point in a coordinate frame.

grid_from_bounding_box(bounding_box[, step, ...])

Create a grid of input points from the WCS bounding_box.

wcs_from_points(xy, world_coordinates, fiducial)

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

gwcs.selector Module

The classes in this module create discontinuous transforms.

The main class is RegionsSelector. It maps inputs to transforms and evaluates the transforms on the corresponding inputs. Regions are well defined spaces in the same frame as the inputs. Regions are assigned unique labels (int or str). The region labels are used as a proxy between inputs and transforms. An example is the location of IFU slices in the detector frame.

RegionsSelector uses two structures:

  • A mapping of inputs to labels - “label_mapper”

  • A mapping of labels to transforms - “transform_selector”

A “label_mapper” is also a transform, a subclass of astropy.modeling.core.Model, which returns the labels corresponding to the inputs.

An instance of a LabelMapper class is passed to RegionsSelector. The labels are used by RegionsSelector to match inputs to transforms. Finally, RegionsSelector evaluates the transforms on the corresponding inputs. Label mappers and transforms take the same inputs as RegionsSelector. The inputs should be filtered appropriately using the inputs_mapping argument which is ian instance of Mapping. The transforms in “transform_selector” should have the same number of inputs and outputs.

This is illustrated below using two regions, labeled 1 and 2

+-----------+
| +-+       |
| | |  +-+  |
| |1|  |2|  |
| | |  +-+  |
| +-+       |
+-----------+

                   +--------------+
                   | label mapper |
                   +--------------+
                     ^       |
                     |       V
           ----------|   +-------+
           |             | label |
         +--------+      +-------+
--->     | inputs |          |
         +--------+          V
              |          +--------------------+
              |          | transform_selector |
              |          +--------------------+
              V                  |
         +-----------+           |
         | transform |<-----------
         +------------+
              |
              V
         +---------+
         | outputs |
         +---------+

The base class _LabelMapper can be subclassed to create other label mappers.

Classes

LabelMapperArray(mapper[, inputs_mapping, name])

Maps array locations to labels.

LabelMapperDict(inputs, mapper[, ...])

Maps a number to a transform, which when evaluated returns a label.

LabelMapperRange(inputs, mapper[, ...])

The structure this class uses maps a range of values to a transform.

RegionsSelector(inputs, outputs, selector, ...)

This model defines discontinuous transforms.

LabelMapper(inputs, mapper[, no_label, ...])

Maps inputs to regions.

Class Inheritance Diagram

Inheritance diagram of gwcs.selector.LabelMapperArray, gwcs.selector.LabelMapperDict, gwcs.selector.LabelMapperRange, gwcs.selector.RegionsSelector, gwcs.selector.LabelMapper

gwcs.spectroscopy Module

Spectroscopy related models.

Classes

WavelengthFromGratingEquation(...)

Solve the Grating Dispersion Law for the wavelength.

AnglesFromGratingEquation3D(groove_density, ...)

Solve the 3D Grating Dispersion Law in Direction Cosine space for the refracted angle.

Snell3D(**kwargs)

Snell model in 3D form.

SellmeierGlass(B_coef, C_coef, **kwargs)

Sellmeier equation for glass.

SellmeierZemax([temperature, ...])

Sellmeier equation used by Zemax.

Class Inheritance Diagram

Inheritance diagram of gwcs.spectroscopy.WavelengthFromGratingEquation, gwcs.spectroscopy.AnglesFromGratingEquation3D, gwcs.spectroscopy.Snell3D, gwcs.spectroscopy.SellmeierGlass, gwcs.spectroscopy.SellmeierZemax

gwcs.geometry Module

Models for general analytical geometry transformations.

Classes

ToDirectionCosines(**kwargs)

Transform a vector to direction cosines.

FromDirectionCosines(**kwargs)

Transform directional cosines to vector.

SphericalToCartesian([wrap_lon_at])

Convert spherical coordinates on a unit sphere to cartesian coordinates.

CartesianToSpherical([wrap_lon_at])

Convert cartesian coordinates to spherical coordinates on a unit sphere.

Class Inheritance Diagram

Inheritance diagram of gwcs.geometry.ToDirectionCosines, gwcs.geometry.FromDirectionCosines, gwcs.geometry.SphericalToCartesian, gwcs.geometry.CartesianToSpherical