我对轨道力学领域很陌生，目前在以下问题上有点挣扎，这个问题应该很容易用Skyfield解决，但我对所有不同的坐标系和它们之间的平移有点不知所措。
我在地球上有一个Topos位置和一个LEO卫星的Topos位置。我正在考虑它们之间的视线。我想确定沿着这条路径的位置的纬度和经度，它与大气的特定层相交。
一个例子是中间层和基于纬度和经度给出的关于其约100公里处属性的现有数据集。交叉点将使我更好地了解这些属性与卫星通信的相互作用。
我尝试直接使用Skyfield执行此操作，但最终只能得到一个Apparent对象，我无法将其转换回地球上的纬度，经度。首先，我三角测量了从地球到达到 100 公里高度的点的距离。
然后，我取了地球上的位置，并使用不变的高程，方位角来保持路径的方向，最后加上计算出的距离到达这个位置。我想我需要得到一个Geocentric对象来使用subpoint()才能获得该位置所需的纬度、经度。
这是我到目前为止所拥有的：
from skyfield.api import load, Distance
from skyfield.toposlib import Topos
import numpy as np
ts = load.timescale()
earth_position = Topos('52.230039 N', '4.842402 E', elevation_m=10)
space_position = Topos('51.526200 N', '5.347795 E', elevation_m=625 * 1000)
difference = (space_position - earth_position).at(ts.now()).altaz()
distance_to_height = 100 / np.sin(difference[0].radians)
position = earth_position.at(ts.now()).from_altaz(alt_degrees=difference[0].degrees, az_degrees=difference[1].degrees, distance=Distance(km=distance_to_height))
我已经多次浏览了文档，并偶然发现了通用 ICRF 对象的frame_latlon(frame)，但不确定如何进一步进行。
完全尝试对纬度和经度进行三角测量也没有产生预期的结果。
不幸的是，我真的没有任何经过验证的结果可以用来更轻松地解决这个问题。再次从三角学上想象，很明显，卫星位置高度的增加会使交叉点的纬度更接近地球上的位置。降低高度将使这个交叉点更接近卫星。
这是一个有趣的问题，Skyfield 的 API 没有提供简单的方法来询问; 如果你能概述通过了解视线与特定高度的交集来解决的更大问题，那么就有可能为未来的用户编写一个解决这个问题的例程来解决同样的问题。
同时：
- 为了让您的脚本运行，我必须从
api导入Distance。 dis这个名字不被认出来，所以我把它换成了distance_to_height，希望这是这个名字。- 呼叫
ts.now()会为您提供每次呼叫的略有不同的日期和时间。虽然脚本运行得如此之快以至于可能无关紧要，但为了清楚起见，我不得不在脚本开始时只调用now()一次，这也比重复调用它略快。(实际上在这种情况下，它要快得多，因为旋转矩阵只计算一次，而不必为每个单独的时间对象重新计算，但这是一个隐藏的细节，不容易看到。 - 我怀疑你的几何形状有问题：我认为，
100 / sin()机动只有在地球是平的时才有效？但也许您总是在处理近头顶的卫星，因此错误是可控的？(或者也许我错误地想象了几何形状;如果数学实际上是正确的，请随时提供图表。 - 为了便于阅读，我将给出
altaz()名称的组成部分而不是数字。
有了这些调整， 我认为答案是您需要手动构建Geocentric通过将观察者的位置相加来定位 和您创建的相对向量 在观察者和沿视线的100公里点之间。 必须采取这样的手动步骤表明了一个可能的区域 天空场可以改进的地方。 下面是它在代码中的外观：
from skyfield.api import load, Distance
from skyfield.positionlib import Geocentric
from skyfield.toposlib import Topos
import numpy as np
ts = load.timescale()
t = ts.now()
earth_position = Topos('52.230039 N', '4.842402 E', elevation_m=10)
space_position = Topos('51.526200 N', '5.347795 E', elevation_m=625 * 1000)
alt, az, distance = (space_position - earth_position).at(t).altaz()
distance_to_height = 100 / np.sin(alt.radians)
e = earth_position.at(t)
p = e.from_altaz(alt_degrees=alt.degrees, az_degrees=az.degrees, distance=Distance(km=distance_to_height))
g = Geocentric(e.position.au + p.position.au, t=t)
s = g.subpoint()
print(s)
print(s.elevation.km, '<- warning: 100/sin() did not produce exactly 100')
我看到的结果是：
Topos 52deg 06' 30.0" N 04deg 55' 51.7" E
100.02752954478532 <- warning: 100/sin() did not produce exactly 100
而对于未来， 我在天空场TODO.rst文件中添加了一些想法 一起可能会解锁一种更惯用的方式 将来进行这种计算—— 尽管我怀疑还需要采取更多步骤：
https://github.com/skyfielders/python-skyfield/commit/ba1172a0ccfef84473436d9d7b8a7d7011344cbd