两个Topos()位置之间的路径:确定到达给定高度的纬度和经度

Question

我对轨道力学领域很陌生，目前正在为下面这个问题而挣扎，这个问题应该很容易用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)，但不知道如何继续下去。

在纬度和经度上完全用三角学方法进行试验，也没有产生预期的结果。

不幸的是，我并没有任何有效的结果可以用来更容易地解决这个问题。从三角学的角度来想象，很明显，卫星位置的高度升高会使交叉口的拉特隆离地球的位置更近。降低海拔会使这个交叉口离卫星更近。

Answer 1

这是一个有趣的问题，Skyfield的API提供了一个很难问的问题；如果你能勾勒出一个更大的问题，可以通过知道视线与特定高度的交集来解决，那么就有可能为将来处理同样问题的用户编写一个解决这个问题的例程。

同时：

为了让您的脚本运行，我不得不从api.

• The name dis导入Distance，因此我将其替换为distance_to_height，希望它是intended.

• Calling ts.now()在每次调用时给您的一个稍微不同的日期和时间。虽然脚本运行得如此之快，可能并不重要，但为了清楚起见，我只在脚本开始时只调用了一次now()，这也比反复调用它稍微快一些。(实际上，在这种情况下，计算速度要快得多，因为自转矩阵只能计算一次，而不必对每个单独的时间物体重新计算一次，但这是一个不易看到的隐藏细节。)

• ，我怀疑你的几何存在一个问题:我认为，只有在地球是平的情况下，100 / sin()策略才能起作用？但是也许你总是在处理几乎架空的卫星，所以这个错误是可控的？(或者我可能是在错误地想象几何；如果数学实际上是correct.)

• For可读性的话，可以提供一个图表，我会给出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

对于未来，我在Skyfield TODO.rst文件中添加了一些想法，这些想法可能在将来一起走向解锁一种更加惯用的方法来执行这种计算--尽管我怀疑还有更多的步骤是必要的：

https://github.com/skyfielders/python-skyfield/commit/ba1172a0ccfef84473436d9d7b8a7d7011344cbd