How to save an Astropy QTable with a galactocentric skycoord to a FITS file?

Question

Creating a SkyCoord object:

c_end = SkyCoord(
    x=xf,
    y=yf,
    z=zf,
    v_x=vxf,
    v_y=vyf,
    v_z=vzf,
    frame='galactocentric'
)

and attempting to put this into an Astropy QTable with other parameters:

tab = QTable([c_end, ...], names = ['skycoord', ...])

tab.write('name.fits', format = 'fits')

Yet I'm getting this error:

yaml.representer.RepresenterError: ('cannot represent an object', <ICRS Coordinate: (ra, dec) in deg (266.4051, -28.936175)>)

which is suspicious as these are the exact coordinates of the galactic centre in the icrs frame. Any thoughts how to overcome this?

I've previously done something similar but with a SkyCoord in the ICRS frame, which worked without issue.

I tried converting to ICRS and saving that but it gives the same error.

To obtain the error you can run this minimum reproducible example.

from astropy.coordinates import SkyCoord
from astropy.table import QTable
import astropy.units as u

xf, yf, zf = 8.3, 0.0, 0.02
vxf, vyf, vzf = -11.1, 232.24, 7.25

c_end = SkyCoord(
    x=xf*u.kpc,
    y=yf*u.kpc,
    z=zf*u.kpc,
    v_x=vxf*u.km/u.s,
    v_y=vyf*u.km/u.s,
    v_z=vzf*u.km/u.s,
    frame='galactocentric'
)

tab = QTable([[c_end]], names=['skycoord'])

try:
    tab.write('name_fail.fits', format='fits', overwrite=True)
    print("Wrote FITS successfully.")
except Exception as e:
    print("Failed to write FITS:", e)

Results in the following output when run with astropy 7.1

Failed to write FITS: ('cannot represent an object', <ICRS Coordinate: (ra, dec) in deg (266.4051, -28.936175)>)

Answer 1

Unfortunately there is a real bug in Astropy, so thank you for posting this to SO. It's worth noting that Astropy can write SkyCoord objects to FITS, so this is not a fundamental limitation of the format.

In this case it turns out that for the galactocentric frame, the SkyCoord object has a galcen_coord attribute which is itself an ICRS Frame object. In order to serialize the column containing your c_end object, astropy needs to represent that galcen_coord object as YAML. That does not currently work in astropy but it will be fixed.

Here is a workaround you can use now:

import astropy.coordinates as coords
import astropy.units as u
from astropy.coordinates import SkyCoord
from astropy.io.misc import yaml
from astropy.table import QTable


def add_frames_to_yaml():
    """Add astropy coordinate frames to classes handled by YAML"""
    for name in coords.frame_transform_graph.get_names():
        cls = coords.frame_transform_graph.lookup_name(name)
        tag = "!astropy.coordinates." + name
        yaml.AstropyDumper.add_multi_representer(cls, yaml._quantity_representer(tag))
        yaml.AstropyLoader.add_constructor(tag, yaml._quantity_constructor(cls))


add_frames_to_yaml()

xf, yf, zf = 8.3, 0.0, 0.02
vxf, vyf, vzf = -11.1, 232.24, 7.25

sc = SkyCoord(
    x=xf * u.kpc,
    y=yf * u.kpc,
    z=zf * u.kpc,
    v_x=vxf * u.km / u.s,
    v_y=vyf * u.km / u.s,
    v_z=vzf * u.km / u.s,
    frame="galactocentric",
)

qt = QTable([[sc]], names=["sc"])
qt.write("test.fits", overwrite=True)
t2 = QTable.read("test.fits")
print(t2)
print(t2["sc"].frame)