python - 随时间传播卫星目录的有效方法

Question

问题陈述

我需要传播 entire catalog of recent TLEs (需要免费帐户才能查看)来自 space-track.org 使用 skyfield或类似的。列表中通常有 15k-16k 个 TLE。我让它工作，但它非常慢。在服务器上使用 46 个内核的小时数。

未与 skyfield 结婚。如果 astropy 或 pyephem 或其他更快，我会很乐意接受显示我正在尝试使用它做什么的答案。

最小示例

对于我的应用程序，我将 TLE 加载到 Pandas 数据框中并在那里进行分析，因此我将在 Pandas 世界中保留我的示例。下面是最小示例。

假设卫星目录保存为catalog.txt，设置环境，然后读取TLEs，生成sf.sgp4lib.EarthSatellite对象，并将所有内容加载到 Pandas 数据框中。我们也offset positions to some observation site .我将选择观察点留给读者(0、0、0 都可以):

import skyfield as sf
import pandas as pd
from skyfield.api import load, Topos
from datetime import datetime, timezone, timedelta

with open('catalog.txt', 'r') as f:
    tle_list = [line.strip() for line in f.read().split('\n')
                if line is not '']
data = []
for i in range(0, len(tle), 2):  # every two lines
    temp = {}
    temp['tle1'] = tle_list[i]
    temp['tle2'] = tle_list[i+1]
    temp['earthsat'] = sf.sgp4lib.EarthSatellite(tle_list[i],
                                                 tle_list[i+1])
    data.append(temp)
df = pd.DataFrame(data=data)
site = Topos(latitude_degrees=site_lat,
             longitude_degrees=site_lon,
             elevation_m=site_alt)
df['earthsat'] = df.earthsat - site  # offset to site location

2.1 s ± 20 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

创建一组可识别时区的日期时间对象，以在其上传播所有卫星。在这里我选择了从我写这篇文章的那天午夜前 4 小时到之后 4 小时每 10 分钟。

ts = load.timescale()
tz = timezone(timedelta(hours=-4))  # Eastern or whatever timezone
midnight = datetime(2018, 4, 4, 0, 0, 0, tzinfo=tz)  # midnight today
start = midnight - timedelta(hours=4)
end = midnight + timedelta(hours=4)
delta_time = timedelta(minutes=10)
# This is ugly, but I had issues using linspace or arange...
times = [start]
now = start
while now <= end:
    now += delta_time
    times.append(now)

189 ms ± 36.9 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

最后，计算每颗卫星每个时间步长的天文位置。这是永远需要的。对我来说再次运行只是为了获得时间太长了，但在服务器上使用 46 个内核时大约需要几个小时。

df['astrometric'] = df.earthsat.apply(lambda x: [x.at(ts.utc(time)) for time in times])

其他详细信息

我找到了 this discussion on date arrays in the documentation这表明一次传递整个数组:x.at(ts.utc(times))。到目前为止，它需要的内核更少，运行速度更快，但仍然需要很长时间。

我通过为天体测量学创建生成器(最初是为什么我放弃了一次传递整个 times 数组的原因)暂时解决了这个问题，但最终我实际上需要评估事物，所以我永远无法避免繁重的工作。

如果最终用例适合某些特定的加速，我最终需要从这些对象中获取站点的观察角度，因此 [x.altaz() for x in row.astrometrics] 类型的东西。

解决方案的想法

我现在的想法是，我正在为目录中的每颗卫星计算整个晚上每个时间步长的观测点位置。我可能是错的，但如果我是对的，那么我想我会通过计算一次然后为每颗新卫星查找它来看到相当不错的加速(可能还不够)。有人知道怎么做吗？

另外，如果有更快的轨道传播器实现，或者有加速天空场实现的方法，我会很乐意接受一个答案，说明如何使用它来做我想做的事情(因此包括 astropy 和 pyephem 标签)。

谢谢。

Answer 1

我最好的建议是使用 NASA SPICE Toolkit对于这类工作。这让您加载两行元素，然后使用 NASA 的 NAIF/SPICE 内核完成其余工作(您可能还会发现一些其他格式的卫星，但 TLE 没问题)。

如果您使用 C 语言工作，您会使用 getelm_c阅读它们的方法和spkpos_c得到的位置。幸运的是，有一个名为 spiceypy 的 python 包装器。 !

getelm_c 方法被包裹在 spiceypy.spiceypy.getelm(frstyr, lineln, lines) 到 load在你的 TLE 中。您希望他们使用 spiceypy.spiceypy.spkpos(targ, et, ref, abcorr, obs) 来获取您的 position相对于引用体。

如果您决定使用 NASA 的 SPICE 内核，我建议使用来自 spiceypy 的文档，特别是卡西尼示例的位置，以确保您拥有所有正确的文件: https://spiceypy.readthedocs.io/en/master/exampleone.html