«СОВЕРШЕНСТВОВАНИЕ МЕТОДИКИ ТОЧНОГО ДИФФЕРЕНЦИАЛЬНОГО ПОЗИЦИОНИРОВАНИЯ С ИСПОЛЬЗОВАНИЕМ ГЛОБАЛЬНЫХ НАВИГАЦИОННЫХ СПУТНИКОВЫХ СИСТЕМ ...»

Министерство образования и наук

и Российской Федерации

Федеральное государственное бюджетное образовательное

учреждение высшего профессионального образования

«СИБИРСКАЯ ГОСУДАРСТВЕННАЯ ГЕОДЕЗИЧЕСКАЯ АКАДЕМИЯ»

(ФГБОУ ВПО «СГГА»)

На правах рукописи

Липатников Леонид Алексеевич

СОВЕРШЕНСТВОВАНИЕ МЕТОДИКИ ТОЧНОГО ДИФФЕРЕНЦИАЛЬНОГО

ПОЗИЦИОНИРОВАНИЯ С ИСПОЛЬЗОВАНИЕМ ГЛОБАЛЬНЫХ

НАВИГАЦИОННЫХ СПУТНИКОВЫХ СИСТЕМ

25.00.32 – «Геодезия»

Диссертация на соискание учёной степени кандидата технических наук

Научный руководитель – кандидат технических наук, доцент Сурнин Юрий Венедиктович Новосибирск –

ОГЛАВЛЕНИЕ

ВВЕДЕНИЕ

1 МЕТОДИКА ТОЧНОГО ДИФФЕРЕНЦИАЛЬНОГО

ПОЗИЦИОНИРОВАНИЯ

1.1 Общие сведения о методике

1.2 Математическая модель динамической системы и алгоритм методики точного дифференциального позиционирования

1.2.1 Измеряемая геометрическая дальность

1.2.2 Модель движения навигационного космического аппарата

1.2.3 Модель движения позиционируемого объекта

1.2.4 Измерительная информация

1.2.5 Система уравнений измерений

1.2.6 Стохастическая модель динамической системы

1.2.7 Рекуррентное оценивание параметров

1.3 Направления развития методики точного дифференциального позиционирования

2 АНАЛИЗ И СОВЕРШЕНСТОВАНИЕ МАТЕМАТИЧЕСКОЙ МОДЕЛИ

ДИНАМИЧЕСКОЙ СИСТЕМЫ И АЛГОРИТМА МЕТОДИКИ ТОЧНОГО

ДИФФЕРЕНЦИАЛЬНОГО ПОЗИЦИОНИРОВАНИЯ

2.1 Статистический анализ погрешностей математической модели ГНССизмерений

2.1.1 Корреляционный анализ остаточных невязок уравнений фазовых измерений

2.1.2 Анализ погрешностей интерполяции поправок бортовых шкал времени навигационных космических аппаратов

2.2 Совершенствование математической модели ГНСС-измерений

2.2.1 Способ уточнения функциональной модели фазовой псевдодальности........ 2.2.2 Способ уточнения априорной стохастической модели погрешностей ГНССизмерений

2.3 Усовершенствованная методика точного дифференциального позиционирования и её программная реализация

3 ЭКСПЕРИМЕНТАЛЬНОЕ ИССЛЕДОВАНИЕ ХАРАКТЕРИСТИК

УСОВЕРШЕНСТВОВАННОЙ МЕТОДИКИ ТОЧНОГО

ДИФФЕРЕНЦИАЛЬНОГО ПОЗИЦИОНИРОВАНИЯ

3.1 Постановка задачи и план экспериментального исследования методики точного дифференциального позиционирования

3.2 Результаты экспериментов и их анализ

3.2.1 Результаты экспериментов

3.2.2 Сравнительный анализ времени сходимости решения по исходному и усовершенствованному вариантам методики точного дифференциального позиционирования

3.2.3 Анализ пространственного распределения погрешностей оценок координат по методике точного дифференциального позиционирования

3.2.4 Анализ результатов применения уточнённой априорной стохастической модели погрешностей измерений

3.2.5 Анализ влияния типа кодовых псевдодальностей на время сходимости решения

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

3.3.1 Выводы по результатам экспериментального исследования

3.3.2 Рекомендации по практическому применению методики

3.3.3 Актуальные направления работы

ЗАКЛЮЧЕНИЕ

СПИСОК СОКРАЩЕНИЙ И УСЛОВНЫХ ОБОЗНАЧЕНИЙ

СПИСОК ЛИТЕРАТУРЫ

ПРИЛОЖЕНИЕ А (справочное) АПРИОРНЫЕ СТОХАСТИЧЕСКИЕ МОДЕЛИ ПРИЛОЖЕНИЕ Б (справочное) СВЕДЕНИЯ О ЗАДЕЙСТВОВАННЫХ ПРИЛОЖЕНИЕ В (обязательное) РЕЗУЛЬТАТЫ ПРИМЕНЕНИЯ ИСХОДНОГО

И УСОВЕРШЕНСТВОВАННОГО ВАРИАНТОВ МЕТОДИКИ

ТОЧНОГО ДИФФЕРЕНЦИАЛЬНОГО

ПОЗИЦИОНИРОВАНИЯ В СТАТИЧЕСКОМ РЕЖИМЕ.... ПРИЛОЖЕНИЕ Г (обязательное) РЕЗУЛЬТАТЫ ПРИМЕНЕНИЯ РАЗЛИЧНЫХ

АПРИОРНЫХ СТОХАСТИЧЕСКИХ МОДЕЛЕЙ

ПОГРЕШНОСТЕЙ ГНСС-ИЗМЕРЕНИЙ В МЕТОДИКЕ

ТОЧНОГО ДИФФЕРЕНЦИАЛЬНОГО

ПОЗИЦИОНИРОВАНИЯ

ВВЕДЕНИЕ

Актуальность темы исследования. Глобальные навигационные спутниковые системы (ГНСС) в настоящее время являются незаменимым средством выполнения многих видов геодезических работ, уникальным по совокупности таких характеристик, как уровень достигаемой точности, обширность геодезических построений и доступность. Методика точного дифференциального позиционирования (ТДП), более известная под англоязычным названием «Precise Point Positioning» (PPP), является одной из наиболее передовых методик математической обработки ГНСС-измерений. К наиболее важным свойствам методики можно отнести следующие:

высокую точность определения местоположения позиционируемого объекта в общеземной системе координат (среднеквадратические погрешности оценок координат от 1 до 3 см в статическом режиме), глобальную доступность и автономность определения координат (отсутствие непосредственной привязки к пунктам геодезических сетей), высокую степень автоматизации обработки ГНСС-измерений.

Актуальность задач развития и внедрения методики ТДП подчёркивается тем фактом, что соответствующие работы ведутся в рамках федеральной целевой программы «Поддержание, развитие и использование системы ГЛОНАСС на 2020 годы» (опытно-конструкторская работа «КФД-В») [17, 18, 26]. Возможность определения местоположения непосредственно в общеземной системе координат также придаёт методике особую практическую значимость в свете перехода к использованию геоцентрических систем координат в качестве государственных [40].

Одним из недостатков исходной методики ТДП, ограничивающих её применение, является необходимость производства продолжительных сеансов измерений (от нескольких часов до суток) для достижения указанного уровня точности.

Сокращение длительности сеансов ГНСС-измерений при сохранении высокой точности определения координат является актуальной задачей.

Степень разработанности проблемы. В России работы по теме точного дифференциального позиционирования ведутся в рамках создания функциональных дополнений Единой системы координатно-временного и навигационного обеспечения (ЕС КВНО) [49]. Значительный вклад в формирование концепции ЕС КВНО внесли Урличич Ю. М., Финкельштейн А. М., Ревнивых С. Г., Тестоедов Н. А., Данилюк А. Ю., Донченко С. И., Долгов Е. И., Демьянов Г. В., Макаренко Н. Л., Пешехонов В. Г., Красовский П. А., Белов С. А., Бутенко В. В. и другие. Авторами работ на русском языке и переводов с английского, посвящённых дифференциальному методу в ГНСС, и в частности методике ТДП, являются Антонович К. М., Виноградов А. В., Войтенко А. В., Дворкин В. В., Глухов П. Б., Голубев А. Н., Жигулин А. Ю., Карутин С. Н., Першин Д. Ю., Подкорытов А. Н., Щербаков А. С. и другие. Значительный вклад в создание и развитие методики ТДП внесли следующие зарубежные авторы: Bisnath S. B., Burgan J., Collins P., Dodson A. H., Geng J., Gerhatova L., Han S. C., Heflin M. B., Hefty J., Hroux P., Jefferson D. C., Jekeli C., Kouba J., Kwon J. H., Liu J., Lahaye F., Langley R. B., Laurichesse D., Meng X., Shi C., Teferle F. N., Watkins M. M., Webb F. H., Zumberge J. F. и другие. Несмотря на возросший интерес к данной теме и активизацию работ по совершенствованию методики ТДП задача повышения оперативности определения местоположения остаётся актуальной.

Цель исследования – совершенствование методики точного дифференциального позиционирования с использованием глобальных навигационных спутниковых систем.

Задачи исследования:

анализ возможностей и перспектив совершенствования методики точного дифференциального позиционирования и выявление факторов, влияющих на оперативность определения координат;

разработка и теоретическое обоснование способа повышения оперативности точного дифференциального позиционирования;

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

экспериментальная проверка состоятельности предложения по совершенствованию методики точного дифференциального позиционирования;

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

Объект и предмет исследования. Объект исследования – динамическая система, включающая спутники – навигационные космические аппараты (НКА) ГНСС, станции – наземные измерительные пункты (НИП), а также источники возмущающих воздействий, проявляющихся в результатах измерений. Предмет исследования – методика точного дифференциального позиционирования.

Научная новизна результатов исследования заключается в следующем: разработан новый способ компенсации погрешностей интерполяции поправок бортовых шкал времени в методике ТДП, особенностями которого являются простота, применимость для различных ГНСС, нулевое математическое ожидание предлагаемых поправок фазовых измерений, за счёт чего, теоретически, обеспечивается совместимость с различными способами разрешения и фиксации фазовых неоднозначностей на последующем этапе решения.

Теоретическая значимость работы заключается в оценке погрешности модели бортовых шкал времени и её связи со временем сходимости решения по методике ТДП, а также в уточнении классификации методов космической геодезии.

Показано, что методику, именуемую «Precise Point Positioning», следует относить к дифференциальному методу космической геодезии, а не к абсолютному, как это принято в настоящее время.

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

Методология и методы исследований. Методологическую базу исследования составили методы математического моделирования и статистической обработки результатов измерений, методы решения задач космической геодезии. Теоретическая база исследования: теория вероятностей и математической статистики, линейная алгебра, теория математической обработки геодезических измерений. Эмпирической базой исследования послужили: совокупность НИП, созвездие НКА GPS, результаты GPS-измерений, геофизические модели, эфемериды и поправки бортовых шкал времени НКА, координаты станций Международной ГНССслужбы (МГС).

Положения, выносимые на защиту:

а) новый способ компенсации погрешностей интерполяции поправок бортовых шкал времени спутников ГНСС и уточнённая априорная стохастическая модель погрешностей ГНСС-измерений в методике точного дифференциального позиционирования обеспечивают двукратное сокращение продолжительности сеансов GPS-измерений по сравнению с ранее достигнутой при сохранении высокого уровня точности определения координат;

б) рекомендации по практическому применению методики точного дифференциального позиционирования позволяют повысить надёжность и точность определения координат.

Соответствие диссертации паспорту научной специальности. Тематика диссертации соответствует паспорту научной специальности 25.00.32 – «Геодезия», разработанному экспертным советом ВАК Минобрнауки РФ, по пунктам: – «Геодезические (глобальные) навигационные спутниковые системы и технологии. Формирование активного координатно-временного пространства на основе навигационной инфраструктуры ГЛОНАСС и др. Геодезические системы наземного, морского и космического базирования для определения местоположения и навигации подвижных объектов геопространства, в том числе транспорта, военной техники, людей и животных»; 11 – «Теория и практика математической обработки результатов геодезических измерений и информационное обеспечение геодезических работ. Автоматизированные технологии создания цифровых трёхмерных моделей технологических объектов, процессов и явлений по геодезическим данным».

Степень достоверности и апробация результатов исследования. Предложения по совершенствованию методики точного дифференциального позиционирования были проверены и применены на реальных объектах: сети станций МГС и геодезической сети активных базовых станций Новосибирской области. Точность решения определялась по внешней сходимости – на основании сравнения оценок координат станций с их эталонными значениями из каталога. Таким образом, представленные результаты обладают высокой степенью достоверности.

Результаты исследований и основные положения диссертации обсуждались и были одобрены на Международных научных конгрессах «Интерэкспо ГЕОСибирь» (2011, 2012, 2013 гг., г. Новосибирск); Международном семинаре аспирантов и молодых учёных «3S-2012 GEOMIR» (27.07.2012, г. Москва), II Международной научно-технической конференции «Навигационные спутниковые системы, их роль и значение в жизни современного человека», посвящённой 30-летию запуска на орбиту первого навигационного космического аппарата ГЛОНАСС (10–14.10.2012, г. Железногорск, Красноярский край).

Основные результаты диссертационного исследования использованы во ФГУП «Сибирский государственный ордена Трудового Красного Знамени научно-исследовательский институт метрологии», в научно-производственном центре по метрологии, стандартизации и сервисному обслуживанию средств измерений ФГБОУ ВПО «СГГА», а также в учебном процессе на кафедре физической геодезии и дистанционного зондирования СГГА при изучении специальных дисциплин студентами направления «Геодезия и дистанционное зондирование» и специальности «Космическая геодезия» и при выполнении курсовых и дипломных работ.

Публикации. Основные результаты диссертационного исследования опубликованы в четырёх научных работах, из них две статьи – в изданиях, входящих в перечень рецензируемых научных журналов и изданий для опубликования основных научных результатов диссертаций.

Объём и структура диссертации. Общий объём диссертации составляет страниц. Диссертация состоит из введения, трёх разделов, заключения, списка литературы, списка сокращений и условных обозначений, четырёх приложений, содержит 13 таблиц и 23 рисунка. Список литературы включает 122 наименования.

1 МЕТОДИКА ТОЧНОГО ДИФФЕРЕНЦИАЛЬНОГО

ПОЗИЦИОНИРОВАНИЯ

1.1 Общие сведения о методике Методика точного дифференциального позиционирования (Precise Point Positioning) была представлена в статьях [94, 120] и получила развитие в работах [52, 60, 62, 69, 74, 75, 97, 99, 100, 105, 116, 117]. Её суть заключается в оценивании положения некоторой точки позиционируемого объекта в общеземной системе координат по результатам фазовых и кодовых ГНСС-измерений с использованием апостериорно уточнённой эфемеридно-временной информации (ЭВИ): эфемерид и поправок бортовых шкал времени навигационных космических аппаратов. Методика точного дифференциального позиционирования (ТДП) является одной из наиболее передовых и перспективных методик математической обработки ГНССизмерений в силу следующих её особенностей:

высокой точности определения местоположения в общеземной системе координат (ОЗСК) – среднеквадратические погрешности (СКП) оценок геодезических координат: широты и долготы в линейной мере порядка 1 см и менее, геодезической высоты – от 1 до 3 см [7, 28, 36, 53];

глобальной доступности и автономности определения координат (потребителю необходимы ГНСС-измерения, выполненные лишь на позиционируемом объекте, без непосредственной привязки к пунктам геодезических сетей);

высокой степени автоматизации обработки ГНСС-измерений.

При оценивании координат позиционируемого объекта по методике ТДП в роли исходной координатной основы выступает созвездие НКА, обеспеченное уточнённой ЭВИ. Уточнение ЭВИ выполняется путём обработки результатов ГНСС-измерений сети опорных станций. При этом сам процесс уточнения ЭВИ не входит в методику ТДП, скрыт от потребителя и выполняется специализированными центрами обработки информации, в частности центрами анализа МГС.

Методика ТДП в работах [7, 28, 36] ассоциируется с абсолютным методом определения координат. В работе [117] справедливо отмечается, что эта методика близка и к абсолютному методу, и к дифференциальному. Необходимо различать методику и алгоритм. Действительно, алгоритм методики ТДП в равной степени может быть использован и при абсолютном, и при дифференциальном позиционировании – в зависимости от того, какая ЭВИ используется: оперативная (прогнозная) или апостериорно уточнённая. Однако методика ТДП как целое предполагает использование именно апостериорно уточнённой ЭВИ, и только в этом случае достигается высокая точность оценивания координат, продемонстрированная в работах [7, 28, 36, 53]. Так как методика ассоциируется с некоторым уровнем точности определения координат, который обеспечивается только с использованием результатов предварительно выполненного сетевого решения (см.

1.2.4), методику ТДП следует относить именно к дифференциальному методу космической геодезии.

В нормативно-технической документации [13, 14, 15, 16, 41], устанавливающей термины и определения в области геодезии и ГНСС-технологий, отсутствует определение методики Precise Point Positioning (PPP). В русскоязычной литературе название методики чаще всего либо употребляют в исходном англоязычном варианте, как в работах [11, 38], либо переводят как «метод точного позиционирования», например в статьях [7, 34]. Применяемый в данной работе термин «методика точного дифференциального позиционирования» однозначно указывает на соответствующий метод космической геодезии.

В настоящее время методика ТДП представлена несколькими вариантами, описанными в статьях [52, 60, 62, 69, 74, 75, 94, 97, 99, 100, 105, 116, 117, 120].

Место методики ТДП в структуре методов космической геодезии показано на рисунке 1: методы приведены от общего к частному в направлении сверху вниз на основе классификаций, представленных в источниках [2, 43, 80].

КНП Т Т ТАП Т

Рисунок 1 – Точное дифференциальное позиционирование (ТДП) в структуре методов космической геодезии Для решения узкой задачи определения местоположения какого-либо объекта, находящегося вблизи Земли, наилучшим образом подходят геометрические методы космической геодезии, к которым относится дифференциальный метод.

Абсолютный геометрический метод является наиболее простым. Он позволяет выполнять позиционирование в режиме реального времени с использованием одного комплекта аппаратуры потребителя с погрешностями порядка нескольких метров, согласно оценке точности в работе [80 с. 437]. Абсолютные координатные определения доступны в любой точке, где наблюдается не менее четырёх НКА.

Точность оценивания координат в этом методе ограничена точностью прогноза таких параметров динамической системы, как поправки бортовых шкал времени (см. 1.2.4) и эфемериды НКА.

Относительный геометрический метод позволяет определять приращения координат между точками по выполненным в этих точках синхронным ГНССизмерениям. Высокая точность оценивания приращений координат достигается в том случае, если достаточно велика корреляция погрешностей синхронных измерений. По оценке, приведённой в работе [80 с. 432], погрешность определения вектора относительным методом в статическом режиме пропорциональна его длине с коэффициентом от 10-7 до 10-6 при длинах до 100 км. Таким образом, погрешности оценивания координат в относительном методе значительно возрастают при увеличении длин определяемых векторов.

Согласно [80 с. 416–417], в дифференциальном методе существует три способа представления поправок: в пространстве решений, в пространстве измерений и в пространстве состояний динамической системы. Примерами представления поправок в пространстве измерений могут служить сетевое дифференциальное позиционирование в реальном времени (Network Real Time Kinematic – NRTK), как указано в статье [117], а также дифференциальное позиционирование по кодовым измерениям с коррекцией по навигационному параметру (КНП), описанное в [2 с. 28–29]. В пространстве решений (которое является подпространством состояний динамической системы) поправки представляются в методике локального дифференциального позиционирования с коррекцией координат (КК) [2 с. 29–30].

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

Методика ТДП относится к дифференциальному методу космической геодезии с представлением поправок в пространстве состояний динамической системы:

поправки учтены в исходных данных – уточнённой ЭВИ. Такой способ позволяет обеспечить адекватность модели динамической системы независимо от расстояния между опорной геодезической сетью и позиционируемым объектом. Поэтому в системах глобальной и широкозонной дифференциальной коррекции поправки представляются в пространстве состояний. Методика ТДП применяется при расчёте координат потребителя, использующего такие системы дифференциальной коррекции, как GDGPS [93], StarFire [67], OmniSTAR XP [8], Trimble CenterPoint RTX [98]. Методики расчёта координат потребителя, использующие поправки в пространстве состояний динамической системы, не позволяющие достигать дециметрового уровня точности, как правило, предполагают использование результатов кодовых измерений и применяются преимущественно для целей навигации в таких системах дифференциальной коррекции ГНСС, как СДКМ, WAAS, EGNOS, MSAS, GAGAN. На рисунке 1 совокупность этих методик обозначена как «стандартное дифференциальное позиционирование» (СДП).

Таким образом, методика точного дифференциального позиционирования характеризуется уникальным сочетанием высокого уровня точности определения координат и глобальной применимости.

Методика ТДП может применяться в постобработке или в режиме реального времени. Согласно определению, приведённому в ГОСТ 15971–90 [12], режим реального времени – это «режим обработки информации, при котором обеспечивается взаимодействие системы обработки информации с внешними по отношению к ней процессами в темпе, соизмеримом со скоростью протекания этих процессов». Оценивание координат по методике ТДП в режиме реального времени может происходить с задержкой порядка нескольких секунд, в течение которых выполняется уточнение ЭВИ и передача её потребителю. При этом для доступа к источнику апостериорно уточнённой ЭВИ требуется наличие дополнительного канала связи в аппаратуре потребителя ГНСС, то есть к аппаратуре предъявляются дополнительные требования по сравнению с абсолютным методом.

Оценивание координат по методике ТДП может выполняться в различных режимах:

в статическом режиме оцениваются постоянные координаты некоторой точки, относительно которой движение принимающей ГНСС-антенны в течение сеанса измерений задаётся известным законом;

в кинематическом режиме координаты некоторой точки позиционируемого объекта оцениваются на каждый момент измерений и формируют случайный процесс;

режим «стой-иди» является комбинацией поочерёдно используемых режимов: кинематического – при движении позиционируемого объекта, статического – в периоды остановок (сеанс ГНСС-измерений при этом не должен прерываться).

По эмпирическим оценкам, представленным в работах [7, 28, 36, 53], в статическом режиме среднеквадратические погрешности оценивания координат по методике ТДП составляют от 1 до 3 см; в кинематическом режиме – порядка дециметра, согласно [53].

Методика ТДП может применяться в указанных режимах для определения местоположения различных объектов: транспортных средств, военной техники, геофизической и другой аппаратуры, конструктивных частей зданий и сооружений, геодезических пунктов и так далее. В силу того, что совместно с координатами позиционируемого объекта в методике ТДП определяется дополнительные параметры, возможно применение этой методики для точной синхронизации шкал времени [66] и получения метеорологической информации [34].

Существует несколько программных продуктов, позволяющих уточнять координаты позиционируемого объекта по методике ТДП. Эта возможность реализована как в некоторых научных программных продуктах, в частности Bernese [56], так и в коммерческом программном обеспечении, к которому относится GRAFNAV/GRAFNET (Waypoint GPS) [7]. Существуют свободно распространяемое программное обеспечение, например GPS Toolkit [115], RTKLIB [107], BNC [59], а также бесплатные или условно бесплатные онлайн-сервисы уточнения координат основанные на методике PPP: ВМ СДКМ (Россия) [39], APPS (США) [113], CSRS-PPP (Канада) [104], GAPS (Канада) [73], Magic GNSS (Испания) [101]. Аналогичный по функциям сервис в России также представлен Информационно-аналитическим центром координатно-временного и навигационного обеспечения (ИАЦ КВНО) [21], однако, в данном случае методика ТДП не используется*.

Информация предоставлена ИАЦ КВНО по запросу.

1.2 Математическая модель динамической системы и алгоритм методики точного дифференциального позиционирования 1.2.1 Измеряемая геометрическая дальность Рассматриваемая динамическая система включает Землю, расположенную на ней или вблизи неё совокупность измерительных средств, созвездие навигационных космических аппаратов, а так же все источники возмущающих воздействий на их движение и результаты измерений в системе. Вероятно, наиболее полная и современная модель этой динамической системы представлена в Соглашениях Международной службы вращения Земли и систем отсчёта (IERS Conventions 2010) [82]. Также детальные модели ГНСС-измерений в динамической системе представлены в работах [1, 2, 10, 80, 118]. Некоторые специфические вопросы моделирования динамической системы, связанные с применяемыми продуктами МГС, рассмотрены в работе [95], а также в материалах на сайте координационного центра МГС [83].

Основой модели динамической системы является фундаментальное уравнение космической геодезии [27 с. 10]. Представим фундаментальное уравнение в инерциальной (небесной) системе координат где r – геоцентрический вектор положения фазового центра антенны НКА в момент излучения радионавигационного сигнала t в небесной системе координат;

R – геоцентрический вектор положения фазового центра антенны позиционируемого объекта в момент приёма сигнала t в небесной системе координат;

– топоцентрический вектор НКА в небесной системе координат.

Моменты излучения t и приёма t радионавигационного сигнала определяются по единой шкале времени (ЕШВ).

В инерциальных системах координат длина топоцентрического вектора НКА (геометрическая дальность) в первом приближении соответствует пути распространения радиосигнала, то есть измеряемой дальности. Во вращающихся системах координат, в частности в ОЗСК, это соответствие нарушается. Исходные данные о положении позиционируемого объекта и навигационных космических аппаратов, а также оценка местоположения позиционируемого объекта (результат) в рассматриваемой задаче представляются в ОЗСК. Связь ОЗСК и небесной геоцентрической системы координат задаётся математической моделью вращения Земли. Согласно [6, 31], в классическом представлении математическая модель вращения Земли определяется произведением девяти элементарных матриц поворота, учитывающих прецессию и нутацию, собственно вращение Земли, неравномерность её вращения и полярное движение мгновенной оси в теле Земли, в соответствии с формулой где – матрица направляющих косинусов осей ОЗСК относительно осей небесной системы координат;

A, A, zA – экваториальные параметры прецессии на интервале [T0, t] – между стандартной эпохой J2000.0 и текущей эпохой t;

A – средний наклон эклиптики к экватору в эпоху t;

, – компоненты астрономической нутации в долготе и наклоне в эпоху t;

S – истинное звёздное время на меридиане Гринвича в момент времени t;

xp, yp – параметры полярного движения оси вращения в теле Земли в эпоху t;

для правой системы координат:

Вектор в небесной системе координат связан с векторами положения фазовых центров антенн НКА и позиционируемого объекта в ОЗСК формулой где r – вектор положения фазового центра антенны НКА в момент излучения радионавигационного сигнала t в ОЗСК;

R – вектор положения фазового центра антенны позиционируемого объекта в момент приёма сигнала t в ОЗСК.

Для удобства решения задачи позиционирования целесообразно ввести мгновенную инерциальную систему координат, совпадающую с ОЗСК в момент достижения сигналом антенны потребителя. Умножим слева левую и правую части уравнения (4) на матрицу W t :

Для мгновенной инерциальной системы координат (МИСК), совпадающей с ОЗСК в момент достижения радионавигационным сигналом фазового центра антенны потребителя, введём индекс «м»:

Введём обозначение для матрицы поворота осей ОЗСК относительно инерциального пространства за время распространения радиосигнала между антенной НКА и антенной потребителя:

Связь векторов положения фазовых центров антенн НКА и потребителя в ОЗСК и топоцентрического вектора в МИСК выражается формулой:

Уравнение (6) позволяет связать измеряемую дальность м c вектором положения позиционируемого объекта в ОЗСК:

В методике ТДП используются результаты фазовых и кодовых ГНССизмерений на двух частотах. Положение фазовых центров антенн НКА и позиционируемого объекта различно для разных частот. ГНСС-измерениям на частоте f i соответствует геометрическая дальность мi, моделируемая следующим образом, согласно [106]:

где is – отклонение фазового центра антенны НКА от среднего положения вдоль линии позиционируемый объект – НКА для частоты с номером i, обусловленное диаграммой направленности антенны, м;

Ai – отклонение фазового центра принимающей антенны позиционируемого объекта для частоты с номером i вдоль линии наблюдаемый НКА – позиционируемый объект, м.

Отклонения фазовых центров (phase center variations) определены путём калибровки. Величины отклонений фазовых центров вместе с параметрами средних фазовых центров (phase center offsets) доступны в файлах ANTEX [106], размещённых на сервере МГС*. В настоящее время применяется модель ГНСС-антенн IGS08. Расчёт положения средних фазовых центров будет рассмотрен в 1.2.2, 1.2.3.

Измеряемая геометрическая дальность между эффективными фазовыми центрами антенны НКА в момент излучения сигнала и антенны позиционируемого объекта в момент приёма сигнала, соответствующая «ионосферно-свободным»

комбинациям результатов ГНСС-измерений, вычисляется по формуле где f1, f 2 – несущие частоты радиосигналов ГНСС.

1.2.2 Модель движения навигационного космического аппарата Апостериорно уточнённая эфемеридно-временная информация ГЛОНАСС и GPS, используемая в методике ТДП, свободно предоставляется Международной URL: ftp://igscb.jpl.nasa.gov/igscb/station/general/ (дата обращения: 17.11.2013).

ГНСС-службой (МГС) [84], а также её отдельными центрами, в частности российским Информационно-аналитическим центром координатно-временного и навигационного обеспечения (ИАЦ КВНО) [20]. Также уточнённая ЭВИ ГНСС предоставляется российской Системой высокоточного определения эфемерид и временных поправок (СВОЭВП) [42].

Апостериорно уточнённая ЭВИ включает так называемые «точные» эфемериды в виде таблиц координат НКА в ОЗСК в моменты, заданные в единой шкале времени (ЕШВ) с интервалом 15 мин. В зависимости от конкретного источника информации, эфемериды могут относиться к различным точкам НКА. Эфемериды СВОЭВП содержат координаты фазового центра антенны НКА, а эфемериды МГС – координаты центра масс НКА. Также эфемериды представляются в различных реализациях общеземной системы координат, в частности эфемериды СВОЭВП представляются в ПЗ-90.02 (переход к ПЗ-90.11 запланирован на год*). В продуктах МГС реализацией ОЗСК является координатная основа IGb [83], в качестве ЕШВ применяется шкала времени МГС (IGS Time Scale), связанная со шкалой времени GPS [110].

Рассмотрим модель движения НКА в продуктах МГС. Положение среднего фазового центра НКА для частоты c номером i в момент излучения радиосигнала в ОЗСК определяется формулой:

где rц.м. – вектор положения центра масс НКА в ОЗСК в момент излучения сигнала антенной;

r ф.ц.i – вектор положения среднего фазового центра НКА относительно его центра масс в ОЗСК в момент излучения сигнала для частоты с номером i.

Аппаратура потребителя ГНСС записывает измерения, выполненные в моменты, определяемые по собственной шкале времени – шкале времени потребиИнформация предоставлена по запросу службой поддержки СВОЭВП.

теля (ШВП). Момент излучения навигационного сигнала по ЕШВ определяется по формуле где t – момент излучения сигнала по единой шкале времени, с;

Pi – результат кодового измерения дальности на частоте с номером i (кодовая псевдодальность), м;

tA – номинальный момент записи сигнала по ШВП, с;

dt s – сдвиг БШВ НКА относительно ЕШВ в момент излучения сигнала, с.

Вектор положения центра масс НКА rц.м. определяется путём интерполяции координат, представленных в табличном виде в ОЗСК с интервалом 15 мин.

Наиболее часто применяются интерполяционные формулы Лагранжа [80 с. 52–53;

118 c. 34–35], а также Ньютона, Чебышёва [70 с. 38]. Согласно исследованию [111], оптимальный порядок интерполяционного многочлена равен 10.

Положение фазовых центров передающей антенны НКА задаётся в системе координат навигационного космического аппарата (СК НКА), изображённой на рисунке 2: начало находится в центре масс НКА; ось z' направлена к центру масс Земли (геоцентру); направление оси y' (оси вращения солнечных панелей) определяется векторным произведением орта оси z' и вектора НКА – Солнце; орт оси x' дополняет правую тройку векторов [95, 106]:

где i, j, k – орты осей x', y',z' системы координат НКА, представленные в ОЗСК;

r НКА-С – орт вектора НКА – Солнце в ОЗСК.

Рисунок 2 – Система координат навигационного космического аппарата При расчёте положения фазового центра НКА орт направления НКА – Солнце может быть заменён ортом вектора положения Солнца в ОЗСК. Для расчёта орта геоцентрического вектора Солнца могут быть применены эфемериды DE\LE421 [91] или EPM2011 [35].

Оси СК НКА жёстко связаны с телом НКА, описанное положение осей – номинальная ориентировка. В описанной системе координат параметры фазовых центров антенн НКА публикуются МГС. В других источниках могут использоваться другие системы координат. Так, согласно [96], направления осей могут отличаться от описанных для НКА некоторых моделей, в частности для НКА GPS блока IIR.

Для сохранения номинальной ориентировки НКА поддерживается его вращение вокруг оси с приблизительно постоянной угловой скоростью кроме тех случаев, когда НКА находится в тени, возвращается к номинальной ориентировке непосредственно после выхода из тени, совершает полуденный или полуночный манёвр. Неучтённые движения НКА в эти периоды могут вносить дециметровые погрешности в результаты измерений, поэтому следует либо исключать эти измерения из обработки, либо применять модели ориентации, учитывающие особенности систем управления различных модификаций НКА [96, 114].

Вектор положения среднего фазового центра антенны НКА относительно его центра масс в ОЗСК для каждой из частот определяется выражением:

где xф.ц.i, yф.ц.i, zф.ц.i – координаты среднего фазового центра передающей антенны НКА для частоты fi в СК НКА, представленной на рисунке 2.

Координаты фазовых центров антенн ГНСС, в том числе антенн НКА, определяются МГС и приводятся в файле параметров антенн в формате ANTEX (atx) [106].

В настоящее время точность определения положения НКА с использованием окончательных эфемерид МГС характеризуются среднеквадратическими погрешностями по одной координате 2,5 см для НКА GPS, 3 см для ГЛОНАСС (приведены средние значения из СКП трёх координат, оценка точности выполнялась по сходимости с результатами лазерной локации); среднеквадратические погрешности прогнозных эфемерид МГС (ultra-rapid predicted half) составляют порядка 5 см по каждой координате для НКА GPS [84]. То есть в результате апостериорного уточнения погрешности эфемерид НКА GPS уменьшаются примерно в 2 раза.

1.2.3 Модель движения позиционируемого объекта Требуется определить местоположение некоторой точки позиционируемого объекта, на котором установлена аппаратура потребителя ГНСС. Связь векторов положения этой точки и среднего фазового центра антенны потребителя в ОЗСК в момент приёма навигационного сигнала может быть представлена формулой:

где Rц0 – вектор приближённого положения позиционируемой точки в ОЗСК;

R – оцениваемый вектор поправок координат позиционируемой точки в ОЗСК;

Rц t – априорная модель движения позиционируемой точки относительно её несмещённого положения в ОЗСК (другие слагаемые в формуле (14) также могут зависеть от времени);

R т.о. – вектор в ОЗСК от позиционируемой точки до точки относимости антенны (точки, относительно которой при калибровке определялись параметры фазового центра антенны [106]);

R ф.ц. i – вектор в ОЗСК от точки относимости антенны до среднего фазового центра антенны потребителя для частоты fi.

С учётом выражений (8), (9), (14) связь точного и приближённого значений геометрической дальности между фазовыми центрами антенн НКА и позиционируемого объекта определяется формулой:

где м 0 – приближённая геометрическая дальность.

Рассмотрим формулу (14) для случая, когда позиционируемым объектом является геодезический пункт. В этом случае чаще всего требуется определить несмещённое положение центра пункта в эпоху ГНСС-измерений Вектор приближённого положения центра Rц 0 может задаваться по результатам навигационного решения [2 c. 9–15] с погрешностью порядка 10 м.

Априорная модель движения пункта в течение сеанса измерений относительно его несмещённого положения определяемого положения Rц в ОЗСК:

где R т t – вектор смещения вследствие прилива в твёрдом теле Земли, вызванный внешним приливообразующим потенциалом (до нескольких дециметров), Rок t – вектор смещения вследствие деформации земной коры из-за океанической приливной нагрузки (до дециметра в прибрежной зоне);

Rатм t – вектор смещения вследствие деформации земной коры из-за колебаний атмосферного давления в течение суток (порядка сантиметра);

Rпол t – вектор смещения в результате деформации земной коры, вызванной изменением положения оси вращения в теле Земли (до 2 см по высоте, несколько миллиметров в плане);

Rпол.ок. t – вектор смещения пункта вследствие перераспределения океанической нагрузки на поверхность твёрдой Земли из-за изменения положения оси вращения в её теле (до 2 мм по высоте, менее 1 мм в плане).

Выражение (17) включает поправки, которые рекомендованы к учёту МСВЗ.

Информационные материалы об этих поправках (в том числе рабочие формулы, таблицы параметров, ссылки на теоретические труды, ссылки на исходные коды программных модулей для вычисления указанных поправок) приведены в Соглашениях МСВЗ 2010 [82 с. 99–122]. Однако не все эффекты, приводящие к смещению центра пункта и антенны, могут быть адекватно отражены представленной моделью. Например, чаще всего неизвестным остаётся характер деформации зданий и сооружений, на которых зачастую располагаются центры пунктов и устанавливаются антенны.

Положение точки относимости антенны (antenna reference point) определяется по результатам промеров относительно центра пункта в топоцентрической горизонтальной системе координат. Для перехода от горизонтальной системы координат пункта к ОЗСК требуется выполнить следующее преобразование:

где – матрица поворота осей системы координат при переходе от горизонгор тальной СК к ОЗСК;

Eт.о., N т.о., U т.о. – координаты точки относимости антенны в горизонтальной СК относительно центра пункта в направлениях: на восток, на север, в зенит – соответственно;

, – астрономические координаты пункта: широта и долгота – могут быть в данном случае заменены на геодезические B, L, которые в соответствии с [ с.190–191] могут быть вычислены по формулам:

где Bj – значение геодезической широты пункта на j-ой итерации;

L – геодезическая долгота пункта;

N j – значение радиуса кривизны эллипсоида в первом вертикале пункта на jой итерации;

a – большая полуось общеземного эллипсоида;

e – эксцентриситет общеземного эллипсоида.

В ОЗСК вектор от точки относимости ГНСС-антенны потребителя до её фазового центра для частоты fi определяется формулой:

где Eф.ц.i, Nф.ц.i,Uф.ц.i – координаты вектора от точки относимости до среднего фазового центра принимающей ГНСС-антенны для частоты c номером i в горизонтальной системе координат (приводятся в файлах ANTEX для различных моделей антенн [106]).

1.2.4 Измерительная информация В рассматриваемом варианте методики точного дифференциального позиционирования измерительная информация представлена «ионосферносвободными»* комбинациями результатов двухчастотных кодовых и фазовых измерений дальности между фазовым центром антенны каждого НКА и фазовым центром антенны потребителя. «Ионосферно-свободная» комбинация кодовых измерений дальности на двух частотах вычисляется по формуле [1 с. 315] В данной работе термин «ионосферно-свободная» употребляется в кавычках, так как эта комбинация содержит остаточное влияние ионосферы. Во избежание путаницы под «ионосферно-свободной» в работе понимается именно двухчастотная комбинация, хотя для НКА, передающих сигнал на трёх частотах, возможно формирование трёхчастотной «ионосферносвободной» комбинации ГНСС-измерений [82 с. 146].

где P, P2 – результаты измерения дальности по дальномерному коду на частотах f1, f 2, соответственно, м.

Наиболее высокая точность кодовых измерений достигается в том случае, когда на обеих частотах доступны измерения по дальномерному коду высокой точности (P-код). Некоторые модели геодезических ГНСС-приёмников не способны выполнять P-кодовые измерения на частоте L1. В этом случае для формирования «ионосферно-свободной» комбинации используются результаты измерений на частоте L1 по коду стандартной точности (С\A-код, обозначаемый в файлах RINEX как «C1»).

«Ионосферно-свободная» комбинация результатов фазовых измерений дальности в метрах [118 с. 97]:

где 1, 2 – длины волн радиосигналов на первой и второй частотах, м;

1, 2 – результаты фазовых ГНСС-измерений на первой и второй частотах, цикл.

В литературе, посвящённой ГНСС, измеренные дальности принято называть «псевдодальностями» (кодовыми и фазовыми) [1], соответствующие понятия закреплены в ГОСТ Р 53607–2009 [16]. Результат измерения фазы в отечественной литературе называют «псевдофазой» [10].

Необходим контроль качества исходной измерительной информации. Целесообразно предварительно выполнять отбраковку грубых результатов измерений, выявлять и по возможности восстанавливать потери счёта фазовых циклов. Алгоритмы выявления и восстановления потерь счёта циклов рассмотрены в работах [54, 65, 68, 119].

Детально модели кодовых и фазовых псевдодальностей рассмотрены в работах [1, 80, 118], на основании которых можно представить обобщённую модель «ионосферно-свободной» комбинации кодовых измерений дальности в методике ТДП где T – влияние тропосферной задержки сигнала на результат измерения, м;

c – скорость света в вакууме, м/с;

dt s – сдвиг шкалы времени НКА, с;

dtA – сдвиг шкалы времени потребителя, с;

IP – остаточное влияние ионосферы на комбинацию результатов кодовых измерений дальности, м;

dg – гравитационная задержка распространения радиосигнала, с;

vP – случайная погрешность функциональной модели кодовых измерений, м;

Модель «ионосферно-свободной» комбинации фазовых измерений:

где I – остаточное влияние ионосферы на комбинацию результатов фазовых измерений дальности, м;

BA – псевдонеоднозначность* фазовых измерений (параметр, включающий фазовую неоднозначность, задержки в аппаратуре и начальную разность фаз генераторов частоты НКА и НИП), м;

– длина волны суммарной комбинации фазовых измерений, м;

Термин «псевдонеоднозначность» введён автором вместо слова «смещение» (bias), применяемого в зарубежной литературе для обозначения указанного параметра.

w – набег фазы вследствие изменения взаимной ориентировки антенн, цикл;

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

Длина волны суммарной комбинации фазовых ГНСС-измерений определяется формулой:

Модели (26), (27) не учитывают многопутность радиосигналов. В рассматриваемом исходном варианте методики ТДП также не учитываются задержки кодовых сигналов в аппаратуре ГНСС. При оценивании параметров по методике ТДП задержки сигналов в аппаратуре неотделимы от параметров шкал времени. По соглашению МГС, задержка кодовых сигналов в аппаратуре НКА для «ионосферносвободной» комбинации высокоточных кодовых псевдодальностей (P1, P2) равна нулю [95 с. 18]. То есть в случае применения поправок БШВ GPS, предоставленных МГС, и использования комбинации высокоточных кодовых измерений P1, P при оценивании координат учёт задержек в аппаратуре не требуется. Однако при использовании комбинации других типов кодовых псевдодальностей, например стандартной точности C1 и высокой точности P2, а также в одночастотном варианте методике ТДП необходим учёт разности задержек кодовых сигналов (differential code bias – DCB), согласно [95 с. 18].

Рассмотрим некоторые детерминированные составляющие моделей (26),(27), значения которых могут быть вычислены по известным формулам.

Сдвиг бортовой шкалы времени НКА dt s, согласно [95 c. 18], включает в себя релятивистскую составляющую и согласуется с определением погрешности, представленной в интерфейсном контрольном документе GPS [89 с. 93]:

где dt s – сдвиг БШВ НКА относительно ЕШВ, с;

dtрел – релятивистское отклонение БШВ, с;

dtт – поправка БШВ по данным МГС, с;

e – эксцентриситет орбиты НКА;

3,986004418 1014 м3 c2 ;

a – большая полуось орбиты НКА, м;

E – эксцентрическая аномалия НКА в момент излучения сигнала, рад;

с – скорость света в вакууме, м/с.

Поправки БШВ НКА dtт представляются в табличном виде в продуктах МГС и отдельных центров анализа данных ГНСС, а также СВОЭВП, с интервалами мин, 5 мин, 30 с, 5 с [42, 77, 84]. Значение dtт на момент излучения радиосигнала НКА определяется путём многочленной интерполяции. При использовании таблиц поправок БШВ с интервалом 30 с и менее, согласно [77], предпочтительна линейная интерполяция.

При использовании поправок БШВ, представленных отдельными центрами анализа, а также СВОЭВП, необходимо уточнять, по отношению к какой эталонной шкале времени определены поправки, и при необходимости приводить их к единой шкале времени. В продуктах МГС в качестве единой шкалы времени применяется шкала времени МГС (IGS Time Scale), связанная со шкалой времени GPS [110].

В настоящее время низкая точность прогнозирования поправки dtт является главным фактором, ограничивающим точность абсолютного метода и определяющим необходимость апостериорного уточнения ЭВИ, применяемой в методике ТДП. Среднеквадратические погрешности прогнозных значений поправок dtт по данным МГС (ultra-rapid predicted) для БШВ НКА GPS составляют около 3 нс, что соответствует погрешности в измеряемой дальности около 0,9 м и сопоставимо с погрешностями модели БШВ, передаваемой в навигационном сообщении. Апостериорное уточнение по данным ГНСС-измерений сети станций МГС позволяет сократить среднеквадратические погрешности модели БШВ до 0,15 нс (0,045 м) в «сверхбыстрых» предварительных продуктах МГС (ultra-rapid observed) и до 0, нс (0,02 м) в окончательных (final).

Тропосферная задержка – суммарная задержка сигнала, возникающая при прохождении радиосигналом нейтральной атмосферы. Модель тропосферной задержки в соответствии с Соглашениями МСВЗ 2010 [82 с. 135–137]:

где Thz – гидростатическая зенитная тропосферная задержка, м;

Twz – влажная зенитная тропосферная задержка, м;

GN, GE – горизонтальные градиенты тропосферной задержки, м;

m – функции отображения, зависящие от угловой высоты h НКА в горизонтальной топоцентрической системе координат позиционируемого объекта;

A – азимут НКА в горизонтальной системе координат позиционируемого объекта, рад.

Наибольший вклад в тропосферную задержку (от 80 % до 100 %) вносит гидростатическая составляющая. Гидростатическая зенитная задержка с высокой точностью может быть вычислена по формуле Саастамойнена где p0 – атмосферное давление вблизи фазового центра ГНСС-антенны потребителя, гПа;

B – геодезическая широта фазового центра ГНСС-антенны потребителя, рад;

H – геодезическая высота фазового центра ГНСС-антенны потребителя, м.

Общий вид функций отображения для гидростатической и влажной составляющих тропосферной задержки:

где h – угловая высота НКА над горизонтом;

a, b, c – коэффициенты модели.

Коэффициенты могут быть получены несколькими способами. Для наиболее точного учёта параметров тропосферы МСВЗ рекомендует применять функцию отображения VMF1 (Vienna Mapping Function 1). В этой модели коэффициенты b, c вычисляются по эмпирическим формулам, коэффициент a публикуется с периодичностью 6 ч для всех станций МГС, а также в узлах глобальной сетки 2,52, [82 с. 136]. Также может быть использована менее точная функция отображения GMF (Global Mapping Function), параметрами которой являются приближённые геодезические координаты станции и номер дня в году.

МСВЗ рекомендует оценивать влажную зенитную задержку Twz совместно с искомыми параметрами, так как в настоящее время нет точных способов априорного расчёта этой составляющей [82 с. 135].

Также рекомендуется учитывать горизонтальный градиент тропосферной задержки. Величины этих параметров порядка 1 мм. При малых угловых высотах НКА влияние горизонтального градиента тропосферной задержки на ГНССизмерения может достигать нескольких сантиметров. Часть горизонтального тропосферного градиента в плоскости меридиана обусловлена несферичностью атмосферы и для средних широт составляет по направлению к экватору около минус 0,5 мм в северном полушарии и плюс 0,5 мм в южном. Функция отображения горизонтального градиента тропосферной задержки задаётся формулой Горизонтальный градиент может оцениваться совместно с другими параметрами по методике ТДП или может быть вычислен с помощью модели, разработанной Европейским центром среднесрочного прогнозирования погоды.

Величина полной тропосферной задержки может достигать 2,3 м в зените и значительно возрастает вблизи горизонта [82 с. 135].

Ионосферная задержка – эффект возникающий при распространении радиосигнала в диспергирующей среде. В соответствии с Соглашениями МСВЗ [82], ионосферная задержка аппроксимируется разложением в ряд Тейлора до величин порядка f i 4. Вклад наиболее значимого первого члена разложения, составляющий около 99,9 % и достигающий нескольких десятков метров, устраняется путём образования «ионосферно-свободных» комбинаций измерений на двух частотах по формулам (24), (25). Остаточное влияние ионосферы на «ионосферносвободные» комбинации результатов ГНСС-измерений, рассмотренное в работах [72, 78, 81, 82], может достигать порядка 3 см (100 пс). В методике ТДП следует учитывать влияние второго и третьего членов в соответствии с [82 с. 144] по формулам:

где I P – остаточное влияние ионосферы на «ионосферно-свободную» комбинацию кодовых псевдодальностей, м;

s2, s3 – второй и третий коэффициенты разложения ионосферной задержки;

f1, f 2 – несущие частоты радиосигналов, Гц.

I – остаточное влияние ионосферы на «ионосферно-свободную» комбинацию фазовых псевдодальностей, м;

При учёте остаточного влияния ионосфера представляется в виде тонкого слоя, расположенного на высоте 450 км [82 с. 144] (в разных моделях высота слоя может отличаться). Точка, в которой сигнал пересекает этот слой, называется «ионосферной точкой прокола» (ionospheric pierce point) [78] или «ионосферной точкой» [1]. Согласно [72, 82], коэффициенты s2, s3 для ГНСС могут быть вычислены по следующим приближённым формулам*:

где Bp – вектор магнитной индукции в ионосферной точке прокола, Тл;

– угол между направлением вектора магнитной индукции и направлением распространения радиосигнала, рад;

S – полное содержание свободных электронов вдоль пути распространения сигнала, аппроксимированного прямой линией НКА – НИП (STEC – Slant Electron Content), м –2;

h – угловая высота НКА, рад.

Для вычисления вектора магнитной индукции в точке прокола может использоваться планетарная модель магнитного поля Земли IGRF (International Geomagnetic Reference Field) [78, 85].

В соответствии с IERS Conventions 2010 [82 c. 140], должна также учитываться дополнительная поправка за искривление пути радиосигнала в ионосфере s3, однако приведённый в этом документе способ учёта данного эффекта противоречит описанному в статье [81], на которую ссылается документ. Согласно [81], поправка должна учитываться с одинаковыми знаками и коэффициентами в результатах фазовых и кодовых измерений.

Полное содержание электронов по направлению НКА – позиционируемый объект для учёта остаточного влияния ионосферы по формулам (35)–(38) может быть приближённо вычислено с применением ионосферной комбинации кодовых и фазовых измерений, свободной от геометрии [82 с. 145]:

где I – ионосферная комбинация фазовых псевдодальностей, м;

1, 2 – длины волн радиосигналов на первой и второй частотах, м;

1, 2 – результаты фазовых ГНСС-измерений на первой и второй частотах, цикл;

PI – ионосферная комбинация кодовых псевдодальностей, м;

I PI – разность ионосферных комбинаций фазовых и кодовых псевдодальностей, осреднённая на непрерывно наблюдаемом участке дуги траектории НКА (без потерь счёта циклов фазы), м;

d DCB – разность задержек кодовых сигналов в аппаратуре НКА на двух частотах (differential code bias – DCB), м;

dA DCB – разность задержек кодовых сигналов в аппаратуре потребителя на двух частотах, м.

Погрешности моделей остаточного влияния ионосферы на «ионосферносвободные» комбинации результатов измерений (35), (36) могут достигать величин более сантиметра при угловой высоте НКА около 10 и высокой солнечной активности, согласно [81].

Гравитационная задержка – эффект, возникающий при распространении сигнала между точками с разным гравитационным потенциалом. Гравитационная задержка может быть вычислена по формуле Формула (42) эквивалентна формуле представленной в [82 с. 164]. В ГНСС гравитационная задержка может достигать 70 пс (2 см), погрешность учёта по формуле (42) порядка 3 пс (1 мм) [82 c. 164].

Набег фазы – эффект, возникающий в процессе изменения взаимной ориентировки антенн НКА и позиционируемого объекта. Эффект рассмотрен в работе [57]. Величина набега фазы может быть вычислена по формулам:

где w – набег фазы, цикл;

D, D – эффективные дипольные векторы антенн НКА и НИП, соответственно;

iм, jм – орты осей x', y' системы координат НКА (рисунок 2) в МИСК;

0 – орт топоцентрического вектора НКА в МИСК;

eм, nм – орты осей горизонтальной топоцентрической системы координат НИП, направленных на восток и на север*, в МИСК.

Набег фазы оказывает влияние только на фазовые измерения. Требуется непрерывно учитывать изменение набега фазы в течение сеанса, поэтому при переходе величины w через ноль необходимо добавлять или отнимать единицу. Для вычисления набега фазы требуется знать ориентировку осей спутниковой системы координат, которая может быстро меняться в периоды, когда НКА находится в тени или выполняет манёвры. При поддержании номинальной ориентировки НКА (рисунок 2) величина набега фазы может достигать половины длины волны несущей [95 c. 11]. Теоретически рассчитанная погрешность учёта эффекта по формуле (43) составляет 0,1 мм [57].

1.2.5 Система уравнений измерений В методике ТДП модели кодовых (26) и фазовых (27) псевдодальностей представляется следующими уравнениями:

где оцениваемыми параметрами являются:

RA – вектор оцениваемых поправок координат позиционируемой точки в ОЗСК (мгновенный или средний – см. 1.2.6), м;

Tzw A – влажная зенитная тропосферная задержка, м;

В данной работе орты векторов топоцентрической системы координат, используемые при расчёте набега фазы представляют правую систему, в соответствии со статьёй [57], а также алгоритмами расчёта набега фазы в программном обеспечении GPS Toolkit и RTKLIB. В комментарии к формуле расчёта набега фазы в работах [94, 95] ошибочно указана левая тройка векторов (на север, на восток, в зенит).

dtA – сдвиг ШВП относительно ЕШВ в момент измерения, с;

BA – псевдонеоднозначность, м;

коэффициенты при оцениваемых параметрах:

0 – орт топоцентрического вектора НКА в мгновенной инерциальной системе координат, совпадающей с ОЗСК в момент достижения радионавигационным сигналом антенны потребителя;

mw – функция отображения влажной тропосферной задержки;

c – скорость света в вакууме, м/с;

правая часть уравнений:

fP, f – свободные члены уравнений кодовых и фазовых измерений, м;

vP, v – случайные погрешности функциональных моделей измерений, м.

Свободные члены вычисляются на основании выражений (26), (27), (31) по следующим формулам:

Свободные члены и коэффициенты при определяемых параметрах уравнений (47), (48) вычисляются в соответствии с моделями, описанными в пунктах 1.2.1 – 1.2.4.

Уравнения измерений (47), (48) формируют систему, параметры которой оцениваются рекуррентно на каждую эпоху измерений с использованием математического аппарата фильтра Калмана. Связь состояния динамической системы на шаге k, соответствующем моменту измерений tk, с измерительной информацией и предыдущим состоянием описывается функциональной моделью:

где xk – вектор оцениваемых параметров системы уравнений (47), (48) в эпоху k;

v k – вектор невязок уравнений измерений на шаге k, характеризуемый ковариационной матрицей Rk ;

Hk – матрица коэффициентов при оцениваемых параметрах уравнений измерений (47), (48) на шаге k;

f k – вектор свободных членов уравнений измерений (47), (48) на шаге k;

wk – вектор случайного изменения состояния системы на шаге k (шум случайного процесса), характеризуемый ковариационной матрицей Qk ;

Fk – матрица эволюции системы, отражающая функциональную связь между векторами состояния системы в последовательные эпохи.

Вектор оцениваемых параметров состояния динамической системы в момент записи измерений tk в рассматриваемом варианте методики ТДП включает следующие параметры:

где XA, YA, ZA – поправки координат позиционируемой точки в ОЗСК (элементы вектора RA );

n – количество НКА, наблюдаемых c объекта А в момент записи измерений tk.

1.2.6 Стохастическая модель динамической системы Стохастическая модель динамической системы характеризует её случайные составляющие, к которым относятся случайные погрешности уравнений измерений и оцениваемые параметры. Стохастическая модель представляется в виде ковариационной матрицы случайных величин.

При рассмотрении погрешности уравнения измерений, принято выделять погрешность модели и погрешность результата измерения, однако, такое деление является достаточно условным. В данной работе погрешность результата измерения отождествляется с погрешностью его функциональной модели. Результат измерения понимается как факт, который сам по себе не содержит погрешностей, а лишь является следствием различных физических эффектов, в том числе связанных с аппаратурой, методикой измерений и так далее. Погрешность возникает вследствие того, что функциональная модель (уравнение) не может идеально представить все эффекты, сказывающиеся на результатах измерений.

В [95 с. 14] указано, что формирование информационных продуктов МГС и их дальнейшее использование, в частности в методике ТДП, должны опираться на единую непротиворечивую модель динамической системы. Иначе результаты оценивания по методике ТДП будут искажены. Стохастические модели погрешностей измерений, используемые при уточнении ЭВИ и при уточнении координат потребителя по методике ТДП также должны быть совместимы. Однако для одних и тех же результатов измерений погрешности могут различаться в зависимости от их функциональной модели и решаемой задачи. Например, при уточнении эфемеридно-временной информации поправки бортовых шкал времени НКА являются определяемыми. В методике ТДП эти параметры считаются известными и учитываются в свободных членах с некоторой погрешностью, которая сказывается на общей погрешности модели. Таким образом, способы задания априорных ковариационных матриц погрешностей измерений, описанные в документах центров анализа МГС, не соответствуют функциональной модели результатов измерений в методике ТДП.

Кроме того, каждый центр анализа МГС применяет свою модель, а продукты МГС формируются как комбинация решений отдельных центров. Поэтому нет единого способа задания априорной ковариационной матрицы погрешностей измерений, строго совместимого с моделью МГС.

Были проанализированы априорные стохастические модели погрешностей ГНСС-измерений, применяющиеся центрами анализа МГС при уточнении эфемерид и поправок БШВ НКА, описанные в документах этих центров [86]. Те же модели применялись в методике ТДП, согласно источникам [53, 62, 74, 109]. Рассмотренные априорные стохастические модели погрешностей ГНСС-измерений представляются диагональными ковариационными матрицами, в которых элементы главной диагонали, соответствующие дисперсиям погрешностей, в общем виде задаются выражениями:

где D – априорная дисперсия погрешности измерений, м2;

n – номинальное стандартное отклонение погрешности, соответствующее данному типу ГНСС-измерений, м;

f h – функция от угловой высоты НКА.

Априорные стохастические модели погрешностей ГНСС-измерений, применяемые центрами анализа МГС, представлены в приложении А. Чаще всего отношение номинальных дисперсий результатов кодовых и фазовых измерений P принимается равным 100. В зависимости от применяемой модели множитель, выражающий зависимость стандартного отклонения от угловой высоты НКА, может определяться формулами:

Так как проанализированные априорные стохастические модели погрешностей ГНСС-измерений, применяемые центрами анализа МГС, изначально не предназначались для использования в методике ТДП, они не учитывают некоторых значимых факторов. В частности, не учитываются:

зависимость разброса погрешностей ЭВИ от типа применяемых продуктов МГС;

зависимость уровня шума кодовой «ионосферно-свободной» комбинации измерений от типа используемых кодовых измерений.

Далее, в 2.2.2, предложена новая априорная стохастическая модель погрешностей результатов измерений, учитывающая данные факторы.

При рекуррентном оценивании параметры x представляются случайными процессами. Каждый элемент случайного процесса – случайная величина, представляющая параметр в момент tk. Характеристики случайных процессов задаются с помощью матриц эволюции Fk и ковариационной матрицы шумов процессов Qk. Основными применяемыми в данном случае моделями случайных процессов являются: постоянный, кусочно-постоянный процесс, белый шум [4 с. 137], случайное блуждание [44 с. 150]. В зависимости от специфики временные ряды параметров (53) системы уравнений (51) моделируются различными типами случайных процессов.

Поправки координат позиционируемой точки могут моделироваться поразному в зависимости от применяемого режима: статического, кинематического, режима «стой-иди». В статическом режиме поправки координат представляются постоянными величинами, то есть происходит уточнение оценок по мере обработки дополнительных порций измерительной информации. В статическом режиме достигается наивысшая точность оценивания координат. В кинематическом режиме временной ряд поправок координат чаще всего представляется моделью белого шума, то есть поправки оцениваются независимо на каждом шаге. Также в кинематическом режиме изменение поправок координат может представляться случайным блужданием в том случае, если известны предельные скорости движения позиционируемого объекта и возможно предсказать разброс его случайных перемещений между моментами измерений. Режим «стой-иди» является комбинацией поочерёдно используемых статического и кинематического режимов.

Параметру ШВП dtA соответствует модель белого шума, значения параметра в различные моменты времени полагаются статистически независимыми.

Влажная составляющая зенитной тропосферной задержки Tzw представляется случайным блужданием согласно формуле где Tk, k 1 – шаг случайного блуждания, случайная величина с нулевым математическим ожиданием.

Дисперсия шага случайного блуждания между эпохами tk 1, tk определяется выражением где параметр q для влажной зенитной тропосферной задержки может быть принят равным 310-8 м2/с, как в программном обеспечении GPS Toolkit [115].

Псевдонеоднозначность BAs представляется кусочно-постоянным процессом:

постоянным на интервале времени, в течение которого НКА s наблюдается с объекта A непрерывно без потерь счёта циклов. При наличии не восстановленной потери счёта циклов на интервале времени между двумя измерениями tk 1, tk значения параметра BAs в эти моменты времени считаются независимыми.

1.2.7 Рекуррентное оценивание параметров Применяемый в методике ТДП для рекуррентного оценивания параметров фильтр Калмана, представленный в статье [92], подробно рассмотрен в работе [5].

Далее приведён модифицированный алгоритм фильтра Калмана, применённый в реализации методики ТДП в программном обеспечении GPS Toolkit со ссылкой на работу [58], обладающий большей численной устойчивостью, согласно [112].

Рекуррентное оценивание параметров системы уравнений (51) на каждом шаге выполняется в 2 этапа. Первый этап – экстраполяция:

где xk – вектор прогноза состояния системы на шаге k, ковариационная матрица которого K k ;

xk 1 – оценка вектора состояния системы на шаге k –1, характеризуемая ковариационной матрицей K k 1.

Для начала рекуррентного процесса оценивания формально задаётся нулевой вектор прогноза x0 и его диагональная ковариационная матрица K 0, диагональные элементы которой – величины порядка 1014.

Второй этап – коррекция, оценивание вектора состояния в эпоху tk по методу максимума апостериорной вероятности:

Оценка точности решения выполняется по формулам:

где 2 – масштабный коэффициент ковариационных матриц на шаге k;

xk i – среднеквадратическая погрешность оценки i-ого элемента вектора состояния системы в эпоху tk ;

nk – общее число кодовых и фазовых измерений в эпоху tk.

Описанный алгоритм представляет прямой ход оценивания и предполагает последовательное уточнение оценок параметров по мере поступления новой измерительной информации. В случае ГНСС-измерений накопление новой полезной информации зависит не только от количества измерений, но и требует некоторого времени: со временем ослабевает корреляция между случайными погрешностями ГНСС-измерений, меняется геометрия созвездия НКА, за счёт чего улучшается обусловленность решаемой задачи. Продолжительность сеанса измерений, при обработке которого достигается заданный уровень точности оценивания параметров, называется временем сходимости (convergence time). При прямом ходе оценивания первые оценки параметров будут иметь низкую точность, что в некоторых случаях неприемлемо. В таких случаях в постобработке оценивание может выполняться дважды с использованием той же измерительной информации: в прямом и обратном направлении. Такой подход позволяет достичь равномерной точности временных рядов оценок параметров.

1.3 Направления развития методики точного дифференциального позиционирования В настоящее время методика ТДП активно совершенствуется, разрабатываются новые её варианты. Основные направления развития – расширение спектра применений методики, повышение надёжности, доступности и оперативности определения координат.

Методика ТДП позволяет оценивать геометрическое положение точки в ОЗСК, которое может быть представлено геодезическими координатами, в том числе геодезической высотой. В то же время для решения многих практических задач необходимо определение нормальной высоты. Расширение сферы применения систем дифференциальной коррекции может быть достигнуто за счёт разработки методики спутникового нивелирования на основе ТДП. Такая методика позволила бы наряду с геометрическим положением позиционируемой точки определять её нормальную высоту и уклонения отвесной линии от нормали к эллипсоиду. В настоящее время спутниковое нивелирование выполняется с использованием моделей квазигеоида [9, 19]. В СГГА разработан теоретически более строгий метод совмещения пространственной геометрической и гравиметрической информации [46, 47, 48]. В будущем, при условии использования в составе аппаратуры потребителя ГНСС генератора частоты с относительной суточной нестабильностью порядка 10-18 появится теоретическая возможность выполнения спутникового нивелирования, основанного на зависимости темпа времени в точке установки генератора частоты от потенциала поля силы тяжести [61, 30]. В этом случае определение нормальной высоты не потребует привлечения модели гравитационного поля. Прототип генератора частоты с указанным уровнем нестабильности уже существует [79].

Одним из недостатков методики ТДП при использовании её для позиционирования подвижных объектов в кинематическом режиме является необходимость повторной инициализации в случае потери видимости спутников. Для решения указанной проблемы существуют способы быстрого восстановления потерь счёта циклов фазовых ГНСС-измерений. Некоторые из этих способов представлены в статьях [54, 65, 68, 119].

Повышение надёжности определения местоположения в кинематическом режиме, а также определение ориентации объекта в пространстве и уточнение характеристик его движения возможно за счёт совместного использования данных аппаратуры потребителя ГНСС и микроэлектромеханических систем (МЭМС):

гироскопов, акселерометров, магнитометрических датчиков (компасов). Примеры представлены в работах [52, 69]. При выборе конфигурации размещения МЭМС и ГНСС-антенн на позиционируемом объекте возможно использование принципа рычага: разнородные данные МЭМС и ГНСС могут быть использованы для взаимного контроля и повышения надёжности определения местоположения, скорости и ориентации объекта [90]. Частным случаем использования инерциальных навигационных систем на основе МЭМС для повышения надёжности решений по методике ТДП является применение их данных при восстановлении потерь счёта циклов в результатах фазовых ГНСС-измерений [68].

Важным направлением совершенствования методики ТДП является обеспечение возможности использования различных ГНСС. В настоящее время существуют примеры реализаций методики ТДП для разных комбинаций систем ГЛОНАСС, GPS, Compass, Galileo [60, 99, 105]. Объединение требует пересмотра модели результатов измерений. В частности, совместная обработка результатов измерений GPS и ГЛОНАСС требует включения дополнительного оцениваемого параметра во все уравнения измерений ГЛОНАСС – относительной поправки шкал системного времени [60]. Совместная обработка результатов измерений нескольких ГНСС позволяет значительно сократить время сходимости [60, 99], более чем на 50 % для GPS/ГЛОНАСС в сравнении с GPS [60]. Авторами [60] также экспериментально продемонстрировано достижение уровня точности оценивания координат порядка нескольких сантиметров с использованием исключительно ГЛОНАСС-измерений. Совместное использование сигналов нескольких ГНСС повышает доступность определения координат в условиях, когда значительная часть небесной сферы закрыта, например в городах.

Удешевление процесса определения координат может быть достигнуто за счёт использования аппаратуры потребителя ГНСС топографического или навигационного класса, в том числе одночастотной, цена которой может быть на один – два порядка ниже цены аппаратуры геодезического класса. При использовании одночастотных измерений возникает проблема коррекции ионосферной задержки.

Поправки, компенсирующие влияние ионосферы, рассчитываются по измерениям региональной или глобальной сети опорных станций. При этом необходимо принимать во внимание, что погрешности глобальных моделей ионосферы составляют от 10 % до 50 % [82 с. 143]. Однако, согласно [116], одночастотный вариант методики ТДП позволяет достигать дециметрового уровня точности позиционирования и демонстрирует даже большую эффективность в течение первых 10 мин сеанса измерений в сравнении с исходным вариантом методики, предполагающим использование «ионосферно-свободных» комбинаций результатов ГНССизмерений. Главным фактором, ограничивающим точность определения координат с помощью малобюджетной ГНСС-аппаратуры, остаётся качество этой аппаратуры, в частности способность к подавлению многопутности радиосигналов [55]. Достижимый уровень точности оценки пространственного положения по методике ТДП в статическом режиме с использованием аппаратуры потребителя навигационного класса, согласно [55], характеризуется СКП порядка полуметра при малом влиянии многопутности.

Важным шагом в направлении повышения доступности и оперативности высокоточного позиционирования на основе методики ТДП является бесплатное предоставление службой МГС уточнённой ЭВИ в реальном времени через сеть Интернет по протоколу NTRIP. В настоящее время предоставляется уточнённая ЭВИ GPS и ГЛОНАСС (в экспериментальном режиме) [87].

Ещё один способ совершенствования методики ТДП – применение разрешения неоднозначностей фазовых ГНСС-измерений или их одинарных разностей.

Существует несколько различных вариантов реализации этого подхода, рассмотренных в работах [37, 62, 74, 75, 97, 100, 117]. Положительный эффект разрешения и фиксации неоднозначностей заключается в повышении точности позиционирования в кинематическом режиме и значительном сокращении времени сходимости как в статическом [62, 75], так и в кинематическом режиме [97]. Разрешение неоднозначностей результатов фазовых измерений требует учёта дополнительных поправок, определяемых при обработке ГНСС-измерений сети опорных станций, и задержек сигналов в аппаратуре.

Повышение точности определения местоположения по методике ТДП и сокращение времени сходимости также могут быть достигнуты за счёт повышения адекватности применяемой модели результатов ГНСС-измерений, и в частности модели БШВ НКА. Данный вопрос подробно рассмотрен во втором разделе.

2 АНАЛИЗ И СОВЕРШЕНСТОВАНИЕ МАТЕМАТИЧЕСКОЙ МОДЕЛИ

ДИНАМИЧЕСКОЙ СИСТЕМЫ И АЛГОРИТМА МЕТОДИКИ ТОЧНОГО

ДИФФЕРЕНЦИАЛЬНОГО ПОЗИЦИОНИРОВАНИЯ

2.1 Статистический анализ погрешностей математической модели ГНСС-измерений 2.1.1 Корреляционный анализ остаточных невязок уравнений фазовых измерений Автором была предпринята попытка реализовать методику ТДП с разрешением неоднозначностей фазовых измерений на основе работы [75]. В ходе тестирования алгоритма не удалось достичь стабильного разрешения неоднозначностей фазовых измерений, так как СКП остаточных невязок уравнений фазовых измерений были сопоставимы с половиной длины волны суммарной комбинации фазовых измерений. Также на ненадёжность оценок параметров по методике ТДП с разрешением неоднозначностей указывали разработчики программного обеспечения RTKLIB: этот алгоритм, включённый в версию 2.4.2, является экспериментальным и нестабильным [108 с. 176]. Ошибочная фиксация набора неоднозначностей приводит к возникновению систематической погрешности фазовых измерений, кратных длине волны суммарной комбинации фазы. Согласно [62], процесс оценивания координат может расходиться, и одной из вероятных причин этого являются не выявленные потери счёта циклов в результатах фазовых измерений. Потери счёта циклов могут маскироваться случайными погрешностями при высоком уровне шума фазовых измерений. Также разброс случайных погрешностей влияет на время сходимости плавающего решения.

С целью поиска возможностей уменьшения шума фазовых измерений проведён анализ остаточных невязок уравнений фазовых измерений дальности до нескольких НКА за 27.05.2008. Использовалась окончательная ЭВИ МГС: эфемериды и таблицы поправок БШВ с пятиминутным и 15-минутным интервалами. Было выполнено оценивание параметров по методике ТДП с использованием результатов GPS-измерений с трёх станций МГС: Мадрид (MADR), Турин (IENG), Веттцель (WTZZ). В результате были получены значения остаточных невязок уравнений фазовых измерений, представленные на графиках (рисунок 3).

Остаточная невязка, м -0, -0, -0, 23:00 23:05 23:10 23:15 23:20 23:25 23:30 23:35 23:40 23:45 23:50 23: Остаточная невязка, м -0, -0, 23:00 23:05 23:10 23:15 23:20 23:25 23:30 23:35 23:40 23:45 23:50 23:

MADR IENG WTZZ

Рисунок 3 – Остаточные невязки уравнений фазовых измерений НКА GPS 21:

а) при пятиминутном интервале задания поправок БШВ;

б) при 15-минутном интервале задания поправок БШВ.

Результат анализа статистической связи остаточных невязок уравнений фазовых ГНСС-измерений на разных НИП представлен в таблице 1.

Таблица 1 – Корреляция остаточных невязок уравнений синхронных фазовых измерений дальности до НКА GPS Коэффициент корреляции остаточных невязок на удалённых станциях близок к единице, что свидетельствует о существовании эффекта, проявляющегося практически одинаково на всех НИП, синхронно наблюдающих определённый НКА. В то же время, коэффициент корреляции между остаточными невязками, полученным при использовании таблиц БШВ с пяти- и 15-минутным интервалами, оказался значительно меньше и равнялся 0,26. За исключением табличного интервала, условия экспериментов были идентичными. Это позволило сделать вывод о том, что проявляющийся в остаточных невязках эффект обусловлен погрешностями интерполяции поправок БШВ.

2.1.2 Анализ погрешностей интерполяции поправок бортовых шкал времени навигационных космических аппаратов Модель БШВ НКА представляется таблицей поправок, заданных с определённым интервалом времени, в составе нескольких типов продуктов МГС: прогнозных (ultra-rapid predicted half), «сверхбыстрых» предварительных (ultra-rapid observed half), предварительных (rapid), окончательных (final). Характеристики продуктов МГС, включающих поправки БШВ, представлены в таблице 2 в соответствии с их описанием [84].

по данным МГС Задержка в реальном публикации времени Примечание – Среднеквадратическое отклонение (СКО) учитывает систематическую погрешность бортовой шкалы времени, в отличие от СКП.

Как указано на сайте МГС [84], служба также предоставляет окончательную ЭВИ ГЛОНАСС с задержкой от 12 до 18 суток. СКП поправок БШВ указываются непосредственно в файлах, содержащих ЭВИ.

Согласно исследованиям [77, 121], краткосрочная нестабильность генераторов частоты НКА ГНСС снижает точность представления поправок БШВ в виде таблиц с интервалом более 30 с. В работе [77] авторы рассматривают влияние интервала задания поправок БШВ на характеристики решения по методике ТДП и показывают, что сокращение интервала с 5 мин до 30 с не оказывает значительного влияния на точность определения местоположения при статическом позиционировании, но позволяет достичь значительного (от 30 % до 50 %) уменьшения погрешностей оценок координат в кинематическом режиме. При этом влияние табличного интервала на время сходимости не рассматривалось.

Для подтверждения сделанного ранее вывода о том, что проявляющийся в остаточных невязках эффект действительно обусловлен погрешностями интерполяции, было проведено экспериментальное исследование. Погрешности интерполяции поправок БШВ были определены по отклонениям интерполированных значений от эталонных. В качестве эталона при вычислении отклонений использовались табличные значения поправок БШВ, заданные с интервалом 30 с (в файлах clk_30s). Интерполированные значения были получены по таблицам с интервалами 5 мин (файлы clk) и 15 мин (файлы sp3). В этом исследовании использовались только окончательные продукты МГС (final). Таким образом, погрешности исходных данных были одинаковы во всех таблицах, а расхождения обусловлены лишь погрешностями интерполяции. Интерполяция выполнялась по методу Лагранжа [51 с. 103–106]:

где y x – интерполяционный многочлен степени n, проходящий через n + 1 узловую точку с известными значениями x j, y j.

Было проведено два эксперимента: интерполяция выполнялась многочленами третьей и девятой степени (девятая степень используется по умолчанию в применявшемся программном обеспечении GPS Toolkit).

На рисунке 4 сопоставлены графики средних значений остаточной невязки по трём станциям (Мадрид, Турин, Веттцель) и погрешностей кубической интерполяции поправок БШВ для НКА GPS 21 за один час 27.05.2008.

Значение, м -0, -0, -0, 22:00 22:05 22:10 22:15 22:20 22:25 22:30 22:35 22:40 22:45 22:50 22: Значение, м -0, -0, 22:00 22:05 22:10 22:15 22:20 22:25 22:30 22:35 22:40 22:45 22:50 22: Среднее из остаточных невязок на пунктах MADR, IENG, WTZZ Рисунок 4 – Остаточные невязки уравнений фазовых измерений и погрешности интерполяции поправок БШВ НКА GPS 21:

а) при пятиминутном интервале задания поправок БШВ;

б) при 15-минутном интервале задания поправок БШВ.

Корреляция погрешностей интерполяции поправок БШВ и средних значений остаточных невязок в примере, представленном на рисунке 4, а достигает 95 %, и 84 % в примере на рисунке 4, б. Таким образом, подтверждён вывод о том, что проявляющийся в остаточных невязках эффект, рассмотренный в предыдущем пункте, действительно обусловлен погрешностями интерполяции поправок БШВ.

Пример, приведённый на рисунке 4, позволяет идентифицировать один из источников погрешностей функциональной модели ГНСС-измерений. Однако оценка разброса погрешности интерполяции поправок БШВ и её вклада в общую погрешность модели измерений требует анализа более репрезентативной выборки.

Были определены среднеквадратические отклонения интерполированных значений поправок БШВ от табличных по данным всех НКА GPS за неделю GPS 1721. Объём выборки для каждого типа продуктов (sp3 и clk) составил более тыс. элементов. Среднеквадратические отклонения значений поправок БШВ, полученных в результате кубической интерполяции, от табличных значений с интервалом 30 c представлены в таблице 3 по дням.

Таблица 3 – Погрешности кубической интерполяции БШВ Среднеквадратические отклонения интерполированных Среднеквадратические отклонения, полученные при использовании интерполяционного многочлена девятой степени, превышают приведённые в таблице 3 на величины порядка миллиметра. Таким образом, из двух вариантов более предпочтительным является интерполяционной полином третьей степени.

В использованных при вычислении СКО данных встречаются разрывы по некоторым НКА, в несколько раз превышающие интервал интерполяции. Уровень погрешностей, продемонстрированный в таблице 3, достигается только в том случае, если не допускать интерполяции в разрывах. В противном случае погрешности могут возрастать до нескольких метров. В ходе эксперимента выявлены пропуски данных в таблицах поправок продолжительностью более часа для НКА GPS № 3, 6, 26.

Полученные значения среднеквадратических отклонений интерполированных значений от табличных в целом согласуются со значениями, представленными в статьях [77, 121].

При заданной степени многочлена погрешности интерполяции поправок БШВ зависят от стабильности генераторов частоты, реализующих БШВ НКА, табличного интервала и моментов, на которые производится интерполяция. Если номинальные моменты записи ГНСС-измерений аппаратурой потребителя совпадают с узлами интерполяции поправок БШВ, то дисперсию погрешностей интерполяции можно считать равной нулю. При проведении оценки отклонений интерполированных значений от табличных не исключались общие узлы, в которых отклонения были нулевыми. Значения СКО, приведённые в таблице 3 представляют общую оценку отклонений. С учётом того, что погрешности интерполяции независимы от погрешностей табличных значений, общую дисперсию учёта поправки БШВ в свободных членах уравнений (47), (48), (73), (74) можно представить формулой где Ddt – дисперсия погрешности учёта поправки БШВ в свободных членах;

DБШВ – дисперсия погрешности табличного значения поправки БШВ;

Dинт – дисперсия погрешности интерполяции, связанная с представленными в таблице 3 среднеквадратическими отклонениями инт общеизвестной формулой Дисперсия DБШВ может задаваться либо в соответствии с заявленной для используемого типа ЭВИ (таблица 2), либо в соответствии с оценкой точности, приведённой в файле ЭВИ, что предпочтительнее.

Погрешность расчёта поправки БШВ не зависит от погрешностей расчёта других поправок, учитываемых при вычислении свободных членов по формулам (49), (50). Таким образом, с учётом формулы (70) дисперсия погрешности модели ГНСС-измерений может быть выражена формулой где Dv – дисперсия общей погрешности модели ГНСС-измерений;

D – дисперсия общей погрешности модели ГНСС-измерений без учёта погрешности модели бортовых шкал времени НКА.

Таким образом, погрешность интерполяции вносит в общую случайную погрешность измерений дополнительный значимый вклад, который целесообразно уменьшить.

2.2 Совершенствование математической модели ГНСС-измерений 2.2.1 Способ уточнения функциональной модели фазовой псевдодальности Погрешности интерполяции поправок БШВ могут быть уменьшены за счёт использовании таблиц с меньшим интервалом, равным, например, 30 с. Однако таблицы поправок БШВ с таким интервалом доступны с задержкой от 12 до суток после выполнения измерений. «Сверхбыстрая» и предварительная ЭВИ могут быть получены потребителем с задержкой от нескольких часов до двух суток.

Однако эти типы продуктов включают поправки с интервалами 5 мин и более. Таким образом, при использовании этих типов ЭВИ погрешности вычисления свободных членов по формулам (49), (50) будут возрастать из-за вклада погрешностей интерполяции.

Автором предлагается способ компенсации погрешностей интерполяции поправок БШВ в уравнениях фазовых измерений. Способ был представлен в статье [29]. Он заключается в учёте дополнительной поправки фазовых измерений потребителя, рассчитываемой по результатам синхронных фазовых измерений, выполненных на опорных станциях. В качестве опорных предлагается использовать станции МГС, региональных сетей [23, 24, 25, 45] или временные станции, координаты которых могут быть определены потребителем. По измерениям на опорных НИП оцениваются параметры модифицированных уравнений ТДП:

Свободные члены уравнений (73), (74) вычисляются с использованием известных с высокой точностью координат опорных НИП. Поправка фазовых измерений t потребителя рассчитывается по остаточным невязкам v уравнений измерений (74) на опорных НИП по формуле где t – поправка фазового измерения дальности до НКА s в момент t, м;

m – количество опорных станций, наблюдавших НКА s в момент t;

vi t – остаточная невязка уравнения фазового измерения дальности от НИП i до НКА s, полученная в результате оценивания.

«Ионосферно-свободные» комбинации синхронных фазовых измерений на опорных станциях полагаются независимыми и равноточными (см. далее в 2.2.2), остаточные невязки уравнений фазовых измерений на этих станциях не коррелированы и характеризуются одинаковыми весами. Поэтому среднеквадратическая погрешность (СКП) s t поправки t вычисляется по формуле Бесселя где v t – среднее арифметическое остаточных невязок уравнений фазовых изs мерений дальности до НКА s в эпоху t.

Отбраковку наиболее грубого значения остаточной невязки vi t предлагаs ется выполнять при условии:

где квантиль предлагается принимать равным 2s t.

После отбраковки грубых значений остаточных невязок поправка t и СКП s t вычисляются заново по формулам (75), (76). Поправки t рассчиs тываются на каждый момент измерений на позиционируемом объекте по каждому НКА и вводятся в свободные члены уравнений фазовых ГНСС-измерений потребителя:

2.2.2 Способ уточнения априорной стохастической модели погрешностей ГНСС-измерений Как уже отмечалось в 1.2.6, существующие априорные стохастические модели погрешностей ГНСС-измерений не учитывают таких значимых факторов, как точность применяемой ЭВИ, величина табличного интервала представления поправок БШВ, тип используемых кодовых псевдодальностей. Кроме того, учёт дополнительной поправки t, предложенной в 2.2.1, должен уменьшать дисперs сию погрешностей фазовых измерений, что должно быть отражено при формировании соответствующей априорной ковариационной матрицы. Необходимо уточнить стохастическую модель таким образом, чтобы она адекватно представляла погрешности ГНСС-измерений.

На основе способов задания априорных дисперсий, рассмотренных в 1.2.6, и формулы дисперсии погрешности учёта поправки БШВ (72) предлагается уточнённая формула априорной дисперсии погрешности модели измерений:

Стандартное отклонение погрешности табличных значений поправок БШВ рекомендуется задавать в соответствии с заявленным значением среднеквадратической погрешности данного типа продуктов МГС (см. таблицу 2) или в соответствии с оценкой точности, приведённой в файлах поправок БШВ.

Стандартное отклонение погрешности интерполяции инт для таблиц поправок БШВ GPS с 15-минутным интервалом задаётся равным 49 мм, для таблиц с пятиминутным интервалом – равным 31 мм. При совпадении интервала задания поправок БШВ и интервала записи измерений аппаратурой потребителя стандартное отклонение инт равно нулю. При расчёте дисперсии погрешностей уравнений фазовых измерений, в которых учтена поправка t, стандартное отклонение инт также принимается равным нулю. Так как измерения на позиционируемом объекте и сети опорных станций выполняются синхронно, предполагается, что учёт поправок t полностью компенсирует вклад погрешности инs терполяции в общую погрешность модели фазовых измерений.

Значения номинального стандартного отклонения n и функции от угловой высоты f h были подобраны эмпирически для различных типов измерительной информации GPS на основе анализа остаточных невязок уравнений измерений из оценивания по методике ТДП. Параметры формулы (79) для четырёх типов измерительной информации GPS, используемой в методике ТДП, представлены в таблице 4.

Таблица 4 – Параметры уточнённой априорной стохастической модели погрешностей «Ионосферно-свободная» комбинация высокоточных кодовых измерений P1, P2 (GPS) «Ионосферно-свободная» комбинация кодовых измерений станsin 2 h дартной и высокой точности С1, P2 (GPS) Исходная «ионосферно-свободная» комбинация фазовых измерений (GPS) «Ионосферно-свободная» комбинация фазовых измерений с учёsin 2 h том поправки t (GPS) Оценивание выполнялось по результатам суточных сеансов GPS-измерений 12 станций МГС за 01.01.2013. Не использовались остаточные невязки, полученные в течение времени сходимости (первых трёх часов каждого суточного сеанса измерений). Критериями выбора функции f h являлись:

независимость масштабного коэффициента ковариационной матрицы погрешностей измерений от угловой высоты НКА;

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

Оценка эффекта от применения уточнённой априорной стохастической модели представлена в 3.2.4.

2.3 Усовершенствованная методика точного дифференциального позиционирования и её программная реализация Усовершенствованная методика ТДП представлена на рисунке 5. К её особенностям можно отнести простоту и возможность сохранения высокого уровня автоматизации процесса обработки ГНСС-измерений. Расчёт поправок t моs жет выполняться как с помощью вычислительных средств потребителя, так и централизованно. Исходные данные для расчёта поправок t доступны из тех же источников, что и для исходной методики ТДП, уже на следующие сутки после выполнения измерений.