Результат интеллектуальной деятельности: СПОСОБ ОПРЕДЕЛЕНИЯ И ПРОГНОЗИРОВАНИЯ ДВИЖЕНИЯ КОСМИЧЕСКОГО АППАРАТА НА НИЗКИХ ОРБИТАХ, ПОДВЕРЖЕННОГО ВЛИЯНИЮ ТОРМОЖЕНИЯ В АТМОСФЕРЕ

Изобретение относится к области космонавтики, а именно к технике определения и прогнозирования параметров движения космических аппаратов (КА), подверженных влиянию торможения в атмосфере, по траекторным измерениям на заданном мерном интервале.

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

Из уровня техники известны традиционные технологии определения и прогнозирования движения космического аппарата на низких орбитах, подверженного торможению вследствие вариаций плотности верхней атмосферы (см. например [1-6]). Однако за последние 30 лет точность прогнозирования торможения низких спутников повысить не удалось. В лучшем случае погрешности находятся на уровне примерно 10% от величины атмосферных возмущений на интервале прогнозирования.

Недостатками известных из уровня техники технических решений являются:

1) пренебрежение данными о возможных вариациях торможения спутников в атмосфере;

2) ограничение применяемого на практике значения мерного интервала, величина которого зависит от уровня торможения КА в атмосфере (высоты КА).

Традиционный подход к определению начальных условий для прогнозирования движения, современные методы интегрирования уравнений движения с учетом разных возмущений, а также методы оценки вариаций плотности атмосферы разрабатывались разными специалистами и в разное время - без учета их взаимосвязи и возможности объединения в одну общую проблему. Заявленное техническое решение направлено на устранение указанных недостатков. В заявленном изобретении решение перечисленных задач выполняется с учетом их взаимосвязи. Это позволяет сделать следующий шаг в решении проблемы повышения точности определения и прогнозирования орбит в атмосфере, т.е. преодолеть упомянутый выше 10%-й уровень погрешностей прогнозирования.

Технический результат заявленного технического решения заключается в повышении точности определения и прогнозирования орбит КА. Методической основой создания способа определения и прогнозирования торможения КА является модификация метода максимального правдоподобия [3, 7, 8], названная методом оптимальной фильтрации измерений.

Технический результат достигается тем, что в способе определения и прогнозирования торможения космического аппарата на низких орбитах вследствие вариаций плотности верхней атмосферы измеряют траекторные параметры спутника на мерном интервале и записывают в их массив исходных измерений на заданном мерном интервале, который является основной исходной информацией для определения и прогнозирования торможения КА; на основе массива данных для адаптации, содержащего результаты предыдущих оценок баллистического коэффициента (S_b) по рассматриваемому спутнику, выполняют адаптивное уточнение среднего значения баллистического коэффициента нормированного среднеквадратического отклонения (k_atm) и интервала корреляции его вариаций (Δ). При недостаточности данных для адаптации используют соответствующие априорные характеристики торможения; массив априорных данных включает в себя также число измерений на мерном интервале, значения среднеквадратических погрешностей компонентов вектора измерений, массив средних значений индексов солнечной и геомагнитной активности, а также массив заданных пользователем интервалов прогнозирования; далее рассчитывают весовую матрицу максимального правдоподобия с учетом статистических характеристик погрешностей прогнозирования движения на мерном интервале, формируют весовую матрицу, которая учитывает погрешности измерений и влияние случайных атмосферных возмущений, вычисляют оценку вектора состояния на момент времени последнего измерения без учета атмосферного торможения, вычисляют оценки текущих и прогнозных значений баллистического коэффициента на основе остаточных невязок и на основе корреляции его вариации; полученные оценки Sbк записывают в массив априорных данных, далее прогнозируют движение космического аппарата с учетом обращения к стандартным модулям интегрирования уравнений движения и расчета статистических характеристик погрешностей при прогнозировании движения.

Признаки и сущность настоящего изобретения поясняются в последующем детальном описании, иллюстрируемом чертежами (фиг.1, 2), где показано следующее:

Фиг.1. Блок-схема алгоритма реализации способа определения и прогнозирования торможения КА на низких орбитах вследствие вариаций плотности верхней атмосферы, на основе процессора, содержащая следующие блоки:

1 - КА;

2 - блок измерения траекторных параметров;

3 - блок априорные и адаптивные данные;

4 - блок адаптивное уточнение характеристик торможения;

5 - блок прогнозирование статистических погрешностей на мерном интервале;

6 - блок расчет весовой матрицы максимального правдоподобия;

7 - блок уточнения значения шестимерного вектора состояния;

8 - блок определение характеристик торможения на основе остаточных невязок;

9 - блок прогнозирование торможения на основе корреляции его вариаций;

10 - блок прогнозирование орбитальных параметров.

Способ реализуется следующим образом. В сеансе связи с КА (поз.1, фиг.1) выполняют измерения траекторных параметров и записывают их в массив (поз.2, фиг.1) исходных измерений на заданном мерном интервале, который является основной исходной информацией для определения и прогнозирования торможения КА. В общем случае состав измерений может быть разным. Испытания заявленной системы по реальной информации были проведены при форме измерений 6-мерного вектора состояния в разные моменты времени. Также выполняют адаптивное уточнение статистических характеристик торможения (поз.4, фиг.1) на основе массива данных для адаптации (поз.3, фиг.1), содержащего результаты предыдущих решений рассматриваемой задачи, а именно оценки баллистического коэффициента (S_b) по рассматриваемому спутнику в различные моменты времени (t), а также значения соответствующих номеров витков (N). Результатами адаптивного уточнения являются: среднее значение баллистического коэффициента нормированное среднеквадратическое отклонение (k_atm) и интервал корреляции его вариаций (Δ), в течение которого нормированная автокорреляционная функция вариаций баллистического коэффициента уменьшается до значений 0.1-0.2. Эти результаты записываются в массив априорных данных. После каждого решения рассматриваемой задачи он пополняется. Адаптация выполняется тогда, когда накоплены данные не менее чем за 1-2 месяца. На начальном этапе работы этот массив может быть недостаточным для адаптации. В этом случае используются априорные данные. Далее осуществляют расчет статистических характеристик погрешностей прогнозирования движения на мерном интервале (поз.5, фиг.1) и расчет весовой матрицы метода максимального правдоподобия (поз.6, фиг.1). Данная операция позволяет увеличить мерный интервал по сравнению с его традиционным значением и на этой основе повысить точность определения и прогнозирования движения КА в атмосфере. В основу оценки вектора состояния спутника на основе измерений z_i=h_i·x(t_i)+v_i, i=1,…k, в моменты времени t_i, где h_i - известные матрицы, v_i - случайные погрешности, положена зависимость

Здесь: Х_k - матрица частных производных вектора измерений Z_k по компонентам вектора состояния х_k, весовая матрица P_k определяется по формуле

Здесь матрица H имеет вид

Выражение (1) имеет традиционную форму. Отличие состоит в формировании весовой матрицы Р_k, которая учитывает не только погрешности измерений (матрица R_Σ), но и влияние случайных атмосферных возмущений (матрица K_wk). В частном случае, при отсутствии случайных возмущений формула (1) принимает форму классического метода наименьших квадратов. Компоненты матрицы K_wk определяются путем обращения к стандартному модулю расчета статистических характеристик погрешностей при прогнозировании движения. Далее осуществляют оценку вектора состояния на момент времени последнего измерения (поз.7, фиг.1). При разработке данного способа учтено важное положение, обоснование которого изложено в статье [8]; в ней показано, что при любом уровне случайных возмущений наилучшую точность обеспечивает применение метода оптимальной фильтрации измерений. Из этого общего положения следует рекомендация - не включать характеристику торможения (баллистический коэффициент) в состав вектора состояния, который уточняется в процессе оптимизации критерия максимального правдоподобия (максимум условного апостериорного распределения p(x(t)|z₁, z₂,…z_k|. Правомерность этой рекомендации является следствием выполнения операций расчета весовой матрицы максимального правдоподобия (поз.6, фиг.1) с учетом случайных атмосферных возмущений. Хотя формула для оценки вектора состояния имеет традиционный вид (1), ее применение обеспечивает минимизацию погрешностей определения вектора x(t_k) в последней точке мерного интервала - с лучшей точностью, чем в методе наименьших квадратов. Таким образом, в стандартном модуле интегрирования уравнений движения используется среднее значение баллистического коэффициента, т.е. оценка характеристик атмосферного торможения не входит в процедуру уточнения шестимерного вектора состояния спутника. Она определяется по остаточным невязкам между измеренными и расчетными значениями параметров орбиты (поз.8, фиг.1). Такая технология является принципиальной особенностью рассматриваемой системы.

Выходными данными уточнения значения шестимерного вектора состояния (поз.7, фиг.1) являются: оценка вектора состояния на момент последнего измерения и соответствующая корреляционная матрица погрешностей, вычисляемая по формуле

Далее осуществляется оценка текущих и прогнозных значений баллистического коэффициента на основе остаточных невязок и на основе корреляции его вариации (поз.8, 9, фиг.1). Содержание выполняемых операций основано на том обстоятельстве, что в заявленном способе результирующие остаточные невязки между измерениями и их расчетными значениями сильно отличаются от результатов применения метода наименьших квадратов. Так как в процессе обработки измерений в модели движения спутника используется некое среднее (постоянное) значения баллистического коэффициента то остаточные невязки отражают влияние не учтенных атмосферных возмущений на мерном интервале. Отличительной особенностью невязок является их уменьшение на несколько порядков по мере приближения к последней точке мерного интервала (моменту начальных условий). Данное обстоятельство является следствием рассмотренных выше особенностей заявленного способа. Отмеченная характерная особенность изменения остаточных невязок на мерном интервале является основой уточнения характеристик торможения спутника. Для каждого из уточнений выполняется полиномиальная аппроксимация остаточных временных невязок - разностей между измеренным и расчетным временем попадания спутника в данную точку. Эти остаточные невязки обозначены ниже как delta z(N_i) где N_i - номер витка i-го измерения. Пример значений остаточных невязок при обработке реальных измерений КА (в каталоге НОРАД) представлен на фигуре 2. Так как величина текущих атмосферных возмущений пропорциональна второй производной временных невязок delta z(N_i), как функции номера витка, то степень полинома должна быть не ниже второй. Используются два вида полиномов, которые аппроксимируют остаточные невязки: полином второй

и третьей степени

Вторая производная этих полиномов в конце мерного интервала используется для оценки текущего значения баллистического коэффициента

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

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

где ΔN - интервал прогнозирования в витках.

Далее вычисляют прогноз орбитальных параметров (поз.10, фиг.1), в котором осуществляется обращение к стандартным модулям интегрирования уравнений движения и расчета статистических характеристик погрешностей при прогнозировании движения. Интервал прогнозирования задается пользователем (поз.3, фиг.1); использование прогнозных значений (8) баллистического коэффициента является еще одной отличительной особенностью заявленного способа, тем самым обеспечивая оптимальные условия взаимодействия между операциями (поз.4, 5 и 6, фиг.1).

Стандартный модуль интегрирования уравнений движения реализован в виде специального программного модуля. Используется численное интегрирование уравнений движения. Применяется вектор состояния (t, R, V) в форме компонентов положения и скорости в геоцентрической инерциальной системе координат. Учитываются следующие возмущающие факторы:

- параметры разложения гравитационного потенциала Земли по полиномам Лежандра до 8-го порядка;

- торможение в атмосфере, плотность которой рассчитывается в соответствии с моделью ГОСТ Р 25645.166-2004.

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

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

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

- возможность управления составом параметров гравитационного поля Земли и, в частности, пренебрежением тессеральными гармониками в случае применения измерений в форме т.н. двух рядных элементов (TLE).

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

Исходными данными для применения этой корреляционной функции являются:

ΔT - изменение периода под действием торможения в атмосфере за виток, которое рассчитывается он основе численного интегрирования при среднем значении баллистического коэффициента;

k_atm - СКО случайных атмосферных возмущений по отношению к их среднему значению.

Δ - интервал корреляции атмосферных возмущений (витки).

Первые две величины используется для вычисления СКО вариаций атмосферного торможения по формуле

Вычисление матриц взаимной корреляции погрешностей прогнозирования вектора состояния на моменты времени (t_i и t_l) выполняется по формуле

Здесь U(t_i,ξ) - т.н. переходная матрица размерностью (6×6), В(ξ) - матрица коэффициентов при атмосферном торможении в дифференциальных уравнениях возмущенного движения. Полная формула (11) при t_i=t₁=t используется только при расчете весовой матрицы. При вычислении матрицы К_wл используется только второе слагаемое:

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

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

Испытания заявленного способа подтвердили выводы о его преимуществе исходя из результатов обработки измерений по сравнению с его аналогами.

Список литературы

1. Э.Л.Аким, Т.М.Энеев. Определение параметров движения космических аппаратов по траекторным измерениям. Космические исследования, т.1, №5, М., 1963.

2. П.Е.Эльясберг. Определение движения по результатам измерений. Изд-во «НАУКА», М., 1976.

3. А.И.Назаренко и Б.С.Скребушевский. Эволюция и устойчивость спутниковых систем, «МАШИНОСТРОЕНИЕ». М., 1981.

4. Патент RU 2150414 С1 (ГУП ЦНИИ «КОМЕТА»), опубл. 01.06.1999, «Способ определения параметров орбиты космического аппарата.

5. D.A.Vallado, Fundamentals of Astrodynamics and Applications. Published jointly by Microcosm Press and Kluwer Academic Publishers, 2004.

6. Патент RU 2391265 C1 (ОАО «Российские космические системы), опубл. 27.08.2010, «Система оценивания точности прогнозирования параметров движения искусственного спутника Земли».

7. А.И.Назаренко. Погрешности прогнозирования движения спутников в гравитационном поле Земли, Москва, Институт космических исследований РАН, 2010.

8. A.I.Nazarenko. Accuracy of Determination and Prediction Orbits in LEO. Estimation Errors Depending on Accuracy and Amount of Measurements, Seventh US/Russian Space Surveillance Workshop, Monterey, October-November 2007.