Координаты спутника по TLE

[Archi]

Здравствуйте.! Начал писать программу, которая считает координаты сателлита по TLE. Решил использовать алгоритм SGP4 http://www.celestrak.com/NORAD/documentation/spacetrk.pdf. Первый вопрос, который возник - есть дефинированы константы q0 и s, но не сказана их величина, и в google не могу найти:
qo = parameter for the SGP4/SGP8 density function
s = parameter for the SGP4/SGP8 density function
Можете подсказать, где их достать?

[Sleepwalker]

Смотрите раздел USERS GUIDE, CONSTANTS, AND SYMBOLS (стр 76) того документа

[Archi]

Да, но дело в том, что там не написано, какую цифру ставить в q0 и s.

[Sleepwalker]

Вот же

[Archi]

Ok, просто у него название не q0 а QOMS2T, поэтому не понял. А следущая вещь "For perigee between 98 kilometers and 156 kilometers, the value of the constant s used in SDP4 is changed to" - как можно узнать перигей? Это из TLE файла Аргумент перицентра в градусах, умножить на 111.12 km, или что-то другое?

[Sleepwalker]

Нет, смотрите реализауию на фортране, сразу все вопросы отпадут.
Например 16 стр.
А кстати зачем кодить методику самому? Все же уже сделано на многих языках.

[Archi]

Fortran не самый удобный язык. А насчет кодирования - я пишу дипломную работу по расчету координат сателлита по TLE файлам, поэтому нельзя брать чужой код. Надо самому написать по формулам из SGP4. Зато и нужна помощь.

PERIGE=(AODP*(1.-EO)-AE)*XKMPER - тут перигей но откуда такая формула?

[Sleepwalker]

формула высоты перигея (А(1-e) - 1)*6378.137 (A в радиусах Земли)

[Archi]

Хмм, перигей получается странным - 3.2387, 5.2247... А данные из TLE - их наверно надо конвертировать, например, частота обращения измеряется в оборотах в день, а для большой полуоси надо радианы в секунду. Какие там велечены надо конвертировать?

[Sleepwalker]

Смотрите раздел Driver (стр. 70)

[Archi]

1) А Bstar надо менять? Если его размерность (радиус земли)^(-1) ?
2) насчет эпохи - я знаю как высчитать сколько дней прошло от эпохи J200 (1 января 2000 года) как можно высчитать сколько дней прошло от 1950 года?
        double dwhole = 367 * year - ((int)7*(int)year+((int)(month+9)/12) / 4) +((int)275 * month / 9 ) + day - 730531.5;
        double dfrac = (hour +  (double)minute/60 + (double)second/3600)/24;
        d = dwhole + dfrac; //дни от 2000 года
3) И как из дня, который читается из TLE, достать все части по отдельности - месяц, день, час, минута, секунда?

[Tau]

1) А Bstar надо менять? Если его размерность (радиус земли)^(-1) ?

Не надо.

2) насчет эпохи - я знаю как высчитать сколько дней прошло от эпохи J200 (1 января 2000 года) как можно высчитать сколько дней прошло от 1950 года?

 Зачем в SGP4 эпоха 1950? Она там нигде не фигурирует.

3) И как из дня, который читается из TLE, достать все части по отдельности - месяц, день, час, минута, секунда?

Проще всего из целого количества дней получить день и месяц, используя цикл, начиная от первого января. Но любители формул сначала могут перейти к юлианской дате (юлианская дата 1 января данного года + количество дней), а потом выделить из юлианской даты день и месяц. Формулы есть хотя бы в Википедии.
Дробная часть дня - это время в долях суток. Умножить на 86400  (= время в секундах) и получившиеся число разложить на часы, минуты, секунды.

[Archi]

Там есть такая функция:
The function subroutine THETAG is passed the epoch time exactly as it appears on the input element cards.1
The routine converts this time to days since 1950 Jan 0.0 UTC, stores this in the COMMON E1,
and returns the right ascension of Greenwich at epoch (in radians).

[Tau]

THETAG вызывается процедурой deep-space, которая учитывает лунно-солнечные возмущения. Используется в модели SDP4 и SDP8. В SGP4 не используется. И, кстати, из-за того, что в SGP4 нет учета лунно-солнечных возмущений, в случае высокорбитальных ИСЗ,  для которых эти лунно-солнечные возмущения сравнимы с возмущениями от несимметричности геопотенциала, результат вычислений будет очень грустный. IMHO.

[Archi]

Тогда если t = 4133.83, а t0 = 4126.75, то там где в формулах (t-t0), будет 7.08, или еще надо конвертировать?

А остальные данные из TLE - так правильно как я сделал - Наклонение, Долгота восходящего узла, Аргумент перицентра, Средняя аномалия - из градусов в радианы, Частота обращения из оборотов в день в радианы в секунду, Эксцентриситет и Коэффициент торможения B* - не трогал.
А и для конвертации в радианы использовал готовую Java функцию, она возвращает радианы в рубежах -2pi ... 2pi, а в Spacetrk файле есть функции ACTAN и FMOD2P, которые возвращают радианы в рубежах  0 ... 2pi. Это обязательно делать?

И еще одна проблема - перигей все не получается. Вот пример:
a = 7253.10 km,
1-e0 = 0.9981,
XKMPER = 6378.135 km.

Perigey = (a*(1-e0)-1)*XKMPER = 7239.02 * 6378.135 = 46171446.8277 - огромное число, и близко нет 156 km

[Sleepwalker]

формула высоты перигея (А(1-e) - 1)*6378.137 (A в радиусах Земли)

[Archi]

Если радиус земли 6371 km, тогда надо делить A / радиус = 7253.10 / 6371 = 1.138 ?

[Sleepwalker]

6378.135
Но большая полуось должна была быть рассчитана уже безразмерной.

[Archi]

Большая полуось не в километрах? Если:
ke = 631.3481,
n0 = 14.0488 обороты в день = 0.0010 радианы в секунду

a1 = (ke/n0)^(2/3) = 7255.07 km, ее надо еще 7255.07/6378.135 = 1.1374. И во всех формулах, где есть большая полуось использвать эту цифру - delta1, a0, delta0, a0''?

[Sleepwalker]

631.3481 это откуда?
стр. 76: XKE=.743...
Среднее движение в радианах в минуту