Система координат местного наблюдателя на спутниковой панели (система координат LVLH)

В небесном поле я хочу рассчитать расстояние по высоте от одного эталонного спутника до следующего (для имитации наведения антенны). Мой текущий обходной путь — создать эталонного наблюдателя на высоте sat1 и вычислить углы alt-az:

sat_observer = sat_ref.subpoint()
ref_PoV = Topos(sat_observer.latitude, sat_observer.longitude, elevation_m=sat_observer.elevation.m)

satellite = Sat(name)
orbit = (satellite - ref_PoV).at(time[0])
el, az, distance = orbit.altaz()

Код может вычислять относительные углы наведения за один раз. Этот метод не работает для абсолютных значений, так как привязка alt-az всегда указывает на север Земли (система координат местного наблюдателя на спутниковой панели будет вращаться вокруг оси alt в следующий момент времени).

Каким будет наиболее эффективный способ определить локальную систему координат (кадр LVLH) на спутнике в поле неба, чтобы получить углы альт-азимут, соответствующие этому кадру LVLH?

Это интересный вопрос! Означает ли "Это работает...", что это не сбой или что он возвращает именно те значения, которые вы ожидаете? Если хотите, вы также можете перейти на страницу Skyfield на github и опубликовать там проблему, а также включить ссылку на свой вопрос здесь. Есть некоторое совпадение с пользователями здесь, в Astronomy SE и в github Skyfield. (см. этот вопрос для ссылок)
Кроме того, существует более одной возможной системы координат локального наблюдателя. Можете ли вы объяснить более подробно, что вы ищете? См., например. Являются ли системы координат LVLH и RSW одним и тем же?
@uhoh Спасибо за ссылку на разные системы координат, я обновил свой вопрос и свяжу его с github Skyfield.
выглядит хорошо, спасибо!

Ответы (1)

Итак, после помощи @uhoh, изучения этого поста и обсуждения здесь , мне удалось создать этот минимальный рабочий пример. Комментарии приветствуются.

from skyfield.api import Loader, EarthSatellite
from skyfield.api import Topos, load
from skyfield.timelib import Time
import skyfield.functions as sf
from sklearn import preprocessing
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
%matplotlib inline

halfpi, pi, twopi = [f*np.pi for f in [0.5, 1, 2]]
degs, rads = 180/pi, pi/180

ts = load.timescale()
line1 = '1 25544U 98067A   14020.93268519  .00009878  00000-0  18200-3 0  5082'
line2 = '2 25544  51.6498 109.4756 0003572  55.9686 274.8005 15.49815350868473'
satellite = EarthSatellite(line1, line2, 'ISS (ZARYA)', ts)
print(satellite)

line1 = '1 43205U 18017A   18038.05572532 +.00020608 -51169-6 +11058-3 0  9993'
line2 = '2 43205 029.0165 287.1006 3403068 180.4827 179.1544 08.75117793000017'
satellite2 = EarthSatellite(line1, line2, 'Roadster', ts)
print(satellite2)

time = ts.utc(2020, 24, 11, np.arange(0, 1, 0.01))

#calculate LVLH reference frame for the reference sat
#Z = - R / ||R||
#Y = Z X V / ||Z X V||
#X = Y X Z
R = satellite.at(time).position.km.T
V = satellite.at(time).velocity.km_per_s.T
Z = -preprocessing.normalize(R, norm='l2')
Y = preprocessing.normalize(np.cross(Z, V), norm='l2')
X = np.cross(Y, Z)
Rpos = satellite2.at(time).position.km.T

#check: LVLH coordinate frame at n events
fig = plt.figure(figsize=[10, 8])  # [12, 10]
ax  = fig.add_subplot(1, 1, 1, projection='3d')
axis_length=20
for i in range(0,5):
    x, y, z = R[i,:]
    u, v, w = X[i,:]
    ax.quiver(x, y, z, u, v, w, length=axis_length, color='red')
    u, v, w = Y[i,:]
    ax.quiver(x, y, z, u, v, w, length=axis_length, color='blue')
    u, v, w = Z[i,:]
    ax.quiver(x, y, z, u, v, w, length=axis_length, color='green')

Координаты LVLH

#construct the rotation matrix at time 0
RM = np.array([X[0,:],Y[0,:],Z[0,:]]).T
#view vector PoV = R_sat - R_ref
PoV = Rpos[:] - R
#rotate PoV to LVLH coordinate system
PoV_LHLV = RM.dot(PoV[0,:])
#go to spherical CS
r1, el1, az1 = sf.to_spherical(PoV_LHLV)

# Plot the view angles in polar plot.
plt.figure()
ax = plt.subplot(111, projection='polar')
ax.set_rlim([-90, 90])
ax.set_theta_zero_location('N')
ax.set_theta_direction(1)
ax.set_title('Visibility of satellite2 form satellite PoV', y=1.1)
ax.plot(az1, el1*degs, 'r+')

Угол обзора

Красивый!!! +n!