Местоположение спутника

[feanorko]

Есть небольшая проблемка, а именно над высчитать координаты положения спутника (x,y,z)
(0 - центр земли, OX - гринвич, OZ - Полюс.).
Но проблема в том что входящие данные никак не мгу связать с формулами расчета (их просто напросто не дали (((((((((().
Может подскажеть сами формулы, иль место где их искать. (ГуГЛ не предлогать )
(Высчитать координаты имея:

     уход шкалы времени спутника
     долгота первого внутри суток восходящего узла
     наклонение орбиты
     эксцентриситет
     аргумент перигея орбиты спутника
     время прохождения восходящего узла орбиты
     драконический период
     скорость изменения драконического периода
)(альманах спутника)
Спутник системы ГЛОНАСС.

Заранее благодарен. Надеюсь поможете, а то я отчался

[Tau]

Имеем:
е - эксцентриситет
i - наклонение орбиты
Om - долгота восходящего узла
w - аргумент перигея
Не хватает большой полуоси и аргумента широты.

1) Вычисляем большую полуось a из периода.
n = 2 * Pi / Period ' собственное движение
fm = 398601.3 '  km^3/sec^2 для Земли
a^3 = fm / n^2
2) В момент прохождения спутником восходящего узла его аргумент широты равен 0.
U = 0
v = U - w '  истинная аномалия
r = a * (1 - e ^ 2) / (1 + e * Cos(v))  ' расстояние до спутника
3) Вычисляем вектор спутника в экваториальной системе координат
r(1) = r * (Cos(U) * Cos(Om) - Sin(U) * Sin(Om) * Cos(i))
r(2) = r * (Cos(U) * Sin(Om) + Sin(U) * Cos(Om) * Cos(i))
r(3) = r * Sin(U) * Sin(i)
4) Экваториальная система координат отличается от гринвической поворотом вокруг оси Z на угол равный гринвическому звездному времени.
Соответственно, вычисляем звездное время, матрицу поворота, умножаем вектор на матрицу и получаем вектор в гринвической системе координат на момент прохождения восходящего узла орбиты (потому что u=0).

Если нужно вычислить положение на любой момент времени, то после пунта 3 вычисляем еще и вектор скорости, медленно и печально интегрируем с учетом возмущений Луны, Солнца, несферичности Земли (или в серьезном случае - с учетом всех гармоник геопотенциала) до нужного момента и только потом переходим в гринвическую систему.
В этом лучае из литературы может пригодится:
"Справочное руководство по небмеху и астродинамике" под ред. Г.Н. Дубошина http://www.astrolib.ru/library/42.html
"Динамика ИСЗ" Т.В.Бордовициной http://www.astro.tsu.ru/ISZ/ISZ.pdf

Что касается таких параметров как уход шкалы времени спутника и скорость изменения драконического периода, то это в данном случае лишняя цифирь, IMHO.

[feanorko]

Огромное спс, уже становится яснее, но все же я нубище в астрономии(хотя после,надеюсь, выполнения данного ТЗ у мя появятся хотя бы какие то познания о движении спутников.) знаюий что "земля круглая" и она крутится вокруг солнца а не наоборот
Tau
тут такое делоя я не сказал что есть еще 1 переменная - номер дня...
почитав немного литературы, вродь умной литературы, считаю что это номер дня в 4-х летнем периоде (там чет о високосном дне.....). и я так понял данные спутник передает каждый день в одно и то же время. думаю это 0-00.

Думаю имея данную точку времени будет немного легче выяснить x-y-z для спутника который передал сию информацию.

есл не сложно - помоги. а то серьезно - завал, я уже совсем начинаю отчаиваться (((

[Tau]

feanorko, предлагаю определиться с целью задачи, если не сложно. Как сформулирована задача? Откуда берутся входные данные? Что должно быть на выходе?

[feanorko]

такс, Входящие данные - Альманах спутников.
Примерно выглядит так:

     уход шкалы времени спутника=0.011444091796875
     долгота первого внутри суток восходящего узла=1.282938241958618164
     наклонение орбиты=1.131630182266235352
     эксцентриситет=0.00077056884765625
     аргумент перигея орбиты спутника=0.261256098747253418
     время прохождения восходящего узла орбиты=10681875
     драконический период=40543996.09375
     скорость изменения драконического периода=0.1220703125
     номер дня=716

необходимо  найти местоположение спутника (координаты x,y,z) (0 - центр земли, OX - гринвич, OZ - Полюс.) в любую точку времени (отесть написал чтот вродь: 19:50:11 10.01.2010 и программа будет выводить координаты. Надеюсь есл в формулах не погрязну )

[Tau]

Время прохождения узла и драконический период в альманахе даны в миллисекундах, причем время скорее всего московское зимнее.
Все углы в радианах.

Здесь http://www.glonass-ianc.rsa.ru/i/glonass/ICD-2002r.pdf расписан алгоритм расчета, начиная со страницы 45. По-моему все подробно разжевано.

[Tau]

Тут новая документация к Глонассу (2008 год)  http://www.aggf.ru/gnss/glon/ikd51ru.pdf , несколько опечаток исправлено. Смотреть со стр. 56, конкретно расчет по данным альманаха - со стр.64.
Кстати, а таблица соответствия номера дня и даты у вас есть?
Еще вопрос: у вас точно в альманахе долгота восходящего узла, а не гринвическая долгота восходящего узла? Если точно, то пункт 2 раздела П.3.2.2. пропустите, все равно у вас двух поправок нет.

[feanorko]

да прочитал, насчет расчета таблицы - придется начальство напрячь завтра пусть выдают

вопрос: у вас точно в альманахе долгота восходящего узла, а не гринвическая долгота восходящего узла?

завтр уточню.

[feanorko]

ну есл верить документации то у меня в альманахе - гринвеческая долгота восх. угл.ъ
а насчет таблиц они напрягли мол в ручную считай, начиная от 1го января високосного года в периоде.

[Tau]

Запрограммила я расчет вектора по альманаху. Без интегрирования, ага. С учетом лишь второй зональной гармоники (J2), это где тригонометрические ряды расписаны в интерфейсном документе. Ну что сказать... Не знаю, какая точность вам нужна, может быть и такая сгодится.

а насчет таблиц они напрягли мол в ручную считай, начиная от 1го января високосного года в периоде.

Вам нужно найти алгоритм расчета юлианской даты  и искать разницу между всеми юлианскими датам четырехлетия и нулевой датой четырехлетия. В текущем четырехлетии Датой=0  будет 31.12.2007, от нее и плясать.
Номер дня = 716 это 16.12.2009, соответственно данные приведенного вами альманаха относятся к 02:58:01.88 ДМВ 16 декабря 2009 года.
Проверим, как вычисляется прогноз.
Вычисляем инерциальный вектор на 12:28:00 UT 16 декабря 2009 г. по алгоритму альманаха для тех начальных данных, что вы привели.
X  = -20795.9907667501
Y =  1030.95151488405
Z =  14709.0506172043
С чем бы его сравнить? В любой задаче должна быть возможность проверки.
Мы знаем, что 16 декабря в 02:58:01.88 ДМВ какой-то Глонасс прошел через восходящий узел. Долгота восходящего узла была равна OM = L + Sg =157.78, что похоже на правду, не так уж много у группировки Глонассов орбитальных плоскостей.  Берем орбиты Глонассов с восходящим узлом ~ 157 градусов из каталога NORAD, какой-нибудь софт для расчета векторов по TLE и ищем. Вашему альманаху соответствует спутник с норадовским номером 28915, правда 16 декабря он прошел через восходящий узел на 2 секунды позже, чем указано в вашем альманахе.
Посчитаем его инерциальный вектор на 12:28:00 ДМВ 16 декабря 2009 г. по данным Норада.
X' = -20795.0529
Y' = 1044.2734
Z' = 14706.6287
Расхождения
X'-X = 0.9379
Y'-Y = 13.3219
Z'-Z = 2.4216
Ошибки эти можно списать на то, что в NORAD-е не особо точные данные, или на то, что в алгоритме  не все возмущения  учтены, а только J2.

Я даже боюсь спросить, для чего нужно использование альманахов Глонассов. Мне казалось, что их вообще для красоты печатают, очень уж формат данных специфичен.

[feanorko]

сори, возможно я все таки чего то не понял но в интерфесном контрольном доке. от 2008 года.
на стр. 65 в самом низу написано - "набор параметров орбиты для каждого НКА задан в гринвичской геоцентрической системе OXYZ на момент времени...."

а расчет (стр. 67) ведется уже в абсолютной геоцентрической системе...

Проблема в том что значения в абсолютной мне  вообще не нужны, а необходимо всего лишь выгрызть местоположение спутника в ГРИНВИЧСКОЙ системе....

я нарное сума сойду сам и Вас Tau  так же сведу ((((((((((

есть ли вариант как нить попроще найти эти  координаты ? нежели пытатся в код переложить 5 стр кода, в которых я и так уже потерялся

[Tau]

У нас есть две системы координат: инерциальная невращающаяся (ось X направлена в точку весеннего равноденствия) и гринвическая вращающаяся (ось X направлена в гринвический меридиан). Эти две системы отличаются поворотом вокруг оси Z, потому что в каждый момент времени угол между осями X обеих систем равен текущему гринвическому звездному времени. Это время в каждый момент разное. Одна система координат вращается в другой системе координат, раз в сутки они совпадают.
Орбита - это, грубо говоря, эллипс, пространственно-ориентированный в невращающейся системе координат. Параметры орбиты принято задавать в невращающейся системе, но в альманахе заданы во вращающейся (гринвической). НО в гринвической системе параметры все те же самые, что и в инерциальной, КРОМЕ долготы восходящего узла. В невращающейся системе долгота - почти постоянная величина, очень медленно изменяется, а во вращающейся системе долгота изменяется со скоростью изменения звездного времени. Связь между долготой в инерциальной СК и гринвической СК такая:
Om = L + Sg
Om - долгота восходящего узла в инерциальной СК
L - географическая долгота восходящего узла в гринвической СК
Sg - текущее гринвическое звездное время

значения в абсолютной мне вообще не нужны, а необходимо всего лишь выгрызть местоположение спутника в ГРИНВИЧСКОЙ системе....

Проблема не в том, что вам нужно . Вычисление положения спутника ВСЕГДА проводится в невращающейся системе координат, даже если вы будете использовать самую примитивную модель - модель невозмущенного движения. И это не блажь, это математически наиболее удобно, потому что спутник с физической точки зрения движется в невращающейся системе координат, какая там под ним планета и с какой скоростью она вращается вокруг своей оси - спутнику по барабану. Вращающаяся система координат введена для удобства наблюдателя на поверхности Земли, ну и еще кое для чего, что в данном случае неважно. Так что, как бы вам ни хотелось, чтобы было попроще, от переходов из одной системы координат в другую, и наоборот, не отвертеться.

нежели пытатся в код переложить 5 стр кода, в которых я и так уже потерялся

А что вы хотели? Задачи небесной механики в тяп-ляп не решаются, если вам серьезно нужно решить, а не какую-нибудь рисовалку закодить, в которой масштаб отрисовки съест все погрешности. И то 5 страниц формул - это упрощенный вариант. По-хорошему нужно все-таки интегрировать .

[feanorko]

Галилей бы мя собственорочно удушил бы , я забыл что "она крутится".

[feanorko]

а как тогда найти истинное звездное время на гринвичскую полночь, даты N0 к которой относится время ti(стр. 68 интерфесном контрольном доке. от 2008 года.)
Там приведена переменная как S0

Тоесть выше написано что Номер дня = 716 это 16.12.2009 вот а как найти на этот день S0 ?

[Tau]

Для вычисления звездного времени нужно сначала вычислить юлианскую дату JD.
Прицепляю три страницы из "Космической геодезии" Крылова, там есть формулы.

[Tau]

Вот еще полезная документина: http://law.rufox.ru/view/standarti/1672.htm

[feanorko]

а есл сделать проще (точность теряется кнешн. но сейчас не точность главное.)
взять из звездного календаря к примеру Nautical Almanac а именно Observers' Calendar for Cambridge for 2008(беру для кембрижа т-как ЖМТ - 0) на январь месяц (а именно 1е января)... там есть юлианская дата и LMST ну а дальше от этой даты и плясать с учетом разницы во времени между двумя системами исчисления (насколько я понял  она равна примерно 9,71-9.75 с/час иль примерно 0,0027с/с) собственно эт и будет примерное звездное время.....?

вот только литл порблем. там в расчетах время в радианах задано.....

[Tau]

Во всей вашей задаче вычисление юлианской даты и звездного времени - самый простой этап, а вы и от него бежите .

1) Вычисление юлианской даты
Пусть дано Year - год, Mon - месяц, Day - день
   

If Mon <= 2 Then
        Mon = Mon + 12
        Year = Year - 1
    End If
    var1 = Year \ 400 - Year \ 100 + Year \ 4
    var2 = 365 * Year - 679004
    MD = var1 + var2 + 306001 * (Mon1 + 1) \ 10000 + Day

\ - деление нацело
MD - модифицированная юлианская дата

2) Вычисление среднего звездного гринвического времени на полночь даты MD
 

A1 = 24110.54841
 A2 = 8640184.812
 A3 = 0.093104
 A4 = 0.0000062
 T0 = (MD - 51544.5) / 36525
 Sg = A1 + A2 * T0 + A3 * T0 ^ 2 - A4 * T0 ^ 3
 Sg - звездное время в секундах времени.
 S0 = Sg *15 / 3600 ' в градусах

Пример: 14 января 2010
MD = 55210
Sg=  891204.250330178 sec
S0=  3713.35104304241 deg = 113.351043042409 deg

3) Для вычисления истинного звездного времени нужно вычислить нутацию по прямому восхождению и знать поправку времени. Это можно сейчас и пропустить (а то вы совсем перепугаетесь), среднее время в данный период отличается от истинного ~ 1 секунду.

время в радианах задано

Время измеряется в часах, минутах, секундах. 24 часа = 86400 секунд.
Но в то же время часовая мера и угловая взаимосвязаны:
24 часа = 360 градусов
1 час = 15 градусов
1 минута = 15 угловых минут
1 секунда = 15 угловых секунд.
Соответственно, если время задано в радианах, то переводите в градусы ( *180/pi), а затем, например, в часы, делением на 15.

[feanorko]

еще сопутсвующий вопросик
в вышеупомянутой документации по глонасу на стр 67 есть кусочек формулы

может я не понимаю чегот но вродь ее не реально посчитать.....

ну разве что e - эт не экспонента....

я все прально понял ?

[Tau]

e - это эксцентриситет, один из параметров орбиты. Дается в альманахе.

эксцентриситет=0.00077056884765625

Это принятые обозначения:
e - эксцентриситет (дальше в документе возмущенный эксцентриситет обозначается как эпсилон греческая)
i - наклонение орбиты
a - большая полуось
Om (Омега большая греческая)- долгота восходящего узла
w (омега малая греческая) - аргумент перигея
u - аргумент широты