Отслеживание спутника Python с помощью spg4, pyephem- позиции не совпадают

Я пытаюсь написать базовый скрипт на Python, который будет отслеживать данный спутник, определенный с помощью tle, из заданного местоположения. Я не асто / орбитальный человек, но пытаюсь стать умнее в этом.

Я сталкиваюсь с проблемой, когда различные модели, которые я использую, дают мне очень разные ответы о положении. Я попытался с помощью: PyEphem SPG4 предсказать (системный вызов exec из сценария)

Спутники, с которыми я тестирую, - это МКС и directv10 (один фиксированный, один движущийся с отслеживанием через Интернет, доступный для проверки):

0 Direct10
1 31862U 07032A   13099.15996183 -.00000126  00000-0  10000-3 0  1194
2 31862 000.0489 046.9646 0000388 001.7833 103.5813 01.00271667 21104
0 ISS
1 25544U 98067A   13112.50724749  .00016717  00000-0  10270-3 0  9148
2 25544  51.6465  24.5919 0009906 171.1474 188.9854 15.52429950 26067

Я изменил источник прогноза, чтобы дать мне местоположение ECI, чтобы я мог использовать его, чтобы узнать реальное местоположение. У меня также есть это, чтобы дать az, el, диапазон для использования, чтобы проверить наблюдения. Я использую spg4 для получения реального местоположения. Для наблюдаемого местоположения я использую PyEphem.

Я получаю позицию ECEF от spg4 с:

def get_real(epoch, sv):
satellite = twoline2rv(sv.tle1, sv.tle2, wgs84)

#epoch = time.time()
obsTime = datetime.datetime.utcfromtimestamp(epoch)
position, velocity = satellite.propagate(obsTime.year,
                                         obsTime.month,
                                         obsTime.day,
                                         obsTime.hour,
                                         obsTime.minute,
                                         obsTime.second)


x = position[0]
y = position[1]
z = position[2]

x *= 1000
y *= 1000
z *= 1000

Мой код для наблюдений на основе пифема:

def get_ob(epoch, sv, obsLoc):
site = ephem.Observer()
site.lon = str(obsLoc.lat)   # +E -104.77 here
site.lat = str(obsLoc.lon)   # +N 38.95   here
site.elevation = obsLoc.alt # meters    0 here
#epoch = time.time()
site.date = datetime.datetime.utcfromtimestamp(epoch)

sat = ephem.readtle(sv.name,sv.tle1,sv.tle2)
sat.compute(site)

az       = degrees(sat.az)
el       = degrees(sat.alt)
#range in m
range    = sat.range
sat_lat  = degrees(sat.sublat)
sat_long = degrees(sat.sublong)
# elevation of sat in m
sat_elev = sat.elevation

#TODO: switch to using az,el,range for observed location calculation
#x, y, z    = aer2ecef(az,el,range,38.95,-104.77,80 / 1000)
x,y,z  = llh2ecef(sat_lat, sat_long, sat_elev)

Преобразование llh2ecef:

def llh2ecef (flati,floni, altkmi ):
#         lat,lon,height to xyz vector
#
# input:
# flat      geodetic latitude in deg
# flon      longitude in deg
# altkm     altitude in km
# output:
# returns vector x 3 long ECEF in km

dtr =  pi/180.0;

flat  = float(flati);
flon  = float(floni);
altkm = float(altkmi);

clat = cos(dtr*flat);
slat = sin(dtr*flat);
clon = cos(dtr*flon);
slon = sin(dtr*flon);

rrnrm  = radcur (flat);
rn     = rrnrm[1];
re     = rrnrm[0];

ecc1    = ecc;
esq1    = ecc1*ecc1

x      =  (rn + altkm) * clat * clon;
y      =  (rn + altkm) * clat * slon;
z      =  ( (1-esq1)*rn + altkm ) * slat;

return x,y,z

aer2ecef:

def aer2ecef(azimuthDeg, elevationDeg, slantRange, obs_lat, obs_long, obs_alt):

#site ecef in meters
sitex, sitey, sitez = llh2ecef(obs_lat,obs_long,obs_alt)

#some needed calculations
slat = sin(radians(obs_lat))
slon = sin(radians(obs_long))
clat = cos(radians(obs_lat))
clon = cos(radians(obs_long))

azRad = radians(azimuthDeg)
elRad = radians(elevationDeg)

# az,el,range to sez convertion
south  = -slantRange * cos(elRad) * cos(azRad)
east   =  slantRange * cos(elRad) * sin(azRad)
zenith =  slantRange * sin(elRad)


x = ( slat * clon * south) + (-slon * east) + (clat * clon * zenith) + sitex
y = ( slat * slon * south) + ( clon * east) + (clat * slon * zenith) + sitey
z = (-clat *        south) + ( slat * zenith) + sitez

return x, y, z

Когда я сравниваю и строю местоположения на трехмерном глобусе (используя позиции ecef), я получаю ответы повсюду. Предсказанная позиция eci (преобразованная в ecef) соответствует тому, что я вижу на веб-сайтах отслеживания ISS (http://www.n2yo.com/?s=25544)

Результат от get_real () далеко от масштаба и местоположения. Результат от get_ob () имеет правильный масштаб, но неправильный адрес на глобусе

пример результатов:

на основе прогноза:

sv: ISS predict observed response    @ epoch: 1365630559.000000 : [111.485527, -69.072949, 12351.471383]
sv: ISS predict aer2ecef position(m) @ epoch: 1365630559.000000 : [4731598.706291642, 1844098.7384999825, -4521102.9225004213]
sv: ISS predict ecef position(m) @ epoch: 1365630559.000000 : [-3207559.6840419229, -3937040.5048992992, -4521102.9110000003]
sv: ISS predict ecef2llh(m)      @ epoch: 1365630559.000000 : [-41.67839724680753, -129.170165912171, 6792829.6884068651]
sv: Direct10 predict observed response    @ epoch: 1365630559.000000 : [39.692138, -49.219935, 46791.914833]
sv: Direct10 predict aer2ecef position(m) @ epoch: 1365630559.000000 : [28401835.38849232, 31161334.784188181, 3419.5400331273049]
sv: Direct10 predict ecef position(m) @ epoch: 1365630559.000000 : [-9348629.6463202238, -41113211.570621684, 3419.8620000000005]
sv: Direct10 predict ecef2llh(m)      @ epoch: 1365630559.000000 : [0.0046473273713214715, -102.81051792373036, 42156319.281573996]

на основе Python:

sv: ISS ephem observed response    @ epoch: 1365630559.000000 : [344.067992722211, -72.38297754053431, 12587123.0][degrees(sat.az), degrees(sat.alt), sat.range]
sv: ISS ephem llh location(m)      @ epoch: 1365630559.000000 : [-41.678271938092195, -129.16682754513502, 421062.90625][degrees(sat.sublat0, degrees(sat.sublong), sat.elevation]
sv: ISS ephem xyz location(m)      @ epoch: 1365630559.000000 :[-201637.5647039332, -247524.53652043006, -284203.56557438202][llh2ecef(lat,long,elev)]
sv: ISS spg84 ecef position(m) @ epoch: 1365630559.000000 : [4031874.0758277094, 3087193.8810081254, -4521293.538866323]
sv: ISS spg84 ecef2llh(m)      @ epoch: 1365630559.000000 : [-41.68067424524357, 37.4411722245808, 6792812.8704163525]
sv: Direct10 ephem observed response    @ epoch: 1365630559.000000 : [320.8276456938389, -19.703680198781303, 43887572.0][degrees(sat.az), degrees(sat.alt), sat.range]
sv: Direct10 ephem llh location(m)      @ epoch: 1365630559.000000 : [0.004647324660923812, -102.8070784813048, 35784688.0][degrees(sat.sublat0, degrees(sat.sublong), sat.elevation]
sv: Direct10 ephem xyz location(m)      @ epoch: 1365630559.000000 :[-7933768.6901137345, -34900655.02490133, 2903.0498773286708][llh2ecef(lat,long,elev)]
sv: Direct10 spg84 ecef position(m) @ epoch: 1365630559.000000 : [18612307.532456037, 37832170.97306267, -14060.29781505302]
sv: Direct10 spg84 ecef2llh(m)      @ epoch: 1365630559.000000 : [-0.019106864351793953, 63.80418030988552, 42156299.077687643]

Az, el и range не совпадают между двумя наблюдениями. Позиции не соответствуют «истинным» местам. (Lat и long делают, но высота не делает после преобразования ecef2llh.

По сравнению с веб-трекером, я заметил, что "истинные" местоположения будут совпадать с сайтом. Для directv10 pyEphem соответствует азимуту и высоте, но не для МКС

Когда я рисую их на земном шаре, предсказуемое «истинное» местоположение находится в правильном месте - соответствует сайту трекера). Положение spg84 ecef (которое, как мне показалось, должно совпадать с прогнозируемым, находится на другой стороне земного шара. Предсказанное «наблюдаемое» местоположение находится близко к местоположению spg84. PyEphem полностью отключен по высоте и не отображается (слишком низко, внутри земли).

Итак, мой вопрос, где я неправильно использую модели python? Насколько я понимаю, вызов spg84 spreadate () должен возвращать точное положение спутника в метрах. Я бы подумал, что должно соответствовать предсказанной позиции после преобразования eci2efec. Я также ожидал бы, что совпадение затем llh2ecef () при использовании sat.sublat, sat.sublong, sat.elevation.

Как я уже сказал, я новичок во всех вещах на орбите, поэтому я уверен, что совершил простую математическую ошибку или что-то в этом роде. Я пробовал гуглить и искать ответы, примеры и учебные пособия, насколько это было возможно, но пока ничего не помогло (я пробовал несколько методов ecef2llh и llh2ecef, чтобы попытаться устранить эти ошибки.

Будем весьма благодарны за любые предложения, советы, указатели в правильном направлении. Я могу опубликовать / отправить полный код, который я использую, если это будет полезно для кого-то. Я постарался сделать так, чтобы я разместил здесь важные части и не хотел делать этот (уже очень) длинный пост.

Спасибо за помощь.

Аарон

ОБНОВИТЬ:

Я нашел хотя бы часть проблемы. spg84.propagate () возвращает местоположение в ECI, а не в ECEF. Быстрый запуск через eci2ecef, и он идеально согласуется с предсказанным ответом.

Кажется, я всегда нахожу решения после публикации справки;)

Теперь нужно выяснить, что происходит с локациями наблюдателя. Это сводится к следующему: Как я могу взять результат из pyEphem.compute () и получить позицию ecef для спутника? Предпочитаю делать это с аз, эл, значениями диапазона, а не широтой, долготой, высотой.

Я предполагаю ошибку в моем вызове aer2ecef.

Благодарю.

ОБНОВЛЕНИЕ 2:

Получил наблюдение, чтобы выровнять с «истинной» позицией. Похоже, у меня была проблема с юнитами. Рабочий код:

az       = degrees(sat.az)
el       = degrees(sat.alt)
#range in km
range    = sat.range
sat_lat  = degrees(sat.sublat)
sat_long = degrees(sat.sublong)
# elevation of sat in km
sat_elev = sat.elevation


#x, y, z    = aer2ecef(az,el,range,obsLoc.lat,obsLoc.long,obsLoc.alt)
x,y,z  = llh2ecef(sat_lat, sat_long, sat_elev / 1000)

x *= 1000
y *= 1000
z *= 1000
return x,y,z

Теперь просто нужен метод aer2ecef (), чтобы вернуть правильную позицию ...

Ответы на вопрос(1)

Ваш ответ на вопрос