Адаптивная обработка данных авиационной гравиметрии тема автореферата и диссертации по механике, 01.02.01 ВАК РФ

МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ ИМЕНИ М.В. ЛОМОНОСОВА Механико-математический факультет

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

ДАННЫХ АВИАЦИОННОЙ ГРАВИМЕТРИИ

01.02.01 - теоретическая механика

АВТОРЕФЕРАТ

диссертации на соискание ученой степени кандидата физико-математических наук

2 2 НОЯ 2012

Москва - 2012

005055569

Работа выполнена на кафедре прикладной механики и управления механико-математического факультета Московского государственного университета имени

М.В. Ломоносова.

Научный руководитель:

Болотин Юрий Владимирович, доктор физико-математических наук, профессор

Официальные оппонепты:

Семенихин Константин Владимирович, доктор физико-математических наук, профессор

Ведущая организация:

Каршаков Евгений Владимирович, кандидат физико-математических наук, старший научный сотрудник

Институт физики Земли им. О.Ю. Шмидта Российской академии наук

Защита состоится 14 декабря 2012 г. в 16 часов 30 минут на заседании совета Д.501.001.22 при Московском государственном университете имени М.В. Ломоносова по адресу: 119991, Москва, Ленинские горы, Главное здание МГУ, механико-математический факультет, ауд. 16-10.

С диссертацией можно ознакомится в читальном зале отдела диссертаций (Ломоносовский просп., 27, Фундаментальная библиотека, сектор А — 8 этаж, к.812)

Автореферат разослан 14 ноября 2012 г.

Ученый секретарь диссертационного совета кандидат физико-математических наук, доцент

Прошкин В.А.

Общая характеристика работы

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

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

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

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

В лаборатории управления и навигации механико-математического факультета МГУ под руководством H.A. Парусникова, A.A. Голована и Ю.В. Болотина начиная с 1995 года ведутся работы по разработке математического и программного обеспечения для решения задачи авиационной гравиметрии. Диссертация является продолжением этих работ.

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

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

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

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

Теоретическая и практическая ценность.

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

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

Обсуждение работы и публикации. Результаты диссертации докладывались на международной конференции "Современные проблемы математики, механики и их приложений"(Москва, 2009 г.), на симпозиуме Международной ассоциации по геодезии (IAG, С.-Петербург, 2010 г.), на конгрессе международной федерации автоматического управления (The 18th World Congress of the International Federation of Automatic Control (IFAC), Италия, 2011 г.), на научном семинаре им. АЛО. Ишлинского (МГУ, 2009, 2010, 2011, 2012 г.), и др.

Основные результаты диссертации опубликованы в пяти работах. Список публикаций приведен в конце автореферата. Работа над диссертацией выполнялась при поддержке РФФИ (проект 10-01-00703-а).

Структура и объем диссертации. Работа состоит из списка обозначений, трех глав, заключения и списка литературы (88 наименований). В работе приведено 24 рисунка, 5 таблиц. Общий объем диссертации составляет 112 страницы.

2 Содержание работы

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

тивной фильтрации. В четвертой части рассматривается вопрос о возможности применения различных адаптивных фильтров к ЗАГТ.

В гравиметрии силу тяжести д принято представлять в виде суммы нормальной силы тяжести до и аиомалии силы тяжести 6д:

д = до +

Нормальное поле силы тяжести до представляет поле модельного эллипсоида вращения, поверхность которого близка к поверхности Земли. В гравиметрии в качестве референц-эллипсоида часто используется эллипсоид WGS-84. Эллипсоид вращения определяет систему координат (p,X,h, где <р — географическая широта, А — географическая долгота, h — высота над эллипсоидом. Нормальная сила тяжести до(ф, h) является известной функцией географической широты <р и высоты h, и вычисляется по формуле Сомильяны1. Аномалия силы тяжести (ACT) определяет отклонение реального поля от нормального. Основное влияние на структуру ACT района съемок оказывает локальное распределение притягивающих масс. ACT является неизвестной непараметрической функцией широты, долготы, и высоты над референц-эллипсоидом 5g(ip, А, К).

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

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

lrIbrge W. Geodesy, 2nd Edition Walter de Gruyter., Berlin, 1991.

1. Определение координат и скоростей антенн бортовых приемников СНС.

2. Определение угловых ошибок горгоонтирования платформы ИНС.

3. Определение координат и скоростей чувствительного элемента гравиметра.

4. Решение основного гравиметрического уравнения — ЗАГТ.

5. Согласование галсов (устранение трендов гравиметра).

6. Построение карты ACT на высоте движения летательного аппарата.

7. Редукция карты на Землю, учет топографических поправок и т.д.

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

Ошибки горизонтирования платформы гравиметра приводят к отклонению вертикального акселерометра от географической вертикали. В процессе решения угловые ошибки горизонтирования платформы оцениваются по данным ИНС. Данная информация позволяет вычислить значение вертикальной проекции удельной силы, действующей на чувствительного элемента, за счет перепроектирования показаний акселерометров.

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

Далее в диссертации рассматривается ЗАГТ. Сила тяжести на траектории рассматривается как функция полетного времени t:

9(t) = 9о{Ф), /»(*)) + 6g{<p{t), A (t), h(t)).

Основное уравнение скалярной авиационной гравиметрии, являющееся уравнением движения чувствительного элемента гравиметра в проекции на геогра-

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

V = f + go + Sg + AgEiv (1)

Здесь / — вертикальная проекция удельной силы, действующей на чувствительный элемент со стороны подвеса; V — вертикальная компонента скорости чувствительного элемента;-до — нормальная сила, тяжести; 5д — ACT; Agsrv — поправка Этвеша. Измерениями в (1) являются оценка V' вертикальной скорости чувствительно элемента, вычисленная по показаниям СНС, и величина вертикальной проекции удельной силы /', действующей на чувствительный элемент. Уравнения измерений вертикальной скорости и удельной силы в частоте СНС можно записать в виде:

f'(h) = f(tk)+p(tk),

V'(tk) = V(tk) + r(tk), (2)

где tk — дискретные моменты времени измерения СНС, tk — tk-\ = At - интервал дискретизации СНС, риг- ошибки измерений.

Члены до и Ад ¡¡ту в (1) с высокой степенью точности определяются по показаниям СНС. ЗАГТ можно сформулировать следующим образом: определить 5g{t) в (1) по измерениям (2). Далее формулируется ЗАГТ в дискретном времени. Используя (2) и известные go(tk), gETv{tk), вычисляются значения z(tk):

z(tk) = ~ f'itk) ~ gait к) - guTvitk),

где V — оператор конечной разности: W'(tk) = V'(tk) — V'{tk-1). Заменив в (1) V на аппроксимацию в дискретном времени V д^, (1), (2) можно переписать в виде:

zitk)=tgitk) + ^fL-p{tk). (з)

В дискретном времени ЗАГТ формулируется следующим образом: но измерениям z{tk), с учетом (3), определить 5g(tk). Из-за наличия шума эта задача математически некорректна без дополнительных предположений о свойствах <jg(tk) и шумов системы.

Далее в диссертации рассматриваются стандартные методы, позволяющие свести ЗАГТ к задаче оптимального стохастического оценивания, основанные на введении для ACT модели стационарного случайного процесса. Рассматриваются различные модели аномалии, применяемые в АГ. Основное внимание уделяется моделям, которые представляются во временной области формирующим уравнением. В частности представлена модель аномалии вида т — го интеграла от белого шума. В работах Ю.В. Болотина, М.Ю. Попеленского показаны сравнительные преимущества использования моделей такого вида2. Отметим, что на практике обычно применяются модели первого и второго порядка.

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

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

В конце этой главы рассматривается вопрос о возможности применения различных адаптивных фильтров для решения ЗАГТ.

Во второй главе представлена разработанная методика адаптивного оценивания силы тяжести с учетом ее пространственной неоднородности. Задача оценивания решалась в стохастической постановке. Путем введения формирующих уравнений для аномалии силы тяжести и погрешностей измерений АГК задача оценивания сводится к задаче оценивания в пространстве состояний. В работе используется следующая система уравнений, связывающая измерения z(tk) и искомое g(tk):

hm6g(tk) = q(tk) \z(tk) = 5g{tk) + ¿Vr(ifc)

Здесь уравнение (3) записано без учета шума гравиметра p(tk), который па практике много меньше шума V . В качестве модели аномалии взято формиру-

-Болотип Ю. В., Попелепский М. Ю. Анализ точпости решения задачи авиагравиметрии при идентификации параметров гравиметра в полете // Фундам. и прикл. матом., 2005. 11. вып. 7. 167-180.

ющее уравнение в дискретном времени вида ш-ого интеграла ог белого шума, где Ут - оператор т-ой конечной разности, ~ гауссовый белый шум в

дискретном времени:

Е №к)} = 0, Е [д(<*М«,)] = Я{Ьк)6к,и

где бщ — символ Кронекера. Дисперсия <5(4) полагается переменной для моделирования пространственной неоднородности. В качестве модели погрешностей СНС в определении скорости используется белый гауссовый шум:

Е [г&)] - 0, Е [г(**)г(«,)] = Щ1к)6ц.

Рассматривается модель шума с переменной дисперсией

Система (4) описывает процесс авторегрессии, измеряемый с шумом. Путем

введения вектора состояния у(1к) = (бд{Ьк). ..., бд^к-т+г), г(4), г(4_ 1)) система (4) сводится к системе в пространстве состояний:

--Ау(гк^) + Ве(1к) = Су{1к)

Здесь А, В, С — матрицы системы. — центрированный двумерный гауссовый шум, составленный из формирующих шумов системы: е(^) = (<?(4)> Г (£*;)) •

ЗАГТ можно рассматривать как задачу оценивания вектора состояния системы (5). В случае, когда <3(4)1 Я{1к) — известные функции времени, задача оптималыюго оценивания вектора состояния у{1к) может быть поставлена как задача оптимального сглаживания - задача минимизации дисперсии ошибки оценки.

В авиационной гравиметрии, как правило, отсутствует точная априорная информация об интенсивности силы тяжести и об уровне ошибок СНС — <3(4)) 11(1 к) неизвестны. Для постановки задачи адаптивного оценивания вводится модель дисперсии шумов системы. Дисперсии <3(4), Я(4) предполагаются кусочно-постоянными функциями, принимающими значения из конечных наборов:

<2(4) е {ЗьСа,...ад е {ДьД2,...,ад.

В каждый момент времени tk дисперсии шумов могут принимать значения из пары {Qi, Rn}. Общее число возможных пар - N = Ari N2 . Для описания дисперсий как функций времени вводится марковская цепь, состоящая из N состояний. Пусть s(tk) € {1,2,..., N} - состояние марковской цепи в момент времени tk, S = [s(io)i s(ti), ■ ■ ■, s(fr-i)] " траектория цепи, где Т - количество отсчетов. Марковская цепь определяется набором переходных аг] и начальных щ вероятностей:

ац = Ps{s{tk) = j|s(ifc_i) = г}, i,j G {1,2,..., N}, щ = Ps{s(to) = г], г G {1,2,..., N}.

Все возможные пары дисперсий {Qi, /?Tt}, нумеруются индексами i € {1,2,..., N} Каждое состояние марковской цепи ассоциируется с соответствующей парой дисперсий. Таким образом, состояние марковской цепи s(tk) задает значения дисперсий шумов в момент tk, траектория цепи S задает дисперсии Q(tk), Щ^к) как функции времени. Представленная комбинация марковской цепи и модели в пространстве состояний является частным случаем модели линейной системы с марковскими скачками. Процесс генерации данных представленной модели выглядит следующим образом:

• реализация в момент tk состояния марковской цепи s(tk) = i',

• состоянию марковской цепи с номером i соответствует наперед заданная пара {Qi,Rn}i, шумы q{tk), r(tk) генерируются с дисперсиями Qi, Rn;

• используя (/(t-k), r(tk) и y(tk-i), формируется вектор состояния y(tk), далее формируется наблюдение z(tk)-

Задача оценивания аномалии на траектории ставится в диссертации как задача оценивания вектора состояния модели по критерию MAB:

y{tk) = argmax/y(y(fib)|Z), (6)

y(tk)

где Z = [z(ta), z(ti),..., z(tT-i)] — набор измерений на галсе. При этом неизвестными являются как некоторые параметры модели 0 = {{Q;, Bn}i, а^. тг, }, так и не наблюдаемая последовательность состояний марковской цепи S. Плотность распределения fy{y{tk)представляет собой модель гауссовой смеси,

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

1. Идентификация.

• Обучение: оценивание всех возможных дисперсий шумов системы, переходных и начальных вероятностей марковской цепи. Критерий оптимизации - метод максимизации правдоподобия (ММП):

в = argmax fz{Z\e).

• Распознавание: оценивание пути марковской цени. Определение дисперсий шумов системы как функций времени. Критерий оптимизации - максимизация апостериорной вероятности (MAB):

S = avgmax fs(S\Z,Q).

2. Фильтрация: оценивание аномалии при известных дисперсиях шумов системы. Критерий оптимизации - MAB:

y(tk) = argmax/y(ij(ifc)iZ,S,e).

VW

Задача идентификации параметров рассматривается в стандартной постановке, принятой для классической СММ.

Особенность данных авиационной гравиметрии позволила свести задачу идентификации к задаче идентификации параметров скрытой марковской модели смеси скользящих средних (СММ-СС). СММ-СС, в отличие от модели линейной системы с марковскими скачками, обладает конечным радиусом корреляции, что позволяет решить задачи обучения и распознавания в оптимальной

3Ghahra,ma.ni Z., Hinton G. Variational Learning for Switching State-Space Models Neural Computation, Vol 12, No. 4., 831-864, April 2000.

постановке. Отметим, что модель СММ-СС и методика решения задачи идентификации параметров СММ-СС были разработаны в ходе работы над диссертацией.

Существенной особенностью данных авиационной гравиметрии является низкое отношение сигнал-шум (ОСШ) — энергия аномалии ög(tk) на несколько порядков меньше энергии шума V д*^. В диссертации разработан алгоритм регуляризации данных авиационной гравиметрии, необходимый для решения задачи идентификации и являющийся ее частью. Алгоритм позволяет сохранить структуру смеси скользящих средних (СС) для регуляризованных данных. Задачи обучения и распознавания решаются для регуляризованных данных.

Задача обучения для СММ-СС ставится как задача ММП. Метод оптимизации — EM-алгоритм (expectation maximization), сходящийся к локальному максимуму функции правдоподобия4. Составной частью EM-алгоритма является алгоритм прямого-обратного хода, разработанный в данной работе на случай процесса вида СММ-СС. Также разработан алгоритм обучения, позволяющий использовать данные со всех галсов съемки.

Задача распознавания для СММ-СС решается путем MAB. Конечный радиус корреляции СММ-СС позволил свести задачу к алгоритму типа Витерби, являющемуся частным случаем алгоритма динамического программирования.

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

Идентификация. Измерения градиента аномалии

Задача идентификации решается, используя измерения ш-ой конечной разности аномалии силы тяжести x(tk) = Vmz(tk), которые представляется в виде смеси СС:

x(tk) = q(tk) + ¿Vm+1r(4). (7)

Здесь q{tu) — тривиальное СС, отвечающее полезному сигналу - аномалии, ^Vm+1r(ifc) — СС, отвечающее шуму СНС. Отметим, что процесс x(£fc) = x(tk)/{VaAt)m, где Va — модуль вектора скорости летательного аппарата, мож-

4Dempster A., Laird N-, Rubin D. Maximum likelihood from incomplete data via the EM algorithm. Royal Statistical Society, Series D, Vol. 39, No 1., 1-38, 1977

но рассматривать как измерения ш-ого пространственного градиента аномалии вдоль траектории. Размерность х(^) — м1_т/с2. В ходе дальнейшего изложения для удобства будем рассматривать процесс (7). Далее в диссертации рассматривается алгоритм регуляризации данных авиационной гравиметрии, который является частью алгоритма идентификации.

Регуляризация данных

Алгоритм регуляризации заключается в нормализации ОСШ путем сведения задачи в область низких частот, где энергия аномалии силы тяжести и шума соизмеримы. В процессе регуляризации данные авиационной гравиметрии сглаживаются оконным фильтром. В результате сглаживания коэффициенты СС (7) сворачиваются с коэффициентами окна, что приводит к существенному увеличению порядков СС. Для устранения избыточности и уменьшения порядков СС данные прореживаются. После прореживания структура смеси СС теряется. Алгоритм регуляризации позволяет сохранить структуру СММ-СС, для чего решается задача аппроксимации модели сглаженных данных на прореженной временной сетке:

Здесь Тк — прореженное время: тк — 1 = пАЬ, х'(тк) — сглаженные и прореженные измерения (7), <1'ь — коэффициенты СС. КЬ' — порядки СС, г' — аппроксимацнонные гауссовые шумы в прореженном времени. Коэффициенты и порядки СС получаются в результате решения задачи ^-аппроксимации. Отметим, что алгоритм аппроксимации позволяет сохранить исходный набор дисперсий {<51, Ип}г для аипроксимационных шумов г/, г'. При дальнейшем рассмотрении (8) штрихи для краткости опускаются.

Модель (8) является частным случаем СММ-СС общего вида, которая представлена далее в диссертации. СММ-СС общего вида состоит из произвольного конечного числа СС. Шумы, соответствующие различным СС, полагаются попарно независимыми. Динамика параметров модели описывается марковской цепыо с конечным числом состояний. В диссертации без потери общности рас-

к'

V

(8)

Обучение

смотрены решения задач обучения и распознавания на примере процесса (8), являющегося смесыо двух СС. Для дальнейшего изложения значение частоты дискретизации в (8) не важно, далее будем пользоваться индексом измерения t G Z:

к L

;=о г=о

Задача обучения, то есть оценивания параметров СММ-СС, поставлена как задача ММП. Неизвестными являются набор параметров В = {а^, 7Г;, Qi, Rj}, траектория марковской цепи S — [s(— р),..., s(0), s(l),..., s(T — 1)], где р — max{if, L}, а s(—p) — начальное состояние марковской цепи. Оценивание происходит при условии известного набора наблюдений X = [ж(0), ж(1),..., х(Т — 1)]:

В = argmax /х(Х|В).

Обучение проводится путем итерационной оптимизации целевой функции ЕМ-алгоритма, доставляющего локальный максимум функции правдоподобия. На каждой тг-ой итерации максимизируется функционал EM-алгоритма по В и проводится замена старых параметров Вп на новые В„+1:

6n+i = argmax (/(В, В„). Функционал EM-алгоритма имеет вид:

i/(e,eil) = ^log(/x,s(X)5|B))Ps{5|X,e„}, (10)

s

где fx,s(X,S|В) — совместная плотность вероятности набора наблюдений X и последовательности состояний S при заданных параметрах В; /^{¿¡ЧХ, 6П} — апостериорная вероятность реализации последовательности состояний S при условии известных параметров, вычисленных на предыдущей итерации. Суммирование в (10) проводится по всем возможным последовательностям состояний.

Введем следующие обозначения для апостериорных вероятностей:

7*(0 = РзМОМ^е»}, (п)

7«(м) = Р3{з(г-1) = ьз(г) = лх,вп}, (12)

т= Р®{3_т|х,еп}. (13)

Ведем обозначения для блоков данных: Х^ = ^(¿1),..., .-г(/2)], = [*(£]), ..., -Используя (11)-(13), целевая функция (10) представлена в следующем виде:

Е/(©,0„) = ^ +

в{-р) т-1

+ Е Е 1оев.(г-1М07«(я(*-1),в(4))+ (14)

Из (14) видно, что оптимизация параметров марковской цепи и параметров СС разделяются. Получены аналитические выражения для переоценивания параметров марковской цепи. Данные выражения по форме совпадают с аналогичными выражениями для переоценивания параметров классической СММ:

п+1 _ Е?=-р+ 1Ъ{ьз) „+1 _ , .

аг,3 ~ ■с-Л-2 ,Л > ~ 1~Р(г)-

Е*=-Р7Нг)

Оптимизации параметров СС в общем виде нелинейная задача. Для ее решения в диссертации применен модифицированный алгоритм покоординатного спуска. Также в диссертации представлен алгоритм обучения, использующий измерения со всех галсов съемки. Для вычисления апостериорных вероятностей (11)-(13) применен разработанный алгоритм прямого-обратного хода.

Алгоритм прямого-обратного хода

Вычисление вероятностей (11)-(13) является одной из главных сложностей при решении задачи обучения модели линейной системы с марковскими скачками. В первую очередь это связано с бесконечным радиусом корреляции данных модели. СММ-СС обладает конечным радиусом корреляции, это позволяет ре-

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

~ /х(х\вп) ■ (15)

Основная идея, используемая в алгоритме прямого-обратного хода для подсчета 7( — представить числитель (15) в виде комбинации частей, зависящих от наблюдений, проведенных до момента времени I, и от наблюдений, проведенных после этого момента. Для этого используются так называемые прямые и обратные вероятности. Прямые вероятности, зависящие от наблюдений до момента определяются следующим образом:

= /хАК 5(-2и-1|вп)- (16)

Обратные вероятности, зависящие от наблюдений после момента í + р, определяются следующим выражением:

ЖЗГ2*"1) = /хАХн?, ^Г'т, ©„). (17)

Используя прямые и обратные вероятности, получено выражение для вычисления (15):

= (18)

1х(Х\вп)

Здесь /х(Х^р-1\Х^,ХЩр-1, 0т) — условная гауссовая плотность, по-

следовательность состояний полностью определяет параметры распре-

деления Х11рР Л Для прямых и обратных вероятностей получены итерационные соотношения:

Начальные значения прямых и обратных вероятностей находятся следующим образом:

р- 1

ар^Б1!.'1) =/х(Хо-1!^1,©»») П а.з{к-1),з(к)1г3{-Р),

к=-Р-1

т-1

А=Г-2р+1

С помощью значений прямых вероятностей, вычисляется значение функции правдоподобия, стоящее в знаменателе (18):

мх|вя) = £ ат-1(^4)-

сТ-1

В рассмотренной форме (16), (17) прямые и обратные вероятности являются не представимыми с помощью машинной арифметики числами. Численно устойчивая модификация алгоритма прямого-обратного хода представлена в отдельном параграфе второй главы диссертации.

Распознавание

Задача распознавания траектории СММ-СС решается методом МАВ:

5 = агёша хР5(5|Х,в). (19)

Здесь § — распознанная последовательность состояний, © — параметры СММ-СС. В качестве © в (19) выбираются параметры, полученные в результате работы алгоритма обучения. Задача сведена к алгоритму типа алгоритма Витерби, который, в свою очередь, является алгоритмом динамического программирования. Оптимизация (19) эквивалентна следующей задаче:

5 = аг§тах1о§/5д(^ Х|0). (20)

Целевая функция (20) представлена в виде:

Г-1

к^/хДХГ1,^1!©) + + log аа{к_1Ык)].

к=р

Алгоритм Витерби — итерационный алгоритм. В данной работе алгоритм реализован в обратном времени. Пусть последовательность состояний Sjz\р фиксирована. Введем обозначения:

т-1

Jt{Stlp) = щьх^^^ЩХЦ, Skk_2p, 0) + logas(k-1]Ak)}. (21)

St k=t

Траектория Sf'1, максимизирующая выражение (21), является оптимальной траекторией, приводящей систему в фиксированную последовательность состояний SfZjp. Для выражения JtiSt l,,) справедливо итерационное соотношение:

MStllp) = max.[loSfx(x(t)\XttZ^Stt.2p,e) + logai(t_i)>i(t) + Jt+iiSU^)].

syt)

Из данного соотношения следует, что оптимальная траектория прихода системы в последовательность состояний Sjl^p может быть найдена по соответствующей траектории, найденной для момента t + 1. Финальная максимизация (19), (20) описана ниже:

maxfjp^1) + log fxAXo'1, I©)]-

Фильтрация

Проведение идентификации позволяет оценить дисперсии Q(tk), R(tk) как функции времени и свести задачу оценивания аномалии силы тяжести к оптимальному сглаживанию. Для сглаживания используется алгоритм BF-сглаживания5, не требующий вспомогательных обращений матриц. Также отметим, что в диссертации разработана методика оценивания градиента аномалии. Задача ставится как задача оценивания полезной компоненты процесса (9) (МАВ) и сводится к фильтру Випера.

Третья глава посвящена проверке разработанной методики адаптивного оценивания. Методика была опробована на данных авиационных съемок района

5Bryson А.Е., Fnizicr М. Smoothing for Linear and Nonlinear Dynamic Systems // Five. Optimum Sys. Synthesis Con}., U. S. Air Forcc Tech. Rept. ASD-TDB, Feb., 1063, pp. 63-119.

реки Омолон6. Данные обладают сильно выраженной пространственной неоднородностью, обусловленной перемежающимися горными хребтами и равнинами в районе съемок. Аномалия силы тяжести моделировалась процессом первого порядка (т = 1). Марковская цепь состояла из двух состояний. Погрешности СНС моделировались стационарным шумом.

Частота съема СНС 10 Гц, скорость самолета (Ан-26) 100 м/с. В процессе регуляризации данные сглаживались фильтром Кайзера с частотой среза ш3 = 0.0237 Гц и частотой подавления сjstop — 0.05 Гц, что эквивалентно длине волны 4.22 км и 2 км соответственно. Далее данные прореживались в п = 100 раз, что эквивалентно пространственному разрешению в 1 км. В результате регуляризации аномалия и шум моделировались СС седьмого порядка.

Алгоритм обучения СММ-СС проводился па всем множестве галсов. Ниже приведены значения дисперсий шумов, полученные в результате обучения:

Qi = 0.079 мГал2, Q2 = 0.443 мГал2, R = 1.36 см2/с2

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

аа = 0.986, ai2 = 0.014, a2i = 0.027, a22 = 0.973,

и значения начальных вероятностей для каждого галса съемки. Алгоритм распознавания позволил определить моменты переключений состояний марковской цепи.

На рис. 1 представлены регуляризованные измерения аномалии, градиента аномалии и соответствующие результаты распознавания для одного галса съемки. На графике выделяются два участки со сравнительно интенсивной аномалией (0—20 км и 120—180 км). Участки интенсивной аномалии обусловлены полетом летательного аппарата над горной местностью. На графике результатов распознавания виден участок с частыми переключениями, которые интерпретировались как ошибки работы алгоритма распознавания.

6Полетные данные были предоставлены компанией ЗАО "ГНПП Аэрогеофизика".

Для борьбы с частыми переключениями происходило сглаживание результатов распознавания, учитывающее результаты с соседних галсов. Сглаживание основано на применении алгоритма нелокальных средних. Алгоритм позволил более четко выделить участки с более интенсивной и менее интенсивной аномалией. Результаты распознавания представлены на рис. 2.

На рис. 3 приведены результаты оценивания аномалии для всех галсов съемки. При построении оценок использовались сглаженные результаты распознавания. Также в диссертации приводится карта аномалии силы тяжести, полу-

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

1-1-1-1-1-1-'-1-1

65.2 65.4 65.6 65.8 66 66.2 66.4 66.6 66.8

1аМис1е

4§ 161 -о 160 -

157-1-1-1-1-1-1-'-'

65.2 65.4 65.6 65.8 66 66.2 66.4 66.6 66.8

1аШис1е

Рис. 2: Карта распознанных состояний (верхний рисунок), сглаженная карта (нижний рисунок). Серым выделено первое состояние (£}1,Я).. черным - второе состояние (<2г, Д)-

о 160

latitude

Рис. 3: Результаты адаптивного оценивания аномалии для всех полетных галсов.

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

На примере исследуемых полетных данных проведено сравнение адаптивного и неадаптивного оценивания аномалии. Ошибка оценивания определялась отдельно для участков с различной интенсивностью аномалии (рис. 2). Под ошибкой понималось ее среднеквадратическое значение (СКО). Для участка с менее интенсивной аномалией и для участка с более интенсивной аномалией проводилась обработка параметрическим семейством неадаптивпых фильтров с заданными постоянными <5, Я. Сравнение проводилось путем варьирования дисперсии формирующего шума аномалии (2 или, что эквивалентно, частотой среза фильтра. За значение дисперсии шума СНС принималось значение, полученное на этапе идентификации. Перед определением точности проводилось согласование галсов.

На рис. 4 представлены результаты оценки ошибки для участка с менее ин-

тенсивпой и более интенсивной аномалией. На дополнительной оси приведены значения длины волны Л, амплитуда которой подавлялась фильтром вдвое. На рисунке также обозначено СКО, соответствующее оптимальным параметрам, полученном в результате идентификации.

8 7 6 ® с

О 5 Е О

I- 4 со

3 2

10 0,1 0,2 0,3 0,4 0,5 0,6 0,7 0,8 0,9

0, тЭа!2

I.............................I..............................[...........................1............................'............................I....................-—|......................1............-~-~1........................|

12 1 10 1 9 2 8,5 8,1 7,7 7,4 7,2 6,9

ХДт

Рис. 4: Оценка ошибки определения аномалии в зависимости от <Э для участка с менее интенсивной аномалией (пунктирная линия) и более интенсивной аномалией (линия с точками).

Полученные результаты позволяют сделать следующие выводы о работе представленного адаптивного алгоритма.

1. Алгоритм позволил получить оценки параметров, близкие к оптимальным.

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

3. На участке с менее интенсивной аномалией получена оценка точности 1.28 мГал. Длина волны аномалии, амплитуда которой подавлена вдвое — 12.77 км. Ожидаемая точность карты аномалии, с учетом осреднения по соседним галсам, не хуже 0.6 мГал. На участке с более интенсивной аномалией

получена оценка точности 2.83 мГал. Длина волны аномалии, амплитуда которой подавлена вдвое — 8.31 км. Ожидаемая точность карты аномалии, с учетом осреднения по соседним галсам, не хуже 1.4 мГал. Данные результаты подтверждаются вычислением невязок на пересекающихся галсах.

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

Задача адаптивного оценивания аномалии силы тяжести на траектории по данным авиагравиметрии поставлена с учетом пространственной неоднородности данных как задача оценивания состояния линейной системы с марковскими скачками.

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

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

Для решения задачи идентификации стандартный для классических скрытых марковских моделей подход, основанный на проведении обучения и распознавания, впервые распространен на случай СММ-СС общего вида. Алгоритм дает максимумы соответствующих функционалов правдоподобия и апостериорной вероятности.

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

• Разработана методика оценки точности определения аномалии силы тяжести на траектории. Показано, что адаптивный алгоритм дает точность на 20% лучше, чем линейный стационарный алгоритм.

Публикации по теме диссертации

1. Дорошин Д.Р. Применение скрытых марковских моделей для адаптивной фильтрации данных аэрогравиметрии // Международная конференция Современные проблемы математики, механики и их приложений. Москва, 2009, С. 277.

2. Болотин Ю.В., Дорошин Д.Р. Адаптивная фильтрация данных авиационной гравиметрии // Симпозиум Международной ассоциации по геодезии (IAG) "Наземная, морская и аэрогравиметрия: измерения на неподвижных и подвижных основаниях"(TG-SMM2010). Санкт-Петербург, 2010, С. 35.

3. Bolotin Yu.V., Doroshin D.R. Adaptive filtering in airborne gravimetry with hidden Markov chains // The 18th World Congress of the International Federation of Automatic Control (IFAC). Milan 2011, pp. 9996-10001.

4. Болотин Ю.В., Дорошин Д.Р. Адаптивная фильтрация данных авиагравиметрии с использованием скрытых марковских моделей // Вестник МГУ, Математика, механика. Москва, 2011, е 3, С. 36-42.

5. Дорошин Д. Р. Методика обработки данных авиационной гравиметрии с учетом пространственной неоднородности гравитационной аномалии // Электронный журнал Труды МАИ, 50, 2012.

В работах, написанных совместно с Ю.В. Болотипым, научному руководителю принадлежит постановка задачи, а основные результаты получены лично автором.

Подписано в печать 9 ноября 2012 г. Формат 60x90/16. Объём 1,25 п.л. Тираж 100 экз. Заказ № 111239

Оттиражировано на ризографе в ООО «УниверПринт» ИНН/КПП 7728572912\772801001 Адрес: 105066, г. Москва, Лефортовский пер., дом 8, корпус 2. Тел. 728-97-17, +7(499)261-78-22. http://www.onlmecopy.ru

Обозначения и сокращения

1 Введение

1.1 Цели работы.

1.2 Авиационная гравиметрия.

1.2.1 Модели гравитационного поля.

1.2.2 Задача построения карт аномального гравитационного поля по данным авиасъемки.

1.2.3 Основные уравнения авиационной гравиметрии.

1.2.4 Обзор методов оценивания в авиационной гравиметрии.

1.3 Обзор методов адаптивной фильтрации.

1.3.1 Модель линейной системы с марковскими скачками.

1.4 Возможности применения методов адаптивной фильтрация для решения задачи авиационной гравиметрии.

2 Адаптивное оценивание силы тяжести

2.1 Задача оценивания силы тяжести как задача оценивания вектора состояния линейной системы с марковскими скачками.

2.2 Идентификация.

2.2.1 Регуляризация данных.

2.2.2 Общий вид скрытой марковской модели смеси скользящих средних (СММ-СС).

2.2.3 Обучение СММ-СС.

2.2.4 Обучение на множестве независимых испытаний.

2.2.5 Обучение. Алгоритм прямого-обратного хода.

2.2.6 Обучение. Оптимизация целевой функции.

2.2.7 Распознавание СММ-СС.

2.3 Оценивание градиента аномалии.

2.4 Оценивание аномалии.

2.5 Определение точности оценивания аномалии на галсе.

2.6 Дополнения.

2.6.1 ЕМ-алгоритм и алгоритм прямого обратного хода в задаче обучения СММ-СС.

2.6.2 Вычислительная устойчивость алгоритма прямого-обратного хода

2.6.3 Аппроксимация смеси скользящих средних в задаче регуляризации данных.

Задачей прикладной гравиметрии является измерение-поля силы тяжести. Детальная информация о поле силы тяжести находит применение в задачах навигации, прогнозирования климата, геодинамики, поиска полезных ископаемых. В геофизических приложениях гравиразведка выделяется в отдельный метод разведочной геофизики.

Силу тяжести в гравиметрии принято представлять в виде суммы нормальной силы тяжести и аномалии силы тяжести (ACT) [77]. Нормальная сила тяжести представляет собой поле модельного референц-эллипсоида вращения, поверхность которого близка к поверхности Земли. Отметим, что под введенным понятием ACT пониматься аномалия в свободном воздухе [77, 76].

В прикладной гравиметрии прибор, измеряющий силу тяжести, называется гравиметром. Статические гравиметры предназначены для проведения измерений на неподвижном основании. Инерциальной гравиметрией называют измерение силы тяжести с борта подвижного носителя. При измерении с борта подвижного носителя используются инерциальные гравиметры. В инерциальной гравиметрии помимо гравиметров используются инерциальная навигационная система (ИНС) и спутниковая навигационная система (СНС), объединенные в систему, называемую аэрогравиметрическим комплексом (АГК). Инерциальная гравиметрия подразделяется на морскую и авиационную. В случае морской гравиметрии АГК устанавливается на борту плавательного средства, в случае авиационной — на борту летательного аппарата. Целью обработки данных авиационной гравиметрии является решение задачи авиационной гравиметрии (ЗАГ) — построение карты ACT района съемок по измерениям АГК. Как правило, аэрогравиметрические съемки состоят из нескольких вылетов летательного аппарата. Траектория каждого вылета летательного аппарата состоит из прямолинейных участков — галсов, и разворотов между галсами. В процессе съемок происходит облет исследуемого участка сетью пересекающихся галсов.

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

Каждая прикладная область предъявляет свои требования к точности и масштабу карт гравитационных аномалий. В настоящее время в авиационной гравиметрии достигнута точность карт около 0.5 мГал (1 мГал= Ю-5 м/с2 ) при разрешении на местности 5 - 10 км [6], [39]. На данный момент целью исследований является доведение точности карт до 0.2 мГал и разрешения до 1 км, что соответствует требованиям геологоразведки для прогнозирования залежей полезных ископаемых [14].

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

1. Обработка данных ИНС и СНС.

2. Определение аномалии на траектории.

3. Построение карты аномалии района съемок на высоте полета.

4. Трансформация карт в соответствии с различными требованиями.

Авиационную гравиметрию принято разделять на скалярную и векторную. Задача скалярной гравиметрии заключается в определении ACT в проекции на вертикаль. В векторной гравиметрии определяются три компоненты ACT. Данная работа посвящена решению задачи скалярной авиационной гравиметрии. Основное внимание уделяется 2-му этапу — задаче авиационной гравиметрии на траектории (ЗАГТ).

ЗАГТ принадлежит классу обратных задач динамики — определение сил, действующих на механическую систему, по траектории системы. Задача относится к классу некорректных обратных задач, и требует регуляризации, основанной на введении дополнительных предположений о структуре ACT и погрешностей датчиков АГК. Одним из методов регуляризации является стохастическая регуляризация. Стандартным подходом здесь является предположение о том, что ACT и погрешности измерений являются стационарными случайными процессами на траектории. Данные предположения позволяют построить стохастическую модель данных авиационной гравиметрии в пространстве состояний, так что ЗАГТ сводится к задаче сглаживания [6], [77].

Предположение о стационарности ACT не всегда адекватно реальности. Для некоторых участков съемок характерна ярко выраженная пространственной неоднородности силы тяжести. Неоднородность силы тяжести в пространстве приводит к нестационарности силы тяжести на траектории летательного аппарата. Подобная структура поля силы тяжести характерна для районов с перемежающейся гористой и равнинной местностью, также неоднородности могут быть вызваны залежами плотных массовых слоев под поверхностью Земли [77]. Пространственная неоднородность притягивающих масс приводит к пространственной неоднородности силы тяжести на высоте полета летательного аппарата. Также для погрешностей датчиков АГК характерна временная неоднородность, выраженная в нестационарности погрешностей в процессе измерений на траектории. Для погрешностей измерений СНС характерно нестационарное поведение, вызванное помехами многолучевости и потерями из видимости спутников [10]. Применение неадаптивного оценивания в данных условиях может привести к существенному ухудшению точности оценивания. В зависимости от выбранных параметров алгоритма оценивания возможны пересглаженность или недосглажснность итоговых оценок ACT [8]. Указанные эффекты приводят к целесообразности применения методов адаптивного оценивания.

В диссертации разработана методика адаптивного оценивания ACT, учитывающая пространственно-временную неоднородность данных авиационной гравиметрии. Разработанная методика позволяет автоматически настраивать характеристики алгоритма оценивания в зависимости от характера измеряемой аномалии и погрешностей датчиков АГК. Далее кратко изложена структура разработанного в диссертации адаптивного алгоритма.

Для данных АГ строится модель в пространстве состояний. Для описания нестационарности используется марковская цепь (МЦ), позволяющая моделировать изменения во времени параметров формирующей линейной модели [7, 41]. Комбинация системы в пространстве состояний и МЦ сводится [16] к модели линейной системы с марковскими скачками (Markov Jump Linear System (MJLS) [48]). Задача оценивания ACT рассматривается как задача оценивания вектора состояния модели MJLS.

Отсутствие априорной информации об аномалии района съемок приводит к необходимости идентификации параметров модели В данной работе, для случая авиационной гравиметрии, задача идентификации MJLS сводится к эквивалентной задаче идентификации для скрытой марковской модели смеси скользящих средних (СММ-СС) [7, 41]. Модель смеси СС описывает измерения градиента аномалии. В свою очередь, задача идентификации СММ-СС редуцируется на задачи обучения и распознавания. Каждая задача является оптимизацией по отдельной группе параметров. Задача обучения - оценивание параметров формирующей системы и параметров МЦ. Задача обучения решается методом максимума правдоподобия (ММП). Оптимизация проводится [41] с помощью ЕМ-алгоритма. Задача распознавания - оценивание траектории МЦ. Распознавание проводится путем максимизации апостериорной вероятности (МАВ). Задача сводится [7, 41] к алгоритму динамического программирования

Важной особенностью данных авиационной гравиметрии является низкое отношение сигнал-шум (ОСШ), что является проблемой для решения задачи идентификации В работе предложена методика регуляризации данных, позволяющая повысить ОСШ за счет сведения задачи в область низких частот [41]. Рассматриваемая методика позволяет сохранить структуру модели измерений градиента аномалии в виде смеси СС. Методика регуляризации данных АГ является частью решения задачи идентификации.

Идентификация позволяет оценить изменение параметров системы во времени, что дает возможность на этапе фильтрации получить оценку ACT с помощью алгоритма нестационарного сглаживания. Отдельно рассмотрена задача оценивания градиента аномалии, как задача оценивания компоненты СММ-СС [7].

Применение предложенной методики рассмотрено на примере обработки полетных данных АГК GT1A [16]. Местность в районе съемок представляла из себя равнинные участки, перемежающиеся горными хребтами и холмистой местностью. На примере полетных данных показана работа алгоритмов регуляризации, обучения, распознавания, нестационарного оценивания. Также приведены сглаженные результаты распознавания, учитывающие результаты распознавания на соседних галсах. Методика сглаживания основана на применении эмпирического алгоритма типа нелокальных средних (NL-means) [44]. Приводится полученная карта аномалии района съемок.

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

Порядок изложения материала следующий: первая глава носит обзорный характер, и содержит информацию, необходимую для постановки ЗАГТ.

• Приводятся модели поля силы тяжести Земли. Рассматривается влияние пространственной неоднородности на точность решения ЗАГ.

• Описана стандартная методика решения ЗАГ, в частности, показано место ЗАГТ в задаче ЗАГ.

• Вводятся уравнения движения и уравнения измерений АГК, использующиеся в данной работе.

• Приводится обзор стандартных методов, позволяющих свести ЗАГТ к задаче оптимального оценивания.

Далее в первой главе представлен обзор методов адаптивной фильтрации. Рассматривается вопрос о возможности применения различных адаптивных фильтров для решения ЗАГТ с учетом пространственно-временной неоднородности данных АГК. Подробно рассмотрена использующаяся в данной работе модель MJLS.

Вторая глава диссертации посвящена постановке и решению ЗАГТ с учетом пространственно-временной неоднородности данных АГК.

• ЗАГТ ставится как задача оценивания вектора состояния модели MJLS. Представлена редукция задачи оценивания на задачи идентификации и фильтрации.

• Вводится модель СММ-СС. Для СММ-СС рассматриваются задачи обучения и распознавания.

• Приводится методика регуляризации данных АГ, являющаяся частью задачи идентификации.

• Задача фильтрации представлена для адаптивного оценивания ACT и для адаптивного оценивания градиента аномалии.

• Представлена методика оценки точности решения ЗАГТ.

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

Применение разработанной методики к полетным данным представлено в третьей главе диссертации.

1.1 Цели работы

Целью работы является постановка и решение задачи скалярной авиационной гравиметрии на траектории с учетом пространственной неоднородности силы тяжести. Разработка данного подхода в первую очередь обусловлена необходимостью повышения точности оценивания для участков с сильно выраженной пространственной неоднородностью силы тяжести. Рассматривается подход к данной задаче, основанный на адаптации параметров аномалии силы тяжести и погрешностей измерений к локальным характеристикам поля силы тяжести и данных АГК.

Основные результаты диссертации

• Поставлена и решена задача адаптивного оценивания аномалии силы тяжести на траектории с учетом пространственно-временной неоднородности данных аэрогра-виметричсского комплекса. Для данных авиационной гравиметрии построена модель в виде линейной системы с марковскими скачками. Задача поставлена как задача оценивания вектора состояния системы.

• Задача оценивания аномалии разделена на задачи идентификации и фильтрации. Удалось свести задачу идентификации к эквивалентной задаче идентификации параметров скрытой марковской модели смеси скользящих средних (СММ-СС). В диссертации решена задача идентификации для СММ-СС общего вида.

• Разработана методика регуляризации, позволяющая нормализовать отношение сигнал-шум данных авиационной гравиметрии. Данная методика включает алгоритм аппроксимации, позволяющий сохранить структуру СММ-СС для регуляри-зованных данных.

• Для идентификации СММ-СС использован стандартный для классических скрытых марковских моделей подход, основанный на проведении обучения и распознавания. СММ-СС обладает конечным радиусом корреляции, что позволило реализовать алгоритмы обучения и распознавания в оптимальной постановке (обучение - метод максимума правдоподобия, распознавание — максимум апостериорной вероятности). Важнейшей частью алгоритма обучения является алгоритм прямого-обратного хода. Задача распознавания сведена к частному случаю алгоритма динамического программирования.

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

• Разработанная методика адаптивного оценивания опробована на реальных данных, обладающих выраженной пространственной неоднородностью. Данные описывались моделью из двух состояний, соответствующих разной интенсивности аномалии. Алгоритм обучения позволил оценить возможные значения формирующих шумов системы, переходные вероятности и начальные вероятности для каждого галса. Проведено сглаживание результатов распознавания, основанное на применении алгоритма нелокальных средних, позволившее выделить пространственные области, соответствующие участкам с различной интенсивностью аномалии. Построены траекторные оценки аномалии и карта аномалии силы силы тяжести района съемок, точность которой от 0 6 до 1.4 мГал.

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