在Python中获取地球上某个位置的J2000 XYZ坐标

Question

编辑：刚刚发现我需要把地球上某个位置的纬度、经度和海拔转换成J2000坐标，与ra/dec或月球无关。很抱歉。你的回答确实给了我很多见解。请看下面编辑过的问题。

问题：如何将纬度、经度和高程转换为J2000坐标（XYZ）。以弗所有皈依吗？我查看了文档，但找不到我需要的东西（或者可能因为我在这个领域缺乏知识而忽略了一些东西）。谢谢

**********************旧（不考虑）*******************

我有月球在赤经（RA）和赤纬（Dec）的位置，我想把它们转换成X-Y-Z坐标。有内置的PyEphem函数吗？还有，它背后的数学原理是什么？谢谢。

编辑：我使用的是J2000坐标系（我认为这是赤道坐标系，这是我第一次使用天文学）。我有到月球的距离。ra/dec值已经在J2000（赤道）坐标系中。

X点北

Y点西

Z指向天空

Answer 1

最佳答案：

我刚刚注意到，2011年6月，海军天文台发布了一个强大的NOVAS参考软件的Python接口，使用该接口可以执行最高精度的天文计算：

http://aa.usno.navy.mil/software/novas/novas_py/novaspy_intro.php

有了这个库，你可以得到你正在寻找的答案，其精度远远高于毕节所提供的：

from novas import compat as novas
jd_tt = novas.julian_date(2012, 9, 8, 12.00)
delta_t = 66.603  # from http://maia.usno.navy.mil/ser7/deltat.preds
lat = 42.3583     # positive is north
lon = -71.0603    # negative is west
observer = novas.make_observer_on_surface(lat, lon, 0, 0, 0)
print novas.geo_posvel(jd_tt, delta_t, observer)

在我的机器上，这给出了答案：

（-3.5081406460494928e-06、3.135277713623258e-05、2.858239912567112e-05）， （-0.00019753847060164693，-2.232994271278055E-05，2.4885824275734915e-07）

你自己试试看这是否能给你所需要的结果！

---

较新的答案：

似乎答案是“不”—令我惊讶的是，皮埃芬并没有给出一个简单的方法来回答这个问题：“在x，y，z坐标系中，波士顿在时间t？“。”

这是一个惊喜，因为“libastro”，即PyEphem后面的库，当然必须在内部计算这个，以便找出其他对象相对于观察者的位置。似乎有两个地方是这样。在parallax.c中，它定义了ta_par()，它只讨论外部的角度，但在内部，您可以看到它临时计算观察者的x，y，z。你甚至可以看到隐藏在里面的重要常数298.257，它测量地球有多平坦，因为地球不是一个完美的球体。

另一个地方在earthsat.c中，它看起来与“libastro”的其余部分完全不同，因此它复制了一些逻辑。它的EarthFlat常数298.25有点不太精确，但做的是相同的工作。它的函数GetSitPosition()实际上公开了x，y，z坐标，而不是隐藏它们。但是它被声明为static，所以没有办法从外部调用这个有用的函数！

所以目前，PyEphem无法直接计算x，y，z。但它确实提供了一个重要的信息：当前的恒星时间，你（我认为）可以用它来计算波士顿（或任何地方）在地球周围经过的时间有多远，这对计算你的坐标很重要。

我将看看能否用Python编写一个快速的解决方案，将PyEphem中的小时角与一些显式三角法结合起来，得到一个答案。但是，就目前而言，不：不幸的是，PyEphem没有直接公开这些信息；我将把它放在未来版本的清单上！

---

较老的答案，从问题是关于月球的位置时起：

遗憾的是，PyEphem没有内置的功能，可以将业余天文学中使用的极坐标转换为x/y/z坐标，这样就可以绘制出地球周围空间中物体的分布情况。但是你自己很容易转换：

import ephem
import math

m = ephem.Mars('2012/8/1')
print m.ra, m.dec

x = math.cos(m.dec) * math.cos(m.ra)
y = math.cos(m.dec) * math.sin(m.ra)
z = math.sin(m.dec)

print x, y, z
print 'sanity check: vector length =', math.sqrt(x*x + y*y + z*z)

此脚本的输出为：

12:58:51.20 -6:24:05.6
-0.961178016954 -0.252399543786 -0.111495695074
sanity check: vector length = 1.0

我在这里使用的随机日期的火星位置是相当合理的值：一个在大圆的一半以上几乎一个小时的RA（因为12小时正好是一半），以及一个将位置推向南边一点的偏斜。因此，我们得到的x、y和z：由于-6°确实在赤道以南，x和y都是负数，因为在一个24小时的圆上绕13小时会使你处于一个正常单位圆的负/负象限。

请注意，尽管J2000有一个北和南-因此我们可以如实地说，稍负的z是一个向南的方向-它确实有一个东和西，因为它下面的地球在不断地向东和西摆动。相反，RA从“Ares的第一点”开始测量，即春分时太阳所在的方向。所以x和y不是东方或西方；它们是指向太阳系的坐标，在一个固定的轴上，由地球每年春天所在的方向来定义。

我创建的这个xyz向量是一个“单位向量”——一个大小为1.0的小向量，我在脚本中验证了这一点，以确保公式正确无误。如果你在计算距离地球的物体的x{}和z坐标，那么你就可以得到一个真正的向量，其大小是距离，而不是1的分数，将三个xy和z乘以到物体的距离。

这对你有帮助吗？从你的描述——还有你关于东西方的问题——我不知道你是想让RA和dec变成xyz，还是你真的想转换方位角和高度（但两者的数学公式都是一样的）。看起来像是：

x = math.cos(m.alt) * math.cos(m.az)
y = math.cos(m.alt) * math.sin(m.az)
z = math.sin(m.alt)

你想用这些坐标来完成什么？这将有助于我们确保以一种有用的格式将它们提供给您。