spherical_offsets_by followed by spherical_offsets_to fails

[thusser]

Description

I'm trying to convert Alt/Az offsets to RA/Dec offsets at a given coordinate on the sky. For this I transform the coordinates to an AltAz frame, then add an offset via spherical_offsets_by, and convert the result back to ICRS. Then I can use spherical_offsets_to to calculate the offset between the two coordinates. See code for details.

Expected behavior

It should work.

Actual behavior

It does not work and I get a:

ValueError: Tried to use spherical_offsets_to with two non-matching frames!

Steps to Reproduce

Assume a coordinate coords that includes location and obstime, then this should work:

altaz = coords.transform_to(frame="altaz")
altaz2 = altaz.spherical_offsets_by(daz, dalt)
coords2 = altaz2.transform_to(frame="icrs")
dra, ddec = coords.spherical_offsets_to(coords2)

But it only works, when stripping coords and coords2 of all other attributes before calling spherical_offsets_to (thanks to maxnoe for this workaround):

coords = SkyCoord(coords.ra, coords.dec)
coords2 = SkyCoord(coords2.ra, coords2.dec)

System Details

Linux-5.13.0-35-generic-x86_64-with-glibc2.34
Python 3.9.7 (default, Sep 10 2021, 14:59:43)
[GCC 11.2.0]
Numpy 1.22.3
pyerfa 2.0.0.1
astropy 5.0.2
Scipy 1.8.0
Matplotlib 3.5.0

[github-actions]

Welcome to Astropy 👋 and thank you for your first issue!

A project member will respond to you as soon as possible; in the meantime, please double-check the guidelines for submitting issues and make sure you've provided the requested details.

GitHub issues in the Astropy repository are used to track bug reports and feature requests; If your issue poses a question about how to use Astropy, please instead raise your question in the Astropy Discourse user forum and close this issue.

If you feel that this issue has not been responded to in a timely manner, please leave a comment mentioning our software support engineer @embray, or send a message directly to the development mailing list. If the issue is urgent or sensitive in nature (e.g., a security vulnerability) please send an e-mail directly to the private e-mail feedback@astropy.org.

[maxnoe]

Here is a self-contained working example:

from astropy.coordinates import SkyCoord, AltAz, EarthLocation
from astropy.time import Time
import astropy.units as u

location = EarthLocation.of_site("Roque de los Muchachos")
obstime = Time.now()

daz = 5 * u.deg
dalt = 1 * u.deg

icrs_coord = SkyCoord.from_name("Crab Nebula")
altaz_frame = AltAz(location=location, obstime=obstime)

altaz_coord = icrs_coord.transform_to(altaz_frame)
altaz_coord_offset = altaz_coord.spherical_offsets_by(d_lon=daz, d_lat=dalt)
icrs_coord_offset = altaz_coord_offset.transform_to("icrs")

# remove additional frame attrs to avoid error this issue is about
icrs_coord_offset = SkyCoord(icrs_coord_offset.ra, icrs_coord_offset.dec)

dra, ddec = icrs_coord.spherical_offsets_to(icrs_coord_offset)

print(dra, ddec)

[ayshih]

The "problem" is the frame attributes that are carried along by SkyCoord even when they don't correspond to the current frame. In the original example, coords contains the frame attributes location and obstime, which have no meaning for ICRS, and coords2 additionally has obswl, pressure, relative_humidity, and temperature because of the transformation path through AltAz. It makes sense for the equivalency check in spherical_offsets_to() to check all frame attributes, even seemingly irrelevant ones, because it could get weird if two SkyCoords are deemed equivalent in one frame, but no longer equivalent if they are both transformed to a different frame that actually uses those attributes.

The simplest workaround is to do:

dra, ddec = coords.spherical_offsets_to(coords2.frame)

The .frame turns coords2 into an ICRS instance, thus dropping the additional frame attributes carried by SkyCoord. The internal code of spherical_offsets_to() then similarly converts coords to an ICRS instance as well before performing the equivalency check.

[ayshih]

Incidentally, for both posts, the simplest way to make a "clean" SkyCoord is usually:

coords = SkyCoord(coords.frame)

which works for any frame, not just ICRS.

[thusser]

dra, ddec = coords.spherical_offsets_to(coords2.frame)

This one actually works, thanks.
I still don't find it very obvious to do this. Maybe it would be good to add a sentence about this to the documentation?

[ayshih]

You are certainly justified at being confused and not knowing how to get things working. At the very least, a note could be added to the docstring for spherical_offsets_to(), but one could also jazz up the equivalency check therein so that the ValueError message is more helpful in the situation when the equivalency check would have succeeded if not for irrelevant frame attributes.