CC BY
192
УДК 681.783
А.Я. Юданин, Б.С. Могильницкий, А.С. Толстиков ФГУП «СНИИМ», Новосибирск
УТОЧНЕНИЕ ОРБИТ ИСЗ СИСТЕМЫ ГЛОНАСС/GPS
A.Ja. Yudanin, B.S. Mogilnizkii, A.S. Tolstikov
Siberian Scientific-Research Institute of Metrology (SSRIM)
4 Dimitrova UI., Novosibirsk, 630004, Russian Federation
ACCURATE DEFINITION ORBIT SATELITE TO SYSTEM GLONASS/GPS
In the work presents method accurate definition orbites artificial satellite Earth to base information pseudorange to navigation satellites.
Точность вычисления эфемерид НС (навигационных спутников) -основной фактор, определяющий эффективность работы систем глобального позиционирования. Поэтому задача уточнения элементов орбит (совместно с ПВЗ, а также, возможно, другими параметрами модели навигационных измерений) является весьма актуальной. В связи с бурным развитием координатно-временных технологий обеспечения системы ГЛОНАСС в последнее время появилась возможность решения этой задачи на основе массивов данных беззапросных измерений псевдодальностей.
НС регулярно излучает сигналы, содержащие информацию о моменте излучения, которые принимаются беззапросными измерительными станциями (БИС), расположенными на поверхности Земли. Псевдодальностью D называется геометрическое расстояние между точкой пространства, в которой находился НС в момент излучения сигнала (W) и положением БИС в момент приема (tnp ). В идеальном случае
D = I ГнсОиЗл) — f^ôlic(tnp) I — с * (tnp — ! (1)
где rHC(tm]I) - положение НС в момент излучения сигнала, R^c(tnp) -положение БИС в момент приема сигнала, с - скорость света.
Однако в реальности уравнение измерений имеет более сложный вид:
D e * (tnp — tH! I + At ) — daTM + v , (2)
где At- расхождение часов БИС и бортовых часов НС; daTM- запаздывание сигнала в атмосфере; v - случайная (немоделируемая) ошибка.
Расхождение часов на небольшом мерном интервале (7-15 суток) можно считать линейной функцией времени, а влияние атмосферы на прохождение сигнала в зависимости от времени суток, температуры воздуха, высоты НС над горизонтом и т. п. определяется по известным общепринятым моделям [1].
Координаты фазового центра антенны БИС в квазиинерциальной системе отсчёта (ИСК 2000) находятся по формуле [2]
R6Hc(t)=P*N*Z*n*Ro, (3)
где Ио - постоянные координаты БИС в земной (гринвичской ) системе;
Р и N - матрицы прецессии и нутации с известной зависимостью от времени.
Матрица звёздного времени 2 и матрица движения полюса П также известны как функции времени и ПВЗ.
Вектор положения БИС в ИСК 2000 на произвольную эпоху t равен:
R п (t ) = P-Ж-Z с
-1 о" ^о о
о о о о
о о, v1 о
о
с
Л
о о о о о -1
0
1 о
>R =
I R ПЕ
= ^П0 (0 + а(0 • Rna (0 + Хр (0 • Япх (0 + Ур (0 • Rm (0 , (4)
где явно выделены малые параметры (а, хр, yp) подлежащие оцениванию. Зависимость ПВЗ от времени целесообразно представлять в виде суммы небольшого числа баззисных функций с минимальными коэффициентами, подлежащими оценке:
a(t) = Xk ak*Pka(t); xp(t) = Ik bk*Pkx(t); yp(t) = Ik ck*Pky(t).
(5)
Движение центра масс НС по орбите описывается векторным дифференциальным уравнением второго порядка [1] г" + ц*г/г3 = w( t, r(t)), (6)
где \i - константа, aw- вектор « возмущающих » сил, обусловленных тяготением Луны и Солнца, нецентральностью гравитационного поля Земли, а также радиационным давлением на панели солнечных батарей НС.
Уравнение (6) можно преобразовать в систему уравнений первого порядка:
г”= у; у' = - ц*г/г3 + w, (7)
и решать численными методами (например, методом Адамса) с достаточно высокой точностью. Однако точность вычисления эфемерид НС можно существенно повысить, перейдя к переменным
х = r/r, М = [ г х у ], X = [ у х М ] - ц *х (8)
и соответствующей системе уравнений:
М' = М2* [ х х w ]/( |i + (А,*х)) г = х*М2/(ц.+(Х*х))
A,' = [wxM] + [vx М1 ] (9)
у = [ М х (X + ц*х) ]/М2
х' = [Мхх]*(ц + (Х*х))2/М4
Правые части первых двух векторных дифференциальных уравнений (9) пропорциональны малому возмущению w, поэтому М и X изменяются медленно, что уменьшает общую погрешность интегрирования. Формально в системе (9) девять скалярных переменных, однако они связаны тремя соотношениями :
(М*А,) = 0;(М*х) = 0;х2 = 1. (10)
В процессе численного интегрирования (тем же методом Адамса) значения М, X и x на каждом шаге можно “ подправлять “ с помощью (10). Предположим, что на начало мерного интервала (t = 0) известны приближённые значения координат и скорости НС. Тогда реальную орбиту r(t) можно представить в виде суммы « опорной » орбиты r0(t), которая рассчитывается по формулам (7) или (8), исходя из неточных начальных условий, и шести “частных производных орбиты по начальным условиям“, 3j(t), с неизвестными весами (8j). Аналоги частных производных орбиты определяются как разности между опорной орбитой и вариативными орбитами, которые вычисляются так же, но с небольшими отклонениями по каждой компоненте начальных координат и скоростей НС. Таким образом, уравнение измерений можно записать в виде:
I r„(t“”) + 1,68^Э;(С) - R°6nc(tinp) - a(t,)*R“6„c(t1np) - x(t1)*Rx6„c(t1np) -
- y(ti)*Ry6„c(tinp) I - Ai - Л, + dm = v(, (11)
где Ai = c^t^-tD - фактически полученные псевдодальности;
А0 - смещение псевдодальностей за счёт расхождения часов.
Предполагая, что случайная ошибка v имеет нормальное несмещённое распределение (возможное смещение ошибки эквивалентно А0), приходим к методу наименьших квадратов:
Н = E,n { I r„(t“) + XiVa/t"”) - R°6i,c(tinp) - a(t1)*Ra6„(t1np) -x(ti)*Rs6„(t,np) -
- y(ti)*Ry6„c(tinp) I - Aj - До + d,,.., r = min (12)
Суммирование производится по всем сигналам, принятым в обработку
(НС должен находиться в зоне видимости БИС и не слишком низко над горизонтом).
Выражение в фигурных скобках преобразуем к виду:
I ro(t“”) + EiVSjft") - R°6„c(tinp) - ... I - д, - До + d„M =
(13)
I r„(t“) + Е,%»Э№“) - R°6„c(t,np) - ■ ■ ■ 12 - (Ai + До - d„M )2
1 r„(t“) + Е,%»Э№“) - R°6„c(t,np) - ... I + Ai + До - d„M
{r„(t") - RWC)}2 + 2*{r„(t“) - RWOiMSi'Sj*^^“) -
a(t,)*R“6„c(t1np) -
- x(ti)*Rx6llc(tinp) - у(1,)*^б„с(1”Р)} + {Е,6^*Э/0 - a(ti)*R“6„c(tinp) - ... }2
(Aj daro)) 2*( Aj dax^i )*Ao (Ao)
Ä/ —————————————————————————————————————————————————————————————————————————————————————————
2 ф (Ai — daTM)
Возводя это выражение в квадрат и располагая слагаемые числителя в порядке возрастания степени параметров минимизации, получаем:
{(ГоаГЬИ^а,4’))2 - (Д, -(I™)2 }2 + 4*{( Го(1 ,")-К0бЖР))2 -(4, -ёа
»{[го(1,“ш)-Кб„Лпр)] »Р,65,*Э,(1 ,“)-а(11)*К“биЛпр)-■ ■ ■ ]-2*(А,Ч™)*До}+...
Н=2....................................................................
4 * (А1 - ёатм)2
«Невязка» Н является полиномом четвертой степени по параметрам минимизации 8}, ак, Ьк, ск, (см. (7)) и А0. Первое приближение 5j0, ак(), ... получаем, учитывая только члены до второй степени включительно, и минимизируя квадратичную форму стандартным образом. Затем представляем искомые параметры в виде 5| = 5/1 + 5/, ак = ак° + ак', ..., подставляем во все слагаемые (14) и образуем новую квадратичную форму по переменным бД ак' и т.д. Двух-трех итераций, как правило, вполне достаточно. Компьютерное моделирование показывает, что неточность орбиты НС данным способом можно сократить с десятков и сотен метров до сантиметров с помощью только одной БИС.
Программа моделирования вышеописанного способа уточнения орбит НС системы ГЛОНАСС на основе беззапросных измерений псевдодальностей включает в себя следующие элементы:
а) Точное аналитическое вычисление положения НС в любой момент времени, как эталон для проверки численного метода;
б) Численный метод (Адамса) интегрирования уравнений движения НС с достаточно большим (50 сек) шагом, а также “ разгонный “ алгоритм;
в) Метод вычисления положений НС в промежуточных точках (шаг две сек) излучения сигнала;
г) Имитация измерений псевдодальностей с учётом вращения Земли и нестабильности ПВЗ;
д) Задание случайных отклонений начальных параметров орбиты от истинных, вычисление опорной орбиты, производных орбиты по начальным параметрам и коэффициентов невязки;
е) Нахождение уточнённых параметров орбиты НС и ПВЗ путём минимизации невязки.
Для проверки состоятельности предлагаемого способа уточнения орбит задача решалась в упрощённом варианте: орбита НС считалась строго эллиптической, без возмущений; влияние атмосферы на прохождение сигнала, прецессия и нутация Земли, а также смещение полюса не учитывались; динамика ухода времени принималась линейной. Случайные ошибки «измерения» псевдо-дальностей имитировались в двух вариантах: равномерно распределённые в интервале (-0.5м, 0.5 м) и нормально распределённые с такой же дисперсией и (СКО = 0.3 м). Имитация измерения псевдодальностей базируется на уравнении
с *Д, = I Щизл + а(ъ) + А,) - г(1Гл) I, (15)
которое мы решаем относительно Л1 с учетом (3) и (4). Далее к Д: добавлялась случайная ошибка и в таком виде псевдодальность запоминалась (вместе с моментом излучения) для каждого принятого сигнала. Аналитический расчёт истинной орбиты, численное интегрирование уравнений движения НС и имитация измерений псевдодальностей проводились в одном цикле по мерному интервалу. Во втором цикле производилось вычисление опорной орбиты, производных орбиты по начальным условиям и насчитывались коэффициенты невязки. Число этих
2 3 4
коэффициентов равно (п + п + п + п ), где п - количество оцениваемых параметров (в нашей модели п = 6+2=8). Затем в каждой итерации минимизировалась квадратичная форма, элементы которой каждый раз пересчитывались в соответствии с предыдущим решением. Результаты моделирования подтверждают возможность существенного уточнения орбит НС системы ГЛОНАСС предлагаемым способом на основе данных измерений зашумлённых псевдодальностей.
Выявлено ограничение этого метода: он эффективен только в тех случаях, когда период обращения НС с неточной орбитой совпадает с истинным периодом. С другой стороны, влияние Луны и Солнца на орбиту НС не меняет период его обращения вокруг Земли, хотя форму и ориентацию орбиты может менять. Все вычисления производились с двойной точностью. Матрица системы линейных уравнений оказалась сильно вырожденной. Моделирование показало, что метод наименьших квадратов рассчитан на нормальный закон распределения ошибок.
При нормальном законе итерации улучшают первое приближение, а при равномерном распределении случайных ошибок результаты итераций были, как правило, существенно хуже первого приближения.
Разработанная оригинальная методика уточнения параметров орбит ИСЗ, на основе использования псевдодальностей, позволяет существенно повысить точность определения орбит с десятков и сотен метров до сантиметров, используя информацию всего от одного НС. Увеличение числа НС приведет к дальнейшему снижению погрешности в определении параметров орбит навигационных систем ГЛОНАСС/GPS.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Абалкин, В.К. и др. Справочное руководство по небесной механике и астродинамике / В.К. Абалкин, Г.Н. Дубошина. - М.: Наука. - 1971.
2. Шебшаевич, В.С. и др. Сетевые спутниковые радионавигационные системы /
В.С. Шебшаевич. - М.: Радио и связь. - 1993.
© А.Я. Юданин, Б.С. Могильницкий, А.С. Толстиков, 2008
