使用Skyfield计算地球到木星的距离

5
我尝试使用Skyfield绘制地球到太阳系行星的距离随时间变化的au图。这非常简单,甚至在包主页的前页中已经给出了。然而,虽然对于水星、金星和火星它可以完美工作,但对于其他行星则无法正常运行。我不熟悉JPL天文数据文件,但看起来例如木星在de421.bsp文件中没有关键项,这可能是问题所在。
以下是一个最小示例(来自主页):
from skyfield.api import load, now

planets = load('de421.bsp')
earth, planet = planets['earth'], planets['jupiter']

jd = now()
position = earth.at(jd).observe(planet)
ra, dec, distance = position.radec()

print(distance)

错误如下。请注意,如果在上面的代码中将“jupiter”替换为“mars”,它就不会崩溃。
---->  earth, planet = planets['earth'], planets['jupiter']
KeyError: "kernel 'de421.bsp' is missing 'JUPITER' - the targets it supports are:
SOLAR SYSTEM BARYCENTER, MERCURY BARYCENTER, VENUS BARYCENTER, EARTH BARYCENTER, 
MARS BARYCENTER, JUPITER BARYCENTER, SATURN BARYCENTER, URANUS BARYCENTER, 
NEPTUNE BARYCENTER, PLUTO BARYCENTER, SUN, MERCURY, VENUS, MOON, EARTH, MARS"

我是否以错误的天心使用了星历文件(错误的重心)?或者这只是 de421.bsp 文件的限制?我阅读了 Skyfield 网站上有关星历文件的描述(此处),但不确定我是否完全理解了它。
有什么建议可以用 Skyfield 进行地球 - 木星距离的简单计算吗?
谢谢!

1
你尝试过 planets['jupiter barycenter'] 吗? - Psytho
2个回答

8

根据错误提示,你需要使用JUPITER BARYCENTER代替jupiter


谢谢您的快速回答,它解决了问题。只是出于好奇,在文件中MARS BARYCENTER和MARS之间有区别吗?还是它只是一个别名?返回的距离值是相同的。 - Fabio
2
@Fabio 实际的JPL星历计算、制表和插值有点复杂,但是在阅读了一些链接之后(见下一个评论),我想地月重心(质心)被视为太阳系计算中的一个“行星”,然后将地球和月球之间的相对运动作为一个有点独立的计算来处理。对于地月来说,它们也被单独拆分出来,但对于木星质心(行星加上所有的卫星),Skyfield没有单独拆分。由于木星非常大,差异非常小。 - uhoh
请参考这个这个,以及这个 - uhoh
冥王星和卡戎的情况特别有趣 - 重心在两个物体之间的空间中。 - uhoh

3

如果有帮助的话,这只是补充 - 接受的答案解决了问题。

我想表明,由于位置是重心坐标,因此['solar system barycenter']将保持在原点。但我失败了,因为它返回一个零值而不是矢量(或None)。无论如何

sun motion in barycentric frame

import matplotlib.pyplot as plt
from skyfield.api import load, JulianDate

data = load('de421.bsp')
sun  = data['sun']
bary = data['solar system barycenter']

years = [1975+i for i in range(51)]
sunpos, barypos = [], []

for year in years:
    jd = JulianDate(utc=(year, 1, 1))
    sunpos.append(sun.at(jd).position.km)
    barypos.append(bary.at(jd).position.km)

plt.figure()
x, y, z = zip(*sunpos)
plt.plot(years, x)
plt.plot(years, y)
plt.plot(years, z)
# x, y, z = zip(*barypos)
# plt.plot(years, x, '-k')
# plt.plot(years, y, '-k')
# plt.plot(years, z, '-k')

plt.title('suns motion in barycentric frame')
plt.savefig('bary one')
plt.show()

下面的两个图表显示了地球和月球相对于地月重心的运动,称为Skyfield中的['earth barycenter']

various motions of earth, moon an barycenter

import matplotlib.pyplot as plt
import numpy as np
from skyfield.api import load, JulianDate

data  = load('de421.bsp')
earth = data['earth']
moon  = data['moon']
bary  = data['earth barycenter']

days = range(0, 366, 5)
earthpos, moonpos, barypos = [], [], []
for day in days:
    jd = JulianDate(utc=(2016, 1, day))  # seems to work
    earthpos.append(earth.at(jd).position.km)
    moonpos.append(moon.at(jd).position.km)
    barypos.append(bary.at(jd).position.km)
ep = np.array(earthpos).T
mp = np.array(moonpos).T
bp = np.array(barypos).T

plt.figure(figsize=[9,9])
plt.subplot(5,1,1)
for thing in ep:
    plt.plot(days, thing)
plt.subplot(5,1,2)
for thing in mp:
    plt.plot(days, thing)
plt.subplot(5,1,3)
for thing in bp:
    plt.plot(days, thing)
plt.subplot(5,1,4)
for thing in (ep-bp):
    plt.plot(days, thing)
plt.subplot(5,1,5)
for thing in (mp-bp):
    plt.plot(days, thing)
plt.savefig('bary two')
plt.show()

1
非常感谢这个详细的洞察!现在我更好地理解了JPL文件中的“重心”概念。 - Fabio

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接