Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

SCSB-111: Improve documentation #483

Merged
merged 10 commits into from
Dec 18, 2023
Merged

SCSB-111: Improve documentation #483

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
43d6e9d
fix mask issue for IFU example
perrygreenfield Nov 9, 2023
9f8c216
remove asdf in fits examples
perrygreenfield Nov 9, 2023
5ee8d21
remove asdf in fits examples
perrygreenfield Nov 9, 2023
ad5f658
final changes
perrygreenfield Nov 14, 2023
275707f
update CHANGE.rst with PR number
perrygreenfield Nov 14, 2023
f205bc5
add missing doctest-skip directives
perrygreenfield Nov 16, 2023
c72f876
address review comments
perrygreenfield Nov 22, 2023
b957e24
Merge branch 'master' into improvedocs
nden Nov 28, 2023
d51a0e0
Merge branch 'master' into improvedocs
nden Dec 14, 2023
fb48cac
move change log; add text
nden Dec 18, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 2 additions & 3 deletions CHANGES.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,12 +1,11 @@
0.21.0 (unreleased)
-------------------

- Improve documentation (part 1) [#483]
nden marked this conversation as resolved.
Show resolved Hide resolved

0.20.0 (2023-11-29)
-------------------

other
^^^^^

- Replace ``pkg_resources`` with ``importlib.metadata``. [#478]

- Serialize and deserialize ``pixel_shape`` with asdf. [#480]
Expand Down
102 changes: 102 additions & 0 deletions docs/gwcs/fits_analog.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,102 @@
.. _fits_equivalent_example:

FITS Equivalent WCS Example
===========================

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


For this example the following imports are needed:

>>> 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 `~gwcs.coordinate_frames.CoordinateFrame`. Although strings are
acceptable as ``coordinate_frames`` it is recommended this is used only in testing/debugging.
nden marked this conversation as resolved.
Show resolved Hide resolved

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

nden marked this conversation as resolved.
Show resolved Hide resolved
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.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why do we need this detector frame? Can't we just bunch all transforms together from detector all the way to the sky? Also, what is a "detector frame"? is it just the value of the name argument in the Frame2D initializer?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Generally the frames are information containers - they hold the units on the axes, the names of the axes, the axes type and physical type. You could create a WCS without an input frame but some functionality will not be available. For example, if the transforms are defined with units, you will need to make sure the inputs are with units, so any time you pass input in the case above, the inputs need to have u.pix attached to them. This can be annoying.
In addition the input frame does not always represent a detector and imo it's good to have the information on what inputs are expected.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nden comment should be added to the docs. My comment was intended to prompt a discussion/explanation about organization of the (gwcs) pipelines. For example, why do some pipelines have multiple steps (JWST WCS comes to mind)?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I wonder if a paragraphs or two before this example outlining the general structure of a WCS object would be helpful. Any objection @nden?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What if the current "Getting Started" is replaced by an earlier section:

Basic Structure of a GWCS Object
--------------------------------

The key concept to be aware of is that a GWCS Object consists of a pipeline
of steps; each step contains a transform (i.e., an Astropy model) that 
converts the input coordinates of the step to the output coordinates of
the step. Furthermore, each step has an optional coordinate frame associated
with the step. The coordinate frame represents the input coordinate frame, not
the output coordinates. Most typically, the first step coordinate frame is
the detector pixel coordinates (the default). Since no step has a coordinate 
frame for the output coordinates, it is necessary to append a step with no
transform to the end of the pipeline to represent the output coordinate frame.
For imaging, this frame typically references one of the Astropy standard
Sky Coordinate Frames of Reference. The GWCS frames also serve to hold the
units on the axes, the names of the axes and the physical type of the axis
(e.g., wavelength).

Since it is often useful to obtain coordinates in an intermediate frame of
reference, GWCS allows the pipeline to consist of more than one transform.
For example, for spectrographs, it is useful to have access to coordinates
in the slit plane, and in such a case, the first step would transform from
the detector to the slit plane, and the second step from the slit plane to
sky coordinates and a wavelength. Constructed this way, it is possible to
extract from the GWCS the needed transforms between identified frames of
reference.  

The GWCS object can be saved to the ASDF format using the
`asdf <https://asdf.readthedocs.io/en/latest/>`__ package and validated 
using `ASDF Standard <https://asdf-standard.readthedocs.io/en/latest/>`__

There are two ways to save the GWCS object to a files:

- `Save a WCS object as a pure ASDF file`_

@mcara does this address your comment?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes! This is a great introduction.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I like this too.

Furthermore, each step has an optional coordinate frame associated with the step.

Is the frame optional? I thought it's required.

a GWCS Object consists of a pipeline of steps

Suggest adding linear -> consists of a linear pipeline of steps
so that it's clear it's not an arbitrary graph.


>>> 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

Now we have a complete WCS object. The next example will use it to convert pixel
coordinates(1, 2) to sky coordinates:

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

The :meth:`~gwcs.wcs.WCS.invert` method evaluates the :meth:`~gwcs.wcs.WCS.backward_transform` to provide a mapping from sky coordinates to pixel coordinates
if available, otherwise it applies an iterative method to calculate the pixel coordinates.

>>> wcsobj.invert(sky)
(<Quantity 1. pix>, <Quantity 2. pix>)
Binary file modified docs/gwcs/ifu-regions.png
View file
Open in desktop
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
26 changes: 14 additions & 12 deletions docs/gwcs/ifu.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,18 @@ First, import the usual packages.
>>> from gwcs import wcs, selector
>>> from gwcs import coordinate_frames as cf

Next, create the appropriate mapper object corresponding to the figure above:

>>> # Ignore the details of how this mask is constructed; they are using
>>> # array operations to generate the mask displayed for this example.
>>> y, x = np.mgrid[:1000, :500]
>>> fmask = (((-x + 0.01 * y + 0.00002 * y**2)/ 500) * 13 - 0.5) + 14
>>> mask = fmask.astype(np.int8)
>>> mask[(mask % 2) == 1] = 0
>>> mask[mask > 13] = 0
>>> mask = mask // 2
>>> labelmapper = selector.LabelMapperArray(mask)

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

Expand All @@ -59,28 +71,18 @@ In this example the mask is an array with the size of the detector where each it
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 <https://asdf-standard.readthedocs.io/en/latest>`__ format, create the mask:

.. doctest-skip-all

>>> import asdf
>>> f = asdf.open('mask.asdf')
>>> data = f.tree['mask']
>>> mask = selector.LabelMapperArray(data)

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,
... label_mapper=labelmapper,
... 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 = labelmapper.mapper.shape
>>> y, x = np.mgrid[:y, :x]
>>> r, d, l = wifu(x, y)

Expand Down
4 changes: 4 additions & 0 deletions docs/gwcs/using_wcs.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -109,6 +109,9 @@ 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.

Note. the ``backward_transform`` attribute is equivalent to
``forward_transform.inverse``.

In this section, for illustration purpose, we will be using the same 2D imaging
WCS from ``imaging_wcs_wdist.asdf`` created in :ref:`imaging_example` whose
forward transformation converts image coordinates to world coordinates and
Expand Down Expand Up @@ -175,3 +178,4 @@ point far away from the image for which numerical inverse fails.
[5.31656943e-06 2.72052603e-10]
[6.81557583e-06 1.06560533e-06]
[3.96365344e-04 6.41822468e-05]]

Loading