改进Skyfield中日出/日落计算的性能

Question

给定坐标(经度/经度)和时间戳，我使用skyfield库为每个坐标/时间戳对计算太阳是否在该位置和时间升起：

import datetime as dt

import numpy as np
from skyfield import almanac
from skyfield.api import Loader
from skyfield.toposlib import wgs84
from skyfield_data import get_skyfield_data_path


load = Loader(get_skyfield_data_path())
eph = load("de421.bsp")
ts = load.timescale()

# Made up times
datetimes = [
    dt.datetime(2000, 1, 1, tzinfo=dt.timezone.utc)
    + dt.timedelta(minutes=m) for m in range(24 * 60)
]
times = ts.from_datetimes(datetimes)

# Made up coordinates
coordinates = np.array([
    (50.0, lon)
    for lon in np.linspace(0.0, 30.0, 1000)
])

results = []
for lat, lon in coordinates:
    pos = wgs84.latlon(lat, lon)
    is_day = almanac.sunrise_sunset(eph, pos)
    results.append(is_day(times))

这段代码可以工作，但速度太慢，不能满足我的需求。有没有办法将Python循环通过坐标移动到NumPy中？Skyfield文档说Skyfield使用NumPy，但据我所知，wgs84.latlon不接受NumPy数组作为输入。

如果无法避免Python循环，那么有没有一种方法可以提高速度呢？

请注意，在我的实际用例中，时间戳和坐标都是用户提供的，并且事先是未知的，所以我不能假设它们有任何特定的结构。

Answer 1

多处理是CPU密集型活动的常见解决方案:看起来就是这样。

由于酸洗错误，不能在多个进程之间共享eph。这意味着每个进程都必须调用load()。然而，这似乎是相对较快的。下面的代码在我的机器上运行大约3秒。这可能会也可能不会比你目前的经验更好，因为你的原始问题中没有提到执行的持续时间。

import time
import datetime as dt
import numpy as np
from skyfield import almanac
from skyfield.api import Loader
from skyfield.toposlib import wgs84
from skyfield_data import get_skyfield_data_path
from concurrent.futures import ProcessPoolExecutor


load = Loader(get_skyfield_data_path())

ts = load.timescale()

# Made up times
datetimes = [
    dt.datetime(2000, 1, 1, tzinfo=dt.timezone.utc)
    + dt.timedelta(minutes=m) for m in range(24 * 60)
]
times = ts.from_datetimes(datetimes)

# Made up coordinates
coordinates = np.array([
    (50.0, lon)
    for lon in np.linspace(0.0, 30.0, 1000)
])


def docalc(lat, lon):
    pos = wgs84.latlon(lat, lon)
    eph = load('de421.bsp')
    return almanac.sunrise_sunset(eph, pos)(times)


def main():
    results = []
    with ProcessPoolExecutor() as executor:
        futures = []
        for lat, lon in coordinates:
            futures.append(executor.submit(docalc, lat, lon))
        for future in futures:
            results.append(future.result())
    return results


if __name__ == '__main__':
    s = time.perf_counter()
    main()
    e = time.perf_counter()
    print(f'Duration = {e-s:.4f}s')