Введение.
Невозможно представить себе современную метеорологию и моделирование климата без спутниковых систем, позволяющих измерять значительное количество параметров атмосферы и океана и имеющих глобальное пространственное покрытие. Спутниковая аппаратура постоянно
ч
совершенствуется, возрастают точность измерений, спектральное и пространственное разрешение.
Для эффективного использования современных многоканальных спутниковых систем необходимы новые методы решения обратных задач обработки полученной информации, учитывающие нелинейные связи между измерениями и восстанавливаемыми параметрами. Поэтому первая цель данной работы это построение высокоточных методов решения обратной задачи восстановления температуры поверхности океана, профилей температуры и влажности атмосферы по спутниковым измерениям в ИК-области спектра. Кроме того, представляет интерес проведение сравнительного анализа различных методов, изучение их чувствительности к ошибкам измерения, обусловленным несовершенством аппаратуры.
Чтобы эффективно решать обратную задачу восстановления характеристик атмосферы и подстилающей поверхности по измерениям современных многоканальных ИК-радиометров, необходима методика оптимального выбора измерительных каналов. В диссертации предлагаются подходы, позволяющие решить эту проблему.
Еще одной целью диссертационной работы является исследование и сведение к минимуму чувствительности решения прямой и обратной задачи для современных ИК-радиометров (например, AIRS) к ошибкам сдвига частоты измерительных каналов. Эта проблема обусловлена возросшим спектральным разрешением современных измерительных систем и, как следствие, значительно возросшей чувствительностью к сдвигу частот измерительных каналов.
В первой главе диссертации рассматриваются основные подходы, позволяющие решать обратные задачи дистанционного зондирования, содержится обзор теоретических работ по решению обратных задач.
Вторая глава содержит физические аспекты решения прямой и обратной задач, здесь приведены теория расчета коэффициентов поглощения, теория решения прямой и обратной задачи, а также рассчитываются линейные операторы, аппроксимирующие прямую задачу. Во второй главе предложена методика выбора измерительных частот.
Третья глава описывает проведенные численные эксперименты по решению обратной задачи восстановления температуры поверхности океана, а также профилей температуры и влажности атмосферы.
Глава 1. Обзор методов решения обратных задач спутниковой
метеорологии.
1.1 Обратные задачи.
Существует два метода проведения измерений: непосредственно (in situ) и дистанционно. Основное достоинство непосредственных измерений в том, что они производятся без каких-либо сложных промежуточных методов обработки данных. Примером прибора, непосредственно измеряющего характеристики атмосферы, является радиозонд - шар, обычно наполненный гелием, измеряющий вдоль траектории своего движения основные параметры атмосферы (температуру, газовый состав, давление). Основная проблема таких измерений возникает, когда измеряемый объект недоступен или труднодоступен, или имеет большие пространственные размеры. Например, верхняя стратосфера труднодоступна для непосредственных измерений. То есть прямые методы обычно имеют жесткие пространственные и временные ограничения.
Этого недостатка практически лишены методы другого типа -дистанционные. Например, спутниковые измерения параметров атмосферы. Основное достоинство дистанционных измерений, это гораздо лучшее временное и пространственное покрытие. К недостаткам можно отнести сложность их проведения и получения конечных параметров среды. Основная проблема таких методов состоит в том, что они не измеряют непосредственно значение физических величин, например, температуры, а измеряют значение другой физической величины, в случае спутников - интенсивности излучения. Эти две физические величины могут быть связаны сложной зависимостью, и восстановление одной величины по другой может быть непростой математической задачей. Обработка дистанционных измерений приводит к классу математических задач, называемому обратными задачами.
Пусть имеются среда и прибор, который измеряет некую физическую величину, характеризующую состояние среды. Тогда задачу определения
показаний прибора по известному состоянию среды называют прямой, а задачу восстановления состояния среды по показаниям прибора - обратной. Именно обратные физические задачи представляют наибольший интерес. Основной проблемой при их решении является их некорректность.
В математической физике хорошо известно, что подавляющее большинство обратных задач являются некорректно поставленными - малым возмущениям исходных данных (измерений) могут соответствовать сколь угодно большие возмущения решения (состояния системы) [5, 8]. Определение корректной задачи по Адамару [15, 16] следующее: задача называется корректно поставленной, если:
1) ее решение существует,
2) решение единственно и
3) решение непрерывно зависит от входных данных, то есть устойчиво по отношению к малым возмущениям (ошибкам) данных наблюдений.
Если хотя бы одно из этих трех условий не выполняется, задача называется некорректно поставленной.
Наиболее часто в случае обратных задач нарушается условие 3, то есть условие устойчивости решения. В этом случае возникает парадоксальная ситуация: несмотря на то, что задача математически поставлена, ее решение невозможно получить обычными методами. Решение, которое испытывает формально бесконечно большие возмущения при малых возмущениях результатов наблюдений, всегда получаемых с некоторой ошибкой, смысла не имеет. Проблема заключается в том, что по существу все задачи обработки и интерпретации данных дистанционного зондирования, астрономических наблюдений и результатов многих физических экспериментов, являются обратными и, как правило, некорректно поставленными.
До появления современных научно обоснованных методов решения обратных задач исследователь либо, используя детальную физическую модель изучаемого явления, сводил обратную задачу к нахождению небольшого числа параметров, либо, основываясь на физической интуиции, отбирал из множества
допустимых решений то, которое лучше всего соответствует здравому смыслу. Однако такие результаты решения обратной задачи можно подвергать критике: в первом случае часто бывает так, что физическая модель, допускающая жесткую параметризацию решения, не отвечает используемым данным наблюдений. Во втором случае выбор решения субъективен, что нехарактерно для научного метода исследований.
Математически под обратной задачей понимается задача отыскания функции z(s) по функции и(х), получаемой из эксперимента или наблюдений, из уравнения вида
u(x) = A[X,z(s)], (1.1.1)
где А есть некоторый оператор, устанавливающий причинно-следственную связь между z(s) и и(х). В уравнении (1.1.1) по наблюдаемым следствиям и(х) процесса мы должны судить о причинах z(s), породивших его.
Во многих случаях обратная задача может быть представлена интегральным уравнением Фредгольма 1-го рода
ъ u(x)=^K(x,s)z(s)ds (1-1.2)
а
где К(х, s) - ядро (непрерывное или квадратично суммируемое по переменным х, s), которое описывает конкретную модель исследуемого процесса.
Математические трудности решения обратных задач связаны с тем, что обратный оператор А'1 , определяемый уравнением (1.1.1), не является непрерывным. Поэтому если данные наблюдений и(х) получены с некоторой ошибкой д (обозначим приближенные данные символом щ(х)), то
соответствующее приближенное решение, полученное стандартным методом,
zs{s) = A-x[us(x)\ (1.1.3)
будет сколь угодно сильно отклоняться от решения, соответствующего идеально точным входным данным и(х).
Предлагаемые ранее методы решения обратных некорректных задач основывались прежде всего на интуиции авторов, и хотя в ряде случаев удавалось получить важную физическую информацию, необходимость в
9
строгой математической постановке и разработке численных методов решения этого важнейшего для современного естествознания круга проблем остро назрела к 60-м годам, особенно в связи с широким внедрением компьютеров в практику научных исследований.
Методы решения некорректных задач получили интенсивное развитие в 60-е годы. Определяющую роль здесь сыграли работы А.Н. Тихонова [17], М.М. Лаврентьева [18], В.К. Иванова [19] и других математиков. Сейчас можно говорить о научных школах, в частности А.Н. Тихонова, которые создали математическую теорию некорректно поставленных задач, разработали методы их решения (регуляризирующие алгоритмы).
Некорректно поставленные задачи рассматриваются как физически недоопределенные. Они плохо поставлены, множества их приближенных решений очень широки, даже неограниченны. Поэтому некорректные задачи нужно доопределить. Для этого необходима дополнительная информация об искомом решении z(s), вытекающая из обширного опыта всесторонних исследований данного процесса. Важно подчеркнуть, что эта дополнительная информация об искомом решении должна быть известна a priori, до решения соответствующей некорректной задачи. Априорная информация позволяет сформулировать критерий отбора приближенного решения из множества приближенных решений уравнения (1.1.1) и построить регуляризирующий алгоритм. Такой информацией могут служить априорные сведения о гладкости искомого решения z(s), его монотонности, выпуклости, неотрицательности, принадлежности к конечно-параметрическому семейству и т. п.
На рисунке 1 приведен пример точного (а) и приближенного (б) решения некорректной задачи - интегрального уравнения Фредгольма 1-го рода (1.1.2) с ядром K(x,s)=V(l+lOO(x-s)2) [8].
Сплошной линией представлено точное решение z(s), которое было задано заранее. Это решение подставлялось под знак интеграла в уравнение (1.1.2) и вычислялась соответствующая ему функция и(х) - идеально точные "входные данные" обратной задачи (1.1.2). Затем в полученную функцию и(х)
10
вносилась погрешность S = 3% от максимального значения и решалась обратная задача: по возмущенной функции щ(х) находилось приближенное решение zs(s). Приближенное решение zs(s) (точки), представленное на рисунке 1а, получено с помощью регуляризирующего алгоритма, использующего априорную информацию о выпуклости искомого решения z(s). При попытке решить эту же задачу без регуляризации (рисунок 16) получаются сколь угодно большие отклонения "приближенного решения" (точки) от истинного. Так проявляется некорректность обратной задачи (1.1.2).
В настоящее время развитая теория решения некорректно поставленных задач успешно применяется для решения многих обратных задач оптики, спектроскопии, астрофизики, оптимального планирования и т.п.
Важно отметить, что регуляризирующие алгоритмы гарантируют сходимость последовательности приближенных решений к точному решению обратной задачи, то есть при стремлении ошибки наблюдений к нулю приближенное решение стремится к точному.
Исследованию задачи восстановления профиля температуры и состава атмосферы посвящено много работ (см., например, [1,4, 6, 7, 9]).
В работе [1] были введены сильные упрощения, что позволило аналитически исследовать задачу восстановления профиля температуры. В следующем разделе приведены результаты этого анализа.
11
О 0.2 0.4 0,6 0.8 1.0 s
Рисунок 1. Результаты решения обратной задачи, описываемой интегральным уравнением (2) с ядром ^x,^2
12
1.2. Решение обратной задачи линеаризацией и обращением оператора.
Наиболее простой подход к решению задачи восстановления профиля температуры это лианиризация уравнения переноса и нахождение прямого решения.
Предположим, согласно Роджерсу [1], что мы имеем измерения сигнала /„, представляющие собой интенсивности излучения атмосферы от уровня 0 до
°°, измеренные на разных частотах v. Путь измерения проводятся в надир, в этом случае отсутствует сигнал, отраженный от поверхности, и пусть известна температура поверхности, что позволяет вычесть из сигнала вклад излучения поверхности. Тогда, для равновесной нерассеивающей атмосферы, измеряемый сигнал может быть представлен в виде:
?л, (1.2.1)
здесь B(y,T{z)) - функция Планка, являющаяся простой функцией температуры и частоты v, r(v,z) - функция пропускания атмосферы от высоты z до бесконечности. Пусть, кроме того, функция пропускания не зависит от температуры и у нас имеется набор измерений на частотах v1,v2,...,vM, расположенных в окреснтости некоторой частоты г(это позволит нам установить однозначное соответствие между функцией Планка B(v,z) и T(z)). Таким образом, у нас имеется Мизмерений:
/,. = IVi = \В(Г,Т(2))К,{2)<Ь i = 1,...,М, (1.2.2)
где Ki(z) = dr(v,z)/dz. Мы получили уравнение линейное по B(v,z), которое теперь будем считать неизвестным. Таким образом, измерения /, это взвешенное среднее от функции Планка с "весовыми функциями" Kt(z).
Задача (1.2.2) недоопределена, поскольку мы пытаемся по конечному числу измерений восстановить непрерывный профиль. Обычный подход в данном случае это разложение B(v,z) по набору М функций:
f, (1.2.3)
13
которые называются "функциями представления". Уравнение (1.2.2), таким образом, запишется:
гдеД,.- элементы квадратной матрицы.
У нас имеется М уравнений и М неизвестных Ьп которые могут быть вычислены после обращения матрицы А. Разрешив уравнение (1.2.4), подставляя bt в (1.2.3), получаем:
где Лу/ ' ji-я компонента обратной матрицы А г. Это уравнение определяет "функции вклада" Д.(2)как вклад в решение конкретного измерения /,.
Для прямого решения, то есть когда сигнал на спутнике точно соответствует вычисленному по формуле (1.2.2), должно выполняться соотношение:
l~Di(z)KJ(z)dz = SiJ, (1.2.6)
где Sy - символ Кроникера.
Если в измерении /, на спутнике присутствует ошибка et, то это
привнесет вклад в ошибку восстановления, который можно найти, используя равенство (1.2.5):
A=*,.A(z). (1.2.7)
Плохая обусловленность задачи будет характеризоваться большими значениями "функции вклада" Д(г), и, следовательно, ошибка р будет
большой.
Чтобы проиллюстрировать подобный рост ошибки, был построен следующий пример. В качестве "весовых функций" Kt(z) был взят набор идеализированных нормированных функций, изображенных на рисунке 2. В качестве "функций представления" W}{z) был выбран набор полиномов по г.
Чувствительность к шуму в измерениях показана графиками функций D,(z)Ha
14
рисунке 3. Значения функции Д(г) особенно велики при z, лежащих вне максимума соответствующей функции K.(z). Компонента ошибки Д. достигает значений в 1ООСЦ.
До этого момента мы выбирали "функции представления" Wj(z)
произвольно. Попытаемся выбрать их таким способом, чтобы чувствительность к шуму (усредненная длинна вектора /?) была минимальной. Потребуем, например, чтобы ^D,2(z) была минимальной. Такой минимум достигается при
Kjiz), (1.2.8)
где (-1) - обозначает обращение матрицы. Минимизация выражения
неявно предполагает, что компоненты ошибки измерения е. в (1.2.7) имеют
одинаковую дисперсию и не коррелируют. Этот результат можно обобщить на случай известной корреляционной матрицы ошибок.
На рисунке 4 видно, как такой выбор "функций представления" Ws(z)
изменил вид "функций вклада" A(z), по сравнению с набором полиномов. Dt(z) уже не достигает таких больших значений, как в предыдущем случае, но чувствительность к шуму, тем не менее, велика и вклад в ошибку в решении достигает 15^,.. Такая чувствительность к шуму в данной обратной задаче это
общая черта всех прямых решений, линейных и нелинейных задач, когда весовые функции значительно перекрываются. Формально решение не может удовлетворять зашумленному измерению, не порождая осцилляции.
В практических задачах прямое решение найти еще тяжелее. Во-первых, излучение, регистрируемое на спутнике, может приходить не только от атмосферы, но и от поверхности Земли, во-вторых, функция пропускания является функцией температуры. В-третьих, измерения прибора на спутнике не монохроматические, а усредненные. Далее, весовые функции Kt(z) имеют
гораздо более сложный вид, чем идеализированные весовые функции на
15
рисунке 2, они перекрываются значительно сильнее, что приводит к # увеличению ошибок.
Таким образом, нет веских оснований пытаться найти прямое решение, и необходимо применять другие методы решения этой обратной задачи.
16
о -I .2 .з : .4 ~.s Normalised Weighting Function
Рисунок 2. Набор идеальных нормализованных весовых функций. Площадь под каждой равна единице.
Рисунок 3. Функции вклада Dt(z) точного решения для функций представления W(z), выбранных набором полиномов по z. На графике
изображена только ближайшая к оси ординат часть графика, значения функций вклада вне максимума соответствующей весовой функции достигают по модулю 1000.
-10 0 10
Contribution Function
Рисунок 4. Функции вклада Dt(z) для точного решения при оптимально выбранных функций представления W} (z).
17 |
