在天空中确定月食
Determining lunar eclipse in skyfield
我得到了 UTC 日期列表,所有时间都转换为 00:00。
我想确定在给定的一天(即过去 24 小时)是否发生了(月食)
考虑 python 片段
from sykfield.api import load
eph = load('de421.bsp')
def eclipticangle(t):
moon, earth = eph['moon'], eph['earth']
e = earth.at(t)
x, y, _ = e.observe(moon).apparent().ecliptic_latlon()
return x.degrees
我假设一个人能够通过
确定日食是否在某个时间 t 的 24 小时内发生
- 检查第一个角度是否足够接近 180(简单)
- 检查二度是否足够接近 0(不是很简单?)
现在,就评论中的答案而言,仅通过测试角度是否接近 0 来解决第二个问题并不是那么简单。
因此,我的问题是
Can someone provide a function to determine if a lunar eclipse occurred on a given day t?
编辑。编辑此问题以反映 Brandon Rhodes 在下方评论中留下的反馈。
我刚刚阅读了 天文年历的解释性补充 的第 11.2.3 节,并尝试将其转换为 Skyfield Python 代码。这是我想出的:
import numpy as np
from skyfield.api import load
from skyfield.constants import ERAD
from skyfield.functions import angle_between, length_of
from skyfield.searchlib import find_maxima
eph = load('de421.bsp')
earth = eph['earth']
moon = eph['moon']
sun = eph['sun']
def f(t):
e = earth.at(t).position.au
s = sun.at(t).position.au
m = moon.at(t).position.au
return angle_between(s - e, m - e)
f.step_days = 5.0
ts = load.timescale()
start_time = ts.utc(2019, 1, 1)
end_time = ts.utc(2020, 1, 1)
t, y = find_maxima(start_time, end_time, f)
e = earth.at(t).position.m
m = moon.at(t).position.m
s = sun.at(t).position.m
solar_radius_m = 696340e3
moon_radius_m = 1.7371e6
pi_m = np.arcsin(ERAD / length_of(m - e))
pi_s = np.arcsin(ERAD / length_of(s - e))
s_s = np.arcsin(solar_radius_m / length_of(s - e))
pi_1 = 0.998340 * pi_m
sigma = angle_between(s - e, e - m)
s_m = np.arcsin(moon_radius_m / length_of(e - m))
penumbral = sigma < 1.02 * (pi_1 + pi_s + s_s) + s_m
partial = sigma < 1.02 * (pi_1 + pi_s - s_s) + s_m
total = sigma < 1.02 * (pi_1 + pi_s - s_s) - s_m
mask = penumbral | partial | total
t = t[mask]
penumbral = penumbral[mask]
partial = partial[mask]
total = total[mask]
print(t.utc_strftime())
print(0 + penumbral + partial + total)
它生成月食发生时间的向量,然后对月食的总程度进行评级:
['2019-01-21 05:12:51 UTC', '2019-07-16 21:31:27 UTC']
[3 2]
它的日食时间与 NASA 巨大的 table 月球历书给出的时间相差 3 秒以内:
我得到了 UTC 日期列表,所有时间都转换为 00:00。
我想确定在给定的一天(即过去 24 小时)是否发生了(月食)
考虑 python 片段
from sykfield.api import load
eph = load('de421.bsp')
def eclipticangle(t):
moon, earth = eph['moon'], eph['earth']
e = earth.at(t)
x, y, _ = e.observe(moon).apparent().ecliptic_latlon()
return x.degrees
我假设一个人能够通过
确定日食是否在某个时间 t 的 24 小时内发生- 检查第一个角度是否足够接近 180(简单)
- 检查二度是否足够接近 0(不是很简单?)
现在,就评论中的答案而言,仅通过测试角度是否接近 0 来解决第二个问题并不是那么简单。
因此,我的问题是
Can someone provide a function to determine if a lunar eclipse occurred on a given day t?
编辑。编辑此问题以反映 Brandon Rhodes 在下方评论中留下的反馈。
我刚刚阅读了 天文年历的解释性补充 的第 11.2.3 节,并尝试将其转换为 Skyfield Python 代码。这是我想出的:
import numpy as np
from skyfield.api import load
from skyfield.constants import ERAD
from skyfield.functions import angle_between, length_of
from skyfield.searchlib import find_maxima
eph = load('de421.bsp')
earth = eph['earth']
moon = eph['moon']
sun = eph['sun']
def f(t):
e = earth.at(t).position.au
s = sun.at(t).position.au
m = moon.at(t).position.au
return angle_between(s - e, m - e)
f.step_days = 5.0
ts = load.timescale()
start_time = ts.utc(2019, 1, 1)
end_time = ts.utc(2020, 1, 1)
t, y = find_maxima(start_time, end_time, f)
e = earth.at(t).position.m
m = moon.at(t).position.m
s = sun.at(t).position.m
solar_radius_m = 696340e3
moon_radius_m = 1.7371e6
pi_m = np.arcsin(ERAD / length_of(m - e))
pi_s = np.arcsin(ERAD / length_of(s - e))
s_s = np.arcsin(solar_radius_m / length_of(s - e))
pi_1 = 0.998340 * pi_m
sigma = angle_between(s - e, e - m)
s_m = np.arcsin(moon_radius_m / length_of(e - m))
penumbral = sigma < 1.02 * (pi_1 + pi_s + s_s) + s_m
partial = sigma < 1.02 * (pi_1 + pi_s - s_s) + s_m
total = sigma < 1.02 * (pi_1 + pi_s - s_s) - s_m
mask = penumbral | partial | total
t = t[mask]
penumbral = penumbral[mask]
partial = partial[mask]
total = total[mask]
print(t.utc_strftime())
print(0 + penumbral + partial + total)
它生成月食发生时间的向量,然后对月食的总程度进行评级:
['2019-01-21 05:12:51 UTC', '2019-07-16 21:31:27 UTC']
[3 2]
它的日食时间与 NASA 巨大的 table 月球历书给出的时间相差 3 秒以内: