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.

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

Support serialization of astropy.wcs.WCS objects to ASDF. #235

Closed
wants to merge 19 commits into from
Closed
Show file tree
Hide file tree
Changes from 12 commits
Commits
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
2 changes: 2 additions & 0 deletions CHANGES.rst
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,8 @@

- replace usages of ``copy_arrays`` with ``memmap`` [#230]

- Add support for astropy.wcs.WCS and astropy.wcs.wcsapi.SlicedLowLevelWCS [#235]

0.6.1 (2024-04-05)
------------------

Expand Down
3 changes: 2 additions & 1 deletion asdf_astropy/converters/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
"FitsConverter",
"AsdfFitsConverter",
"AstropyFitsConverter",
"FitsWCSConverter",
"ColumnConverter",
"AstropyTableConverter",
"AsdfTableConverter",
Expand Down Expand Up @@ -51,7 +52,7 @@
SkyCoordConverter,
SpectralCoordConverter,
)
from .fits import AsdfFitsConverter, AstropyFitsConverter, FitsConverter
from .fits import AsdfFitsConverter, AstropyFitsConverter, FitsConverter, FitsWCSConverter
from .table import AsdfTableConverter, AstropyTableConverter, ColumnConverter, NdarrayMixinConverter
from .time import TimeConverter, TimeDeltaConverter
from .transform import (
Expand Down
2 changes: 2 additions & 0 deletions asdf_astropy/converters/fits/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,8 @@
"FitsConverter",
"AsdfFitsConverter",
"AstropyFitsConverter",
"FitsWCSConverter",
]

from .fits import AsdfFitsConverter, AstropyFitsConverter, FitsConverter
from .fitswcs import FitsWCSConverter
17 changes: 17 additions & 0 deletions asdf_astropy/converters/fits/fitswcs.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,17 @@
from asdf.extension import Converter


class FitsWCSConverter(Converter):
tags = ("tag:astropy.org:astropy/fits/fitswcs-*",)
types = ("astropy.wcs.wcs.WCS",)

def from_yaml_tree(self, node, tag, ctx):
from astropy.wcs import WCS

primary_hdu = node["hdu"][0]
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
primary_hdu = node["hdu"][0]
primary_hdu = node["hdu"][0]

Doesn't this throw out the distortion (and other non-primary) data that is generated with to_fits? Would you add a test with a wcs that generates more than 1 hdu for to_fits?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

@braingram, Do you have any suggestions for effectively preserving distortion data?

Copy link
Contributor

@braingram braingram Aug 26, 2024

Choose a reason for hiding this comment

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

I believe that's already done with to_fits on write (since it puts the table into additional hdus that are stored in the asdf file with this PR).

This is all new to me so I could be getting this wrong but I think:

astropy.wcs.WCS(node["hdu"][0].header, fobj=node["hdu"])

would work for reading. Do any of the example files have distortion tables?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I am currently experimenting and trying to determine the appropriate serialization logic for this WCS. If you have any suggestions, please let me know.

return WCS(primary_hdu.header)

def to_yaml_tree(self, wcs, tag, ctx):
node = {}
node["hdu"] = wcs.to_fits()
return node
31 changes: 31 additions & 0 deletions asdf_astropy/converters/fits/tests/test_fitswcs.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,31 @@
import warnings

import asdf
import pytest
from astropy.io import fits
from astropy.utils.data import download_file
from astropy.wcs import WCS, FITSFixedWarning


def create_wcs():
Copy link
Member

Choose a reason for hiding this comment

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

I wonder if there's a set of WCSes we can use to test in astropy somewhere? There has to already be one or many fixtures in astropy we can import riiighghtttt?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I have extracted the WCS from the FITS file. I thought it would be better to test directly using a FITS file.

urls = [
"http://data.astropy.org/tutorials/FITS-cubes/reduced_TAN_C14.fits",
Copy link
Contributor

Choose a reason for hiding this comment

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

This would introduce downloading data files as part of the test runs. This is not something done by other tests and I don't think this belongs in unit tests.

Would you add a fixture that creates suitable wcs objects?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Please see this

"http://data.astropy.org/tutorials/FITS-images/HorseHead.fits",
]

with warnings.catch_warnings():
warnings.simplefilter("ignore", FITSFixedWarning)
return [WCS(fits.open(download_file(url, cache=True))[0].header) for url in urls]


@pytest.mark.filterwarnings("ignore::astropy.wcs.wcs.FITSFixedWarning")
@pytest.mark.parametrize("wcs", create_wcs())
def test_wcs_serialization(wcs, tmp_path):
file_path = tmp_path / "test_wcs.asdf"
with asdf.AsdfFile() as af:
af["wcs"] = wcs
af.write_to(file_path)

with asdf.open(file_path) as af:
loaded_wcs = af["wcs"]
assert wcs.to_header() == loaded_wcs.to_header()
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
assert wcs.to_header() == loaded_wcs.to_header()
assert wcs.to_header() == loaded_wcs.to_header()

This isn't testing the distortion information for the wcses.
I tried changing this to wcs == loaded_wcs and this is failing in part due to the sip attribute not round-tripping.

Copy link
Contributor

Choose a reason for hiding this comment

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

The wcs == loaded_wcs appears to be defaulting to object.__eq__ (so will fail since id(wcs) != id(loaded_wcs)). This helper function could be used instead:

def assert_hdu_list_equal(a, b):

and reports a failure for this test since the cards don't match:

-> assert tuple(card_a) == tuple(card_b)
(Pdb) p card_a
('A_0_0', -2.0413459666638e-06, 'SIP distortion coefficient')
(Pdb) p card_b
('A_0_0', -2.0413459666638133e-06, 'SIP distortion coefficient')

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Since there isn't a helper function available to check if two WCS objects are equivalent, I thought of checking the round-tripping of the WCS through its header.

5 changes: 5 additions & 0 deletions asdf_astropy/converters/slicedwcs/__init__.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
__all__ = [
"SlicedWCSConverter",
]

from .slicedwcs import SlicedWCSConverter
34 changes: 34 additions & 0 deletions asdf_astropy/converters/slicedwcs/slicedwcs.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,34 @@
from asdf.extension import Converter
Copy link
Contributor

Choose a reason for hiding this comment

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

I believe the general layout for this package is to put converters in submodules similar to the layout in astropy. Would you move this and fitswcs.py (and the corresponding tests) into a wcs submodule? Having them both contained in the new wcs submodule is sufficient (so the SlicedWCSConverter doesn't have to be further nested in a wcsapi.wrappers.sliced_wcs submodule).

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Agreed, moving them to a wcs submodule makes sense and will keep things organized.



class SlicedWCSConverter(Converter):
tags = ("tag:astropy.org:astropy/slicedwcs/slicedwcs-*",)
types = ("astropy.wcs.wcsapi.wrappers.sliced_wcs.SlicedLowLevelWCS",)

def from_yaml_tree(self, node, tag, ctx):
from astropy.wcs.wcsapi.wrappers.sliced_wcs import SlicedLowLevelWCS

wcs = node["wcs"]
slice_array = []
slice_array = [
s if isinstance(s, int) else slice(s["start"], s["stop"], s["step"]) for s in node["slices_array"]
]
return SlicedLowLevelWCS(wcs, slice_array)

def to_yaml_tree(self, sl, tag, ctx):
node = {}
node["wcs"] = sl._wcs
node["slices_array"] = []

for s in sl._slices_array:
if isinstance(s, slice):
braingram marked this conversation as resolved.
Show resolved Hide resolved
node["slices_array"].append(
{
"start": s.start if s.start is not None else None,
"stop": s.stop if s.stop is not None else None,
"step": s.step if s.step is not None else None,
ViciousEagle03 marked this conversation as resolved.
Show resolved Hide resolved
},
)
else:
node["slices_array"].append(s)
return node
Empty file.
35 changes: 35 additions & 0 deletions asdf_astropy/converters/slicedwcs/test/test_slicedwcs.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,35 @@
import asdf
import numpy as np
import pytest
from astropy.wcs import WCS
from astropy.wcs.wcsapi.wrappers.sliced_wcs import SlicedLowLevelWCS


def create_wcs():
wcs = WCS(naxis=4)
wcs.wcs.ctype = "RA---CAR", "DEC--CAR", "FREQ", "TIME"
wcs.wcs.cunit = "deg", "deg", "Hz", "s"
wcs.wcs.cdelt = -2.0, 2.0, 3.0e9, 1
wcs.wcs.crval = 4.0, 0.0, 4.0e9, 3
wcs.wcs.crpix = 6.0, 7.0, 11.0, 11.0
wcs.wcs.cname = "Right Ascension", "Declination", "Frequency", "Time"
wcs0 = SlicedLowLevelWCS(wcs, 1)
wcs1 = SlicedLowLevelWCS(wcs, [slice(None), slice(None), slice(None), 10])
wcs3 = SlicedLowLevelWCS(SlicedLowLevelWCS(wcs, slice(None)), [slice(3), slice(None), slice(None), 10])
wcs_ellipsis = SlicedLowLevelWCS(wcs, [Ellipsis, slice(5, 10)])
wcs2 = SlicedLowLevelWCS(wcs, np.s_[:, 2, 3, :])
return [wcs0, wcs1, wcs2, wcs_ellipsis, wcs3]


@pytest.mark.filterwarnings("ignore::astropy.wcs.wcs.FITSFixedWarning")
Copy link
Contributor

Choose a reason for hiding this comment

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

What's triggering this warning?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Same as the above reason

@pytest.mark.parametrize("sl_wcs", create_wcs())
def test_sliced_wcs_serialization(sl_wcs, tmp_path):
file_path = tmp_path / "test_slicedwcs.asdf"
with asdf.AsdfFile() as af:
af["sl_wcs"] = sl_wcs
af.write_to(file_path)

with asdf.open(file_path) as af:
loaded_sl_wcs = af["sl_wcs"]
assert sl_wcs._wcs.to_header() == loaded_sl_wcs._wcs.to_header()
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
assert sl_wcs._wcs.to_header() == loaded_sl_wcs._wcs.to_header()
assert sl_wcs._wcs.to_header() == loaded_sl_wcs._wcs.to_header()

If I try to use sl_wcs == loaded_sl_wcs it fails. Any idea why the wcs is not round-tripping?

assert sl_wcs._slices_array == loaded_sl_wcs._slices_array
4 changes: 4 additions & 0 deletions asdf_astropy/extensions.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,8 @@
from .converters.coordinates.sky_coord import SkyCoordConverter
from .converters.coordinates.spectral_coord import SpectralCoordConverter
from .converters.fits.fits import AsdfFitsConverter, AstropyFitsConverter
from .converters.fits.fitswcs import FitsWCSConverter
from .converters.slicedwcs.slicedwcs import SlicedWCSConverter
from .converters.table.table import AsdfTableConverter, AstropyTableConverter, ColumnConverter, NdarrayMixinConverter
from .converters.time.time import TimeConverter
from .converters.time.time_delta import TimeDeltaConverter
Expand Down Expand Up @@ -482,6 +484,8 @@
AstropyTableConverter(),
AstropyFitsConverter(),
NdarrayMixinConverter(),
FitsWCSConverter(),
SlicedWCSConverter(),
]

_COORDINATES_MANIFEST_URIS = [
Expand Down
17 changes: 17 additions & 0 deletions asdf_astropy/resources/manifests/astropy-1.1.0.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -54,3 +54,20 @@ tags:
title: NdarrayMixin column.
description: |-
Represents an astropy.table.NdarrayMixin instance.
- tag_uri: tag:astropy.org:astropy/slicedwcs/slicedwcs-1.0.0
schema_uri: http://astropy.org/schemas/astropy/slicedwcs/slicedwcs-1.0.0
title: Represents an instance of SlicedLowLevelWCS
description: |-
ViciousEagle03 marked this conversation as resolved.
Show resolved Hide resolved
The SlicedLowLevelWCS class is a wrapper class for WCS that applies slices
to the WCS, allowing certain pixel and world dimensions to be retained or
dropped.

It manages the slicing and coordinate transformations while preserving
the underlying WCS object.
- tag_uri: tag:astropy.org:astropy/fits/fitswcs-1.0.0
schema_uri: http://astropy.org/schemas/astropy/fits/fitswcs-1.0.0
title: FITS WCS (World Coordinate System) Converter
description: |-
Represents the FITS WCS object, the HDUlist of the FITS header is preserved
during serialization and during deserialization the WCS object is recreated
from the HDUlist.
18 changes: 18 additions & 0 deletions asdf_astropy/resources/schemas/fits/fitswcs-1.0.0.yaml
Original file line number Diff line number Diff line change
@@ -0,0 +1,18 @@
%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://astropy.org/schemas/astropy/fits/fitswcs-1.0.0"

title:
Represents the fits object
ViciousEagle03 marked this conversation as resolved.
Show resolved Hide resolved

description:
Represents the fits object
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
description:
Represents the fits object
description: >-
Represents the fits object

Same as title would you add a more descriptive description?


allOf:
- type: object
properties:
hdu:
tag: "tag:astropy.org:astropy/fits/fits-*"

required: ["hdu"]
24 changes: 24 additions & 0 deletions asdf_astropy/resources/schemas/slicedwcs/slicedwcs-1.0.0.yaml
Original file line number Diff line number Diff line change
@@ -0,0 +1,24 @@
%YAML 1.1
---
$schema: "http://stsci.edu/schemas/yaml-schema/draft-01"
id: "http://astropy.org/schemas/astropy/slicedwcs/slicedwcs-1.0.0"

title:
Represents the SlicedLowLevelWCS object

description:
Represents the SlicedLowLevelWCS object

allOf:
- type: object
properties:
wcs:
tag: "tag:astropy.org:astropy/fits/fitswcs-1*"
slices_array:
type: array
items:
- oneOf:
- type: integer
- type: object
braingram marked this conversation as resolved.
Show resolved Hide resolved

required: ["wcs", "slices_array"]
Loading