有效地将 Astropy SkyCoord 数组导出到 array/DataFrame/file

Efficiently export Astropy SkyCoord array to array/DataFrame/file

A​​stropy 坐标已被证明是一种在 ICRS 和银河坐标系之间转换的非常快速的方法。我的问题是关于将 SkyCoord 数组写入文本文件,可能使用 Numpy 数组或 Pandas DataFrame 作为中介。我的阵列相当大(~100K 坐标)所以效率是可取的。

SkyCoord class 中似乎没有任何方法可以将 RA/Dec 导出为任何 Python 内置数据结构。我探索过的一种可能的方法是为坐标创建一个 Astropy QTable(将它们添加为字符串),但这很慢,并且需要额外的处理来分离和清理坐标字符串。

例如

import numpy as np
import pandas as pd
from astropy import units as u
from astropy.coordinates import SkyCoord
from astropy.io import ascii
from astropy.table import QTable

if __name__ == "__main__":
    RA = np.array([0.819853, 0.0920091, 1.04977])
    DEC = np.array([-30.8484, -30.671, -31.0921])

    c = SkyCoord(ra=RA * u.degree, dec=DEC * u.degree, frame="icrs")
    c_galactic = c.galactic

    table = QTable([c_galactic], names=["skycoord"])
    ascii.write(
        table, "gal_coords.csv", format="csv", fast_writer=False, overwrite=True
    )

我们非常欢迎任何有关实现此目标的最佳方法的建议。如果解决方案可以推广到将 SkyCoord 数组转换为其他数据结构,我希望它对更广泛的社区有用。

我可能对这个 SO 问题感到有点高兴;通过更多的挖掘,我找到了我想要的东西。希望解决方案对其他人有帮助:

import numpy as np
import pandas as pd
from astropy import units as u
from astropy.coordinates import SkyCoord

if __name__ == "__main__":
    RA = np.array([0.819853, 0.0920091, 1.04977])
    DEC = np.array([-30.8484, -30.671, -31.0921])

    c = SkyCoord(ra=RA * u.degree, dec=DEC * u.degree, frame="icrs")

    c_galactic = c.transform_to("galactic")

    coord_df = pd.DataFrame({"l": c_galactic.l.degree, "b": c_galactic.b.degree})
    coord_df.to_csv("gal_coords.csv", index=False)

Astropy 具有使用 ECSV format 处理将 Table 中的坐标列导出到 CSV 的内置功能,请参见下面的示例。这会自动生成与具有相关坐标属性的原始列相对应的列名称,在本例中为 Galactic lb

ECSV 是 CSV 的扩展,它在文件开头的注释文本中包含元数据,以完全指定列数据并允许数据的无损往返。但重要的是,您可以通过指定 # 是评论(例如 pandas.read_csv(..., comment='#').

来使用任何 CSV reader(不仅仅是 astropy)阅读此内容
>>> RA = np.array([0.819853, 0.0920091, 1.04977])
... DEC = np.array([-30.8484, -30.671, -31.0921])
... 
... c = SkyCoord(ra=RA * u.degree, dec=DEC * u.degree, frame="icrs")
... c_galactic = c.transform_to("galactic")
...
>>> t = QTable([c_galactic], names=["skycoord"])
>>> t.write('skycoord.ecsv')
>>> cat skycoord.ecsv
# %ECSV 0.9
# ---
# datatype:
# - {name: skycoord.l, unit: deg, datatype: float64}
# - {name: skycoord.b, unit: deg, datatype: float64}
# meta: !!omap
# - __serialized_columns__:
#     skycoord:
#       __class__: astropy.coordinates.sky_coordinate.SkyCoord
# ...
# schema: astropy-2.0
skycoord.l skycoord.b
10.618564158206215 -78.83862258651922
12.320325970465513 -78.28290679292496
9.106328561289311 -78.95458548581124

>>> pd.read_csv('skycoord.ecsv', comment='#')
                   skycoord.l skycoord.b
0  10.618564158206215 -78.83862258651922
1  12.320325970465513 -78.28290679292496
2   9.106328561289311 -78.95458548581124

另一种不使用 ECSV 的方法很简单:

>>> RA = np.array([0.819853, 0.0920091, 1.04977])
... DEC = np.array([-30.8484, -30.671, -31.0921])
... 
... c = SkyCoord(ra=RA * u.degree, dec=DEC * u.degree, frame="icrs")
... c_galactic = c.transform_to("galactic")
...
>>> t = QTable([c_galactic], names=["skycoord"])
>>> t.to_pandas().to_csv('skycoord.csv', index=False)
>>> cat skycoord.csv
skycoord.l,skycoord.b
10.618564158206215,-78.83862258651922
12.320325970465513,-78.28290679292496
9.106328561289311,-78.95458548581124

这种方法会丢弃任何元数据,但很容易工作。