CC BY
28
АЛГОРИТМ КОРРЕКЦИИ ВЛИЯНИЯ ПОГЛОЩЕНИЯ СВЕТА В АТМОСФЕРНЫХ ГАЗАХ НА ГИПЕРСПЕКТРАЛЬНЫЕ ДАННЫЕ ДИСТАНЦИОННОГО ЗОНДИРОВАНИЯ
О.В. Николаева1 1 Институт прикладной математики РАН, Москва, Россия
Аннотация
Предложен экономичный алгоритм коррекции влияния поглощения света в атмосферных газах на коэффициент яркости солнечного света, отраженного системой атмосфера -земная поверхность. Алгоритм не требует предварительного задания оптических параметров аэрозоля, предположений о подстилающей поверхности и общем содержании газов. В алгоритме коррекция выполняется домножением коэффициента яркости на корректирующий множитель, полученный из анализа спектральной зависимости коэффициента яркости; алгоритм применим только к гиперспектральным данным. Приведены результаты тестирования алгоритма на модельных задачах.
Ключевые слова: атмосферная коррекция, коэффициент яркости, влияние газового поглощения.
Цитирование: Николаева, О.В. Алгоритм коррекции влияния поглощения света в атмосферных газах на гиперспектральные данные дистанционного зондирования / О.В. Николаева // Компьютерная оптика. - 2018. - Т. 42, № 2. - С. 328-337. - Б01: 10.18287/2412-61792018-42-2-328-337.
Введение
По данным зондирования Земли из космоса определяются значения коэффициента яркости света, отраженного системой атмосфера - подстилающая поверхность. Задача атмосферной коррекции состоит в исключении из коэффициента яркости атмосферных эффектов (газового поглощения, рэлеевского рассеяния, а также аэрозольного поглощения и рассеяния).
При атмосферной коррекции данных для спектральных каналов, расположенных вне полос существенного газового поглощения, выполняется коррекция только эффектов, возникших из-за взаимодействия света с аэрозолями и при рэлеевском рассеянии [1, 2].
Коррекция влияния газового поглощения выполняется в широко используемых в настоящее время пакетах атмосферной коррекции ATREM (Atmospheric REMoval algorithm), ATCOR (Atmospheric CORrection), FLAASH (Fast Line-of-sight Atmospheric Analysis of Spectral Hypercubes), ISDAS (Imaging Spectrometer Data Analysis System), HATCH (High-accuracy ATmosphere Correction for Hyperspectral data), ACORN (Atmospheric CORrection Now) [3, 4]. К данной группе могут быть отнесены и недавно разработанные пакеты SIERRA (Spectral reflectance Image Extraction from Radiance with Relief and Atmospheric correction) [5], SART (Spherical Atmosphere Radiation Transfer) [6], а также методики, опирающиеся на простые малопараметрические модели коэффициента яркости [7]. Реализованные в этих пакетах алгоритмы определяют альбедо поверхности, полагая известными оптические параметры атмосферного аэрозоля и величину общего содержания (ОС) поглощающих газов. Параметры аэрозоля предварительно определяются по значениям коэффициента яркости в «темных» пикселях (где альбедо поверхности предполагается близким нулю) для длин волн из видимой области спектра, лежащих вне полос газового поглощения. ОС каждого газа находится по относительным изменени-
ям коэффициента яркости в полосе поглощения этого газа.
Большую проблему представляет определение ОС водяного пара, обладающего сильными пространственными вариациями. Если для таких газов, как 02, СО2, СН4, СО, N20, N02 , величина ОС определяется для всего снимка, то ОС водяного пара приходится находить отдельно для каждого пикселя.
Алгоритмы определения ОС газов постоянно совершенствуются. Например, первоначально в пакете АТЫЕМ ОС водяного пара определялось по соотношениям коэффициента яркости в трех спектральных каналах, что приводило к большим погрешностям [8]. В настоящее время используются не только данные во многих спектральных каналах, но и совместное восстановление ОС водяного пара, жидкой воды и льда [9].
В работе [10] предложен критерий качества определения ОС газов. Критерий основан на том факте, что ошибка в определении ОС приводит к большим немо-нотонностям в определяемом альбедо поверхности. Поэтому выбирается такое значение ОС, при котором отклонение восстановленного альбедо поверхности от соответствующего ему отрезка ряда Фурье минимально. В этом случае минимальной будет и ошибка определения альбедо. Использовать дополнительное условие минимальной корреляции альбедо и коэффициентов пропускания атмосферы предложено в [11]. Это условие позволяет существенно повысить точность восстановления альбедо, однако приводит к необходимости решать нелинейную оптимизационную задачу.
Тем не менее сравнение восстановленных альбедо поверхности с наземными эталонами показывает, что рассматриваемые пакеты не всегда позволяют правильно восстановить форму альбедо, см. [12].
Поэтому в [13] предложен другой подход к восстановлению альбедо. Альбедо ищется как линейная комбинация спектральных сигнатур, отвечающих определенным материалам. Коэффициенты этой линейной комбинации, а также коэффициенты коррек-
ции, зависящие от оптических параметров атмосферы, определяются из условия минимального отклонения измеренных и модельных коэффициентов яркости. Здесь также возникает нелинейная оптимизационная задача, для решения которой используется градиентный метод. Информация об аэрозоле и газах в алгоритме не требуется, но необходимы спектральные сигнатуры материалов, входящих в состав поверхности в рассматриваемой местности.
Другой способ использования альбедо, измеренных наземными сканерами, предложен для атмосферной коррекции данных в конкретном регионе океана в коротковолновой области спектра [14]. Альбедо ищется как линейная комбинация наземных альбедо; коэффициенты комбинации находятся из условия наилучшего соответствия наземных альбедо с альбедо, восстановленными для длинноволновой области спектра, где точность восстановления выше.
В настоящей работе предложен алгоритм коррекции влияния газового поглощения на коэффициент яркости. Результатом коррекции является коэффициент яркости для «безгазовой» атмосферы. Далее этот коэффициент может быть подвергнут другой процедуре атмосферной коррекции для получения альбедо поверхности. Также этот коэффициент может быть использован для определения оптических параметров аэрозоля.
Предлагаемый алгоритм опирается на формулу Ламберта - Бэра и находит мультипликативную корректирующую поправку к коэффициенту яркости. Поправка определяется из условия гладкости корректированного коэффициента яркости (здесь мы следуем идеям из работ [10], [11]). Однако, в отличие от [11], предлагаемый алгоритм не требует решения сложной оптимизационной задачи и является явным.
1. Алгоритм коррекции влияния газового поглощения
Коэффициент яркости [15] отраженного в зенит солнечного излучения может быть представлен с помощью известной формулы Ламберта - Бэра [16]
R(X) = e-(a(X)+b(X) °(X)). (1)
Здесь a(X) [см2] - сечение газового поглощения для длины волны X. Величина b(X) = pc(X) есть произведение ОС газа р [мол/см2] и коэффициента воздушной массы c(X); произведение Ь(^)о(^)описывает ослабление в атмосфере излучения за счет газового поглощения. С помощью величины a(X) описывается взаимодействие излучения с аэрозолями и поверхностью, а также процесс рэлеевского рассеяния.
Формула (1) используется для восстановления ОС газов; при этом величины c(X) и a(X) определяются с использованием некоторых предположений о составе аэрозоля и поверхности, см. подробный обзор методов в [16].
Рассмотрим отношение коэффициентов яркости на соседних длинах волн Xj, Xj-1:
-(a(X,)-а(Х,_1)) -(b(Xj)g(Xj )-Ь(Х,_1)o(X,_1))
R(X j )/ R(X j4)=e-
(2)
Поскольку при гиперспектральных измерениях альбедо поверхности и оптические параметры аэрозоля незначительно меняются на соседних длинах волн, особенно по сравнению с сечением поглощения а(А), разность а(А-) - а(А--1) аппроксимируем кубической функцией длины волны
2(Xj) - а(Хj_,) = J us (Xj -Xj-,)1
(3)
Здесь и1, и2, и3 - пока не известные величины. Отметим, что кубическая аппроксимация вводится для разности значений функции а (А) на соседних длинах волн, а не для самой а (А,).
Зависимость от длины волны величины Ь(А) - это зависимость от длины волны коэффициента воздушной массы с(А). Этот коэффициент определяется как отношение пути фотона в слое атмосферы к высоте этого слоя. Для слабо поглощающего газа этот коэффициент будем считать линейной функцией длины волны
Ь(А) = у + у2 А .
При увеличении сечения поглощения возникает зависимость Ь(А) от этой величины. Опираясь на соотношение (48) в Приложении, будем использовать следующее представление
b(X) = J (V2t-, + X) a(k-1)/2(X).
(4)
Здесь предполагается, что коэффициенты разложения в (48) слабо зависят от длины волны. Это предположение справедливо, поскольку эти коэффициенты зависят от сечений рэлеевского рассеяния и взаимодействия излучения с аэрозолем, но не от сечения газового поглощения.
Отметим, что все неизвестные пока коэффициенты их и ук являются размерными величинами. В частности, коэффициенты ук включают неизвестную величину ОС газа р.
Подставляя представления (4) и (з) в равенство (2) и логарифмируя получившееся соотношение, приходим к выражению:
1п [ Я(А --О/Я(А -)] = А(и) + Б- (у) - Б-_х(у).
Здесь вектора и и у содержат неизвестные параметры
и = {ир ^ из} , у = {Ур У2К } ,
а функции А--(и), Б-(у) определяются равенствами
Aj (u) = J us (X j-X j-,)s,
s=1 K
Bj (v) = J (v2t-, + vm Xj) a(k+1)/2(Xj).
(5)
(6)
Коэффициенты U1, U2, из и Vk, k = 1, ...,2K, будем искать как решение задачи
J i A(u) + Bj (v) - Bj-1(v) - ln ^^} ^min, (7)
j=1 I R(Xj) I
s=1
k=1
k=1
где ] - число длин волн в рассматриваемом спектральном диапазоне. Решение задачи (5) - (7) может быть найдено явно методом наименьших квадратов.
Далее по формуле (4) получаем функцию ¿(X). В соответствии с законом Ламберта - Бэра (1) коррекция влияния газового поглощения выполняется дом-ножением коэффициента яркости на мультипликативную поправку
Ь(Х,) а(Х,)
R^H t) = R(A, )eb
(8)
Полученная величина ЕС0ГГ(Х/) - коэффициент яркости для атмосферы, из которой удален поглощающий газ.
Отметим, что предложенный метод коррекции требует знания сечения поглощения а(Х). Это сечение зависит от давления и температуры, то есть меняется с высотой в слое атмосферы.
Чтобы учесть эту зависимость, разобьем слой атмосферы на Ь высотных зон, в 1-ой зоне зададим среднее для нее сечение поглощения (X). Вместо вектора V введем вектор
V = {М1,1, У2,1У2К,1,..., V , У2,ЬУ2К,Ь } .
Тогда формулы (1), (4) - (8) меняются следующим образом
-a(A)- £b, (A) a, (A)
R(A) = e ,
K
bt (A) = £ (vlki + vM, i A) a,k -1)/2(A),
(9) (10)
£ I A(u) + ñj (v) - BBj-i(v) - ln ^^ I ^min, (11) j=i l R(Aj)
K L
Bj(V) = £ £(^-1,i + i Aj)af+1)/2(Aj). (12)
Rcorr(^ j)=Ra j) с (a j),
с a,) = e
£b, (Aj) a, (Aj)
(13)
(14)
Сначала решается задача (11), (12), (5) и определяются величины ык, мщ. Затем находятся величины Ь{ (X) с помощью равенств (10), наконец, осуществляется коррекция по формулам (13), (14).
Отметим, что число неизвестных коэффициентов Ы1, Ы2, Ы3, м>к,и к = 1,...,2К, I = 1,...,Ь, в задаче (11), (12), (5) равно 3+2КЬ. Эта величина должна быть существенно меньше, чем число используемых длин волн /. Таким образом, предложенный алгоритм применим только для данных гиперспектральных измерений.
Выполненную коррекцию можно считать качественной, только если выполнено условие
С(Х,) > 1, , = 1,...,] . (15)
В качестве параметра качества коррекции предлагается использовать величину вариации, определяющей максимальное относительное изменение решения на соседних длинах волн
V = max
2< j <i
Ra
r(A j) - Rc
Rcorr (A.) + Rcorr(A M)
(16)
Выбор вариации V как параметра качества опирается на тот факт, что коэффициент яркости атмосферы, лишенной поглощающего газа, определяется только взаимодействием с аэрозолем, подстилающей поверхностью и рэлеевским рассеянием. Поэтому коэффициент Rcorr(Aj) должен быть медленно меняющейся функцией длины волны. Вариация V в этом случае мала. Если же влияние газового поглощения из коэффициента Rcorr(Aj) удалено не полностью, то величина Rcorr(Aj) содержит сильные немонотонности и величина V велика.
Отметим, что представленный алгоритм не требует задания параметров аэрозоля и ОС поглощающего газа. Алгоритм имеет следующие входные параметры H - высота слоя атмосферы (км), L - число высотных зон, K - порядок разложения в (10), J - количество спектральных каналов, at (A) - сечения поглощения в каждой высотной зоне (см2),
R(Aj) - значения коэффициентов яркости. Выходными параметрами алгоритма являются значения корректированных коэффициентов яркости Rcorr(Aj).
Точность алгоритма коррекции зависит от следующих факторов.
1. Точность полуэмпирических формул (9), (10).
2. Точность задания сечений поглощения a, (A).
3. Корректность разбиения атмосферы на высотные зоны.
В следующем параграфе приведены результаты тестирования алгоритма на модельных задачах.
2. Численные результаты
Рассмотрим спектральный интервал [752, 770] нм, где существенно поглощение света в кислороде. Этому интервалу отвечает J = 41 спектральных каналов прибора типа Maya (http://www.oemoptic.ru/spec maya2000t.php) со спектральным разрешением 0,4 нм.
Рассмотрим слой атмосферы H*= 50 км с морским аэрозолем [17]. Высотные профили температуры и давления, а также концентрации кислорода зададим в соответствии со стандартной моделью атмосферы [18] с шагом 1 км по высоте.
Примем сечения поглощения кислорода как функции температуры и давления, полученные с использованием базы данных HITRAN при помощи специальной программы [19].
Зададим следующие наборы параметров. N0 = 6 зенитных углов солнца: 0 = 5°, 15°, 25°, 45, 55°.
Nt = 5 оптических толщин аэрозоля на длине волны 752 нм:
т = 0,05, 0,175, 0,3, 0,425, 0,55.
L
k =1
k=1 l=1
1=1
N = 10 альбедо поверхности, измеренных на экспериментальных площадках; в набор входит одно альбедо почвы и 9 альбедо травянистой растительности [20].
Для каждого набора ге, 1т, ге = 1,. ,N0, 1т = 1,. ,Дт, = 1,...,Дл, численным решением уравнения переноса находим коэффициенты яркости К;е,;т,;5. Для слоя атмосферы без кислорода находим коэффициенты яркости Я"0 . .
ге .т Л
В алгоритме коррекции используем следующие значения входных параметров.
Высота слоя Н = 10, 20, 30, 40 км.
Число зон Ь = 1, 2, 3, 4.
Порядок разложения К = 1, 2, 3, 4.
Подчеркнем, что корректируемые коэффициенты яркости К;е,;т,и соответствующие им коэффициенты
для «безгазовой» атмосферы Я"0 г находятся для
слоя атмосферы высотой Н*= 50 км. Высота слоя атмосферы Н - параметр алгоритма коррекции -задается независимо от величины Н, причем Н < Н. Приведенные ниже результаты подтверждают возможность такого выбора параметра алгоритма Н.
Погрешность алгоритма определим формулой
Е. . . = тах 1 - Я™. (А.)/Я"° . (А,)| . (17)
ге .. 1<-<ге . .л - I ге. Л - | у 7
Величины (17) суть максимальные относительные погрешности.
Также найдем параметры качества - вариации -по формуле (16)
V г г = тах
10 . . 2<-И
ЯГч (А-) - Я™ л (А--1)
К7< (А-) + Я^г (А--1)
(18)
Результаты тестирования алгоритма коррекции приведены на рис. 1 - 5.
Отметим, что для всех наборов необходимое условие качества коррекции (15) выполнено.
На рис. 1, 2 представлены результаты, полученные алгоритмом коррекции в задаче с параметрами: 0 = 550, т = 0,55, поверхность - поле подсолнечника. Показаны исходный Я, искомый ЯП0 и корректированные ЯС0ГГ коэффициенты яркости. Коррекция выполнена при разбиении слоя атмосферы высоты Н = 30 км на Ь зон при порядке разложения К. Из рис. 1а можно видеть, сколь существенные искажения, вызванные газовым поглощением, корректируются предложенным алгоритмом. На рис. 16 и рис. 2 показано, как корректированный коэффициент ЯС0ГГ сближается с искомым ЯП0 при увеличении параметров Ь и К. Одновременно уменьшаются осцилляции в коэффициенте ЯС0ГГ.
Отметим, что, хотя наилучшее восстановление достигается при К = 4, приемлемыми можно считать и результаты, полученные уже при К = 3, и даже при К = 2 после дополнительного сглаживания величины ЯС0ГГ с помощью полинома низкого порядка.
На рис. 3 приведены погрешности коррекции (17), максимальные по альбедо поверхности, оптическим
толщинам т и углам солнца 0 для различных значений параметров алгоритма Н и Ь при фиксированном значении параметра К. Можно видеть, что при увеличении Н и Ь погрешность алгоритма быстро уменьшается и хорошая точность достигается при достаточно грубом высотном разбиении слоя атмосферы. При этом восстановление ЯП0 с хорошей точностью выполняется и при Н < Н. Это связано с тем, что основная масса кислорода находится в нижних слоях атмосферы.
а)
75Я 762 766 770 длина волны, нм
б)
758Л 762 766 770 длина волны, нм
Рис. 1. Коэффициенты яркости в задаче с параметрами 0 = 55°, т= 0,55, поверхность - подсолнечник. Исходный коэффициент Я, искомый Яп°, корректированные Яс°гг. Параметры алгоритма коррекции: Н = 30 км, К = 4, Ь = 1, 2, 3,4 0,29121
2 Г\1
0,32: о,зо-0,280,260,240,220,20
• ДЯ0(Х)
—2—К-2 —з—К-3 —4—К=4
0,2908
0,2904
а)
758 762 766 770 длина волны, нм
0,2900
б)
758 762 766 770 длина волны, нм
Рис. 2. Коэффициенты яркости в задаче с параметрами в = 55°, т= 0,55, поверхность - подсолнечник. Исходный коэффициент Я, искомый Яп°, корректированные Яс°гг.
Параметры алгоритма коррекции: Н = 30 км, Ь = 4, К = 1, 2, 3,4
На рис. 4 приведены погрешности коррекции (17) максимальные по типам альбедо поверхности, при фиксированных значениях параметров алгоритма Н = 40 км, Ь = 4, К = 4, для различных оптических толщин аэрозоля т и угла солнца 0:
Е(т, 0) = тах е.
1<1, <Д 10
(19)
Из рис. 4 можно видеть, что относительная погрешность восстановления не превышает 10-35 « 0,0003, то есть 0,03%. При этом погрешность слабо зависит от оптической толщины аэрозоля т и зенитного угла солнца 0 (погрешность лишь несколько уменьшается с уменьшением т и 0).
На рис. 5 показана зависимость вариаций Vв•iтi от погрешностей Е^., при фиксированных значениях параметров алгоритма Н = 40 км, Ь = 2, К = 1, 2, 3, 4. Можно видеть, что малым погрешностям Е действительно отвечают малые вариации V, тогда как большие вариации V - большим погрешностям Е. Таким образом, величина вариации, найденная по значениям корректированного коэффициента яркости ЯС0ГГ, может быть использована как параметр качества коррекции.
LoglO(E)
Н, км
2,5
20 ^ШШ'^>^3,0 ь 15^^3,5 10 4,0
Рис. 3. Погрешности коррекции Е, см. (17), максимальные по всем оптическим толщинам аэрозоля Т, альбедо поверхности и зенитным углам солнца в. Параметры алгоритма Н = 10, 20, 30, 40 км, Ь = 1, 2, 3, 4, К = 4
LoglO{E)
Рис. 4. Погрешности коррекции Е, см. (17), для различных оптических толщин аэрозоля Ти зенитных углах солнца в, максимальные по альбедо поверхности. Параметры алгоритма: Н = 40 км, Ь = 4, К = 4
V 0,1-
0,01-.
0,001--
0,0001-
К=1
К=2
К=4
0,0001 0,001 0,01 0,1 1
Рис. 5. Зависимость вариации У1в1Т15 от погрешностей коррекции Е1в1Т^ при различных оптических толщинах аэрозоля Т, зенитных углах солнца в, альбедо поверхности.
Параметры алгоритма Н = 40 км, Ь = 2, к = 1, 2, 3, 4
Наконец, приведем результаты коррекции с помощью предложенного алгоритма данных гиперспектральных измерений, выполненных комплексом Фиалка с борта МКС [21]. Используем измерения, выполненные 11.12.2014 вдоль трассы с начальной точкой ^ 14.38, 8 6.03) и конечной точкой (Е 13.1, N 28.91). Число измерений - 7773. Рассмат-
ривается спектральный интервал [751, 772] нм, в котором содержится 48 спектральных каналов.
На рис. 6 представлены измеренный R и корректированные Rcorr коэффициенты яркости для начальной точки трассы. Коррекция выполнена с параметрами H = 30 км, L = 3, K = 1, 2, 3, 4.
0,20'
750 755 760 765 770
Рис. 6. Результаты обработки измеренных данных в начальной точке трассы. Измеренный коэффициент Я, корректированные Ясогг. Параметры алгоритма коррекции: Н = 30 км, Ь = 4, К = 1, 2, 3,4
Можно видеть, что коррекция исключает из коэффициента яркости минимумы, связанные с газовым поглощением. При увеличении параметра К корректированный коэффициент яркости Ясогг становится все менее осциллирующим. Значение вариации здесь V = 0,009 для К = 4.
На рис. 7 представлено распределение значений вариации V для всех точек трассы. Можно видеть, что для всех точек параметр качества V мал. При этом вариация V < 0,01 для 36 % точек, V < 0,02 для 88 % точек.
доля измерений
0,25-
0,20-
0,15-
0,10-
0,05-
ТТтг^
0 0,01 0,02 0,03 0,04 0,05 0,06
Рис. 7. Результаты обработки измеренных данных. Распределение значений вариации V для всех точек трассы.
Параметры алгоритма коррекции: Н = 30 км, Ь = 4, К = 4
Выводы
Предложенный алгоритм коррекции влияния газового поглощения опирается на представление коэффициента яркости формулой Ламберта - Бэра. Эта формула позволяет выделить из коэффициента яркости сомножитель, отвечающий за ослабление излучения в слое атмосферы за счет газового поглощения.
Коррекция в алгоритме выполняется домножением коэффициента яркости на мультипликативную поправку, обратную этому сомножителю.
Важной частью алгоритма является нахождение данного сомножителя. Здесь основным моментом является использование полуэмпирического выражения для коэффициента яркости (9), (10), опирающегося на формулу Ламберта - Бэра, и разбиение слоя атмосферы на несколько зон. Представление (10) опирается на предположение о слабой зависимости всех оптических параметров атмосферы от длины волны по сравнению со спектральной зависимостью сечения газового поглощения.
Параметры модельного коэффициента яркости (9) находятся из условия его наименьшего отклонения от измеренного коэффициента. При этом используется явный метод наименьших квадратов.
Предложенный алгоритм коррекции обладает следующими свойствами.
1. Является явным и тем самым экономичным.
2. Не требует информации об аэрозоле и поверхности.
3. Не требует значения ОС поглощающего газа.
4. Применим только к гиперспектральным измерениям (поскольку в алгоритме требуется определение большого количества неизвестных параметров в модельном коэффициенте яркости). Отметим, что ОС поглощающего газа, а также оптические параметры аэрозоля в алгоритме не определяются.
Тестирование алгоритма коррекции выполнено для спектрального интервала [752, 770] нм, где существенно поглощение света кислородом. При тестировании использовались модельные задачи с реалистичными параметрами: альбедо поверхности, набором спектральных каналов, моделью стандартной атмосферы для задания зависимости давления и температуры от высоты.
Тестирование показало, что при грубом разбиении атмосферы на 4 высотные зоны алгоритм имеет хорошую точность: максимальная полученная ошибка восстановления не превышает 1 %.
Результаты коррекции данных, полученных комплексом Фиалка с борта МКС, показывают возможность применения представленного алгоритма к реальным измерениям.
Применению алгоритма к результатам реальных измерений должно предшествовать его тестирование на модельных задачах с различными моделями атмосферы отдельно для каждого спектрального диапазона. В тестировании должны быть определены оптимальные значения параметров алгоритма Ь (число высотных зон) и К (порядок разложения в формуле (10)) в каждом диапазоне. Это возможно сделать благодаря гладкой зависимости сечения поглощения от давления и температуры.
Сечения газового поглощения могут быть получены с помощью специальной программы [19].
Заключение
Предложен алгоритм коррекции влияния газового поглощения на коэффициент яркости света, отраженного от системы атмосфера - подстилающая поверхность.
Результатом алгоритма коррекции является коэффициент яркости для атмосферы, из которой удален поглощающий газ. Этот коэффициент может быть использован, например, при определении оптических параметров атмосферного аэрозоля. Кроме того, после дальнейшей коррекции в этом коэффициенте искажений, связанных с взаимодействием света с аэрозолями и рэлеевским рассеянием, может быть найдено альбедо поверхности.
Отметим, что другие алгоритмы выполняют коррекцию всех атмосферных искажений (и газового поглощения, и аэрозольного и рэлеевского рассеяния) одновременно.
Предложен параметр качества восстановления -вариация корректированного коэффициента яркости.
Предметом дальнейший исследований является обобщение алгоритма на случай нескольких поглощающих, а также излучающих газов.
Приложение
Рассмотрим задачу для уравнения переноса излучения в однородном плоском слое высоты H км
цЭ^dz + KtI = SI, (20)
I(0, ц, ф) = I* 5(ц-ц0) 5(ф) при ц> 0, (21)
I(H,ц,ф) = AI при ц<0. (22)
Здесь I(z, ц, ф) - интенсивность излучения в точке z слоя в направлении (ц = cos 0, ф), I* - поток внеатмосферного излучения, ц0 - косинус зенитного угла солнца, к - коэффициент экстинкции. Оператор отражения AI имеет вид
1 1 2п
AI = - J dц' J dф'ц' I(H, ц', ф'), (23)
П 0 0
где A - альбедо поверхности. Здесь предполагается, что отражение от поверхности происходит по закону Ламберта. Интеграл рассеяния SI имеет вид
1 2п
SI = J dц' J dф' I (z, ц', ф') KS (ц, ц', ф - ф'), (24)
-1 0
где к/ц, ц', ф-ф') - индикатриса рассеяния, причем выполнено условие
1 2п
Jdц' J dф' к(ц,ц',ф-ф') = Ks0 <K. (25)
-1 0
Здесь Ks0 коэффициент рассеяния.
Коэффициент экстинкции к является суммой следующих величин
к =Ks0 + ка + к. (26)
Здесь к и к суть коэффициенты поглощения на частицах аэрозоля и в газах соответственно.
Коэффициент яркости R [15] определяется формулой
R = п1 (0, -1,0)/I* ц0. (27)
Эта величина есть отношение интенсивности отраженного в зенит излучения к интенсивности излучения, отраженного белым ортотропным экраном, расположенным горизонтально на верхней границе атмосферы [15].
Теперь подставим (26) в уравнение (20):
цЭ//dz + (kj0 + ка + к) I = SI.
Аналогично имеем уравнение для коэффициента поглощения к' ^ к:
цЭ/73z + (kj0 + ка + к') I'=SI'.
Вычитая одно уравнение из другого, получаем уравнение для разности решений SI=I-I ':
цЭ8^dz + (kj0+ка+к) SI = S 5I - (к-к')I'. (28)
Краевые условия следуют из (21), (22)
SI(0, ц, ф) = 0 при ц > 0 ,
SI (Н, ц, ф) = A SI при ц<0 .
Из принципа максимума [22] следует, что знак функции SI обратен знаку разности к-к', то есть интенсивность излучения является монотонно убывающей функцией коэффициента к
Поэтому коэффициент яркости R также является монотонно убывающей функцией коэффициента газового поглощения к, а также оптической толщины т = кН. Поэтому обратный логарифм lh(1/R) является монотонно возрастающей функцией т.
Оценим зависимость lh(1/R) от величины т. Решение задачи (20) - (26) может быть представлено в виде суммы трех компонент
I (z, ц, ф) = I„ (z, ц, ф) + Il( z, ц, ф) + Iu (z, ц, ф).
Интенсивность нерассеянного излучения I0 является решением следующей задачи
цЭ!0/ dz + KtI0 = 0, (29)
I0 (0, ц, ф) = I* S^ - ц0) S(ф) при ц > 0, (30)
I0(H, ц, ф) = AI0 при ц< 0. (31)
Интенсивность один раз рассеянного излучения I1 определяется задачей
цд^/ dz + KtI1 = SI0, (32)
I1(0, ц, ф) = 0 при ц> 0, (33)
Ij(H, ц, ф)=0 при ц< 0. (34)
Задача для интенсивности многократно рассеянного излучения записывается в виде
^Iu/dz+KtIu = S Iu + S I,, (35)
Iu (0, ц, ф) = 0 при ц> 0, (36)
Iu (Н, ц, ф) = A [ Iu +1 ] при ц< 0. (37)
Решение задачи (29) - (31) находится аналитически
10 = I,8(ц - ц0) 8(ф) е-к-г/ц° при ц > 0 ,
10 = I,А (ц0/л)е-к'н/ц° ек'(н-г)/ц при ц<0. (38)
Подставляя эти выражения в равенство (24) и используя теорему о среднем, получаем правую часть уравнения (32)
£ I _ у е-*1гМ) + z е-^н/ц» е*1 (н-г>/ц<) (39)
где
Y = /. к, (ц, ц,, ф) < I* к,о (max к, (ц, ц, ф) /к,о) ,
ц>0,ф /
ц(z, кг) е (-1,0),
Ц 0 ц
Z = I. А — í йц' í dф'к,(ц,ц',ф-ф') < I. А —к0.
П J J п
-1 0
П
Из равенства (39) следует оценка
SI0 < Smax I0 = Y е-к'zц0 + Z e
-к,Н/ц, e-к, (Н-z)
Используя правую часть £>шах 10 вместо £ 10 в задаче (32)-(34), в соответствии с принципом максимума [22] находим оценку сверху для интенсивности 11:
11 < Уgl( z• ц) + Zg2( z• ц) при ц> 0, 11 < У gз (z• ц) + Zg4(z• ц) при ц< 0.
Здесь
&( z, ц) =
e-к, zj ц0 - e-к, z ц
к,(1 -цЮ
(zl ^)e-**z/ ц0
при ц^ц0,
p4z _р ч
g2( z, ц)=e-^H (1+^)e e
при Ц = Цо,
к, zl ц
g3( z, ц) =
к, (1 + ц)
e-Wцо - e-к^Hцо eк'(Hц
к, (1 + И/ц,) '
к, z „к, (H-zVц+к, H -к,Н (1+1 цо) e e_
g4 (z, ц) =
к, (1 + ц)
(Н - z)e-,c'H(Ш/цо) ezк'
при ц^-1, при ц = -1.
Для функций gi(z, ц) имеем оценки
gt (z, ц) < 1/ кг при i = 1,3,
gt (z, ц) < е-к'H/ц" Iкг при i = 2,4.
Объединяя эти оценки, приходим к соотношению
gi (z, ц) < 1/кг при i = 1,2,3,4.
Отсюда следует оценка
I1 < I* К0/кг(40) где безразмерная величина Z1 не зависит от к z = Aц0/ п+ max к, (ц,ц0,ф) Iк,0 .
ц>0, ф /
Для решения задачи (35) - (37) - интенсивности многократно рассеянного излучения IM - выполнена оценка [22]
1 - A
5 Ч + AlllJI
Здесь используется равномерная норма II II (максимальное значение по всем переменным). Из определения (24) и свойства (25) находим
I* 4 ^„1 ш.
Поэтому имеем соотношение
Im <
1
1 - A
- + A
Отсюда находим оценку для интенсивности рассеянного излучения
ii + im <
1
1 - A
+1
Ii
Используя здесь неравенство (40), получаем
I1 + IM <
1
1 - A
+1
olKt )ZV
Учитывая определение (27) и равенство (26), получаем оценку сверху для коэффициента яркости рассеянного излучения
R1 + RM <
vs0
К+К
■Z 2,
(41)
где безразмерная величина
Z2 =
nZ1
(1 - A) Ц
не зависит от коэффициента газового поглощения к.
В соответствии с (27), (38) имеем выражение для коэффициента яркости нерассеянного излучения
Ro = Z3 е-КHZ4
(42)
где безразмерные величины ^з = Ае (ка+к»о)Н(1+1/^°), Сэ = 1+1/^0 не зависят от коэффициента к.
Объединяя соотношения (41) и (42), приходим к следующим оценкам для коэффициента яркости Я = Я0+Я1+ЯМ
Z3 е
-К H Z4
< R <-
^so
К+К
■Z 2 + Z3 е
-КH Z4
(43)
Перейдем к оптическим толщинам: т,о = к,„Н - рассеяния, Та =каН - аэрозольного поглощения, Т = к Н - газового поглощения. Тогда выражение (43) принимает вид
Т
С3е-Т?4 <Я ^ + £3 е-Т?4 .
Т+Та
Преобразуем выражение справа
Тд0^2
Т + Т„
(44)
-+Z3e"
тZ4 <TsoZ2 +_ZVZ^ <2max(TsoZ2, WZ4) T + Ta V Z4 + T min(Ta Л/ Z4 ) + T '
Тогда оценки (44) принимают вид
-т Z4 — D — Z5
Z3 е < R <
T + Z6
(45)
где безразмерные величины Zs = 2max (т5о^2, Z3 / Z4) и Z6 = min (та, 1 / Z4) не зависят от оптической толщины газового поглощения т.
Из соотношений (45) следуют оценки для обратного логарифма коэффициента яркости ln(1 / R)
1п(т+Z6) - In Z5 < ln(1 R) < т Z4 - In Z3.
(46)
Выше было показано, что величина 1п(1 / Я) является монотонно возрастающей функцией т. В соответствии с оценками (46) величина 1п(1 / Я) может быть аппроксимирована степенной функцией оптической толщины т. Однако, поскольку вариации величины т в полосе газового поглощения велики, точность аппроксимации может быть повышена переходом к разложению по корням >/Т :
1п(1/R) = £ шп (yR)" .
(47)
Здесь оп - безразмерные коэффициенты, не зависящие от величины т.
Перейдем теперь от оптической толщины т к сечению газового поглощения о = т / р, где р - содержание газа в столбе атмосферы (молекул/см2). Из выражения (47) получаем представление коэффициента яркости как функции сечения о
R = exp [ П" (>/ö)"
(48)
Здесь размерные величины n« = ß"P"
не зависят
от сечения о, но зависят от остальных сечений os0, оа, индикатрисы рассеяния os(|, ф - ф'), альбедо поверхности A, высоты слоя H, косинуса зенитного угла солнца |о. Размерность величин п« обеспечивает без-размерность величины в аргументе экспоненты в соотношении (48).
Литература
1. Тарасенков, М.В. Комплекс программ восстановления отражательных свойств земной поверхности в видимом и УФ-диапазонах / М.В. Тарасенков, В.В. Белов // Оптика атмосферы и океана. - 2014.- Т. 27. № 7. - С. 622-626.
2. Shanmugam, P. CAAS: an atmospheric correction algorithm for the remote sensing of complex water / P. Shanmugam // Ann. Geophys. - 2012. - V. 30. - P. 2-3-220.
3. Gao, B.C. Atmospheric correction algorithms for hyperspectral remote sensing data of land and ocean / B.G. Gao, M.J. Montes, C.O. Davis, A.F.H. Goetz // Remote Sensing of Environment. - 2009. - V. 113. - P. S17-S24.
4. Minu, S. Atmospheric Correction Algorithms for Hyperspectral Imageries: A Review / S. Minu, A. Shetty // International Research Journal of Earth Sciences. - 2015. V. 3, №. 5. - P. 14-18.
5. Lenit, X. SIERRA: A new approach to atmospheric and topographic corrections for hyperspectral imagery / X. Lenot, V. Achard, L. Poutier // Remote Sensing of Environment. - 2009. - V. 113. - P. 1664-1677.
1
1
К -К*о
К
s0
К -К*о
К
sO
К -К*о
К
К -Кso
6. Hadjit, H. Atmospheric correction of Earth-observation remote sensing images by Monte Carlo method / H. Hadjit, A. Oukebdane, A.H. Belbachir // J. Earth Syst. Sci. - 2013. -V. 122, №. 5. - P. 1219-1235.
7. Катковский, Л.В. Параметризация уходящего излучения для быстрой атмосферной коррекции гиперспектральных изображений / Л.В. Катковский // Оптика атмосферы и океана. - 2016. Т. 29, № 9. - С. 778-784.
8. Griffin, M.K. Compensation of Hyperspectral Data for Atmospheric Effects / M.K. Griffin, H.K. Burke // Lincoln Laboratory Journal. - 2003. V. 1, №. 1. - P. 29-54.
9. Thompson, D.R. Atmospheric Correction for Global Mapping Spectroscopy: ATREM Advances for the HysplRI Preparatory Campaign / D.R. Thompson, B.C Gao, R.O. Green, D.A. Roberts, P.E. Dennison, S.R.Lundeen // Remote Sens. Environ. - 2015. - V. 167. - P. 64-77/ doi: 10.1016/j.rse.2015.02.010
10. Qu, Z. The High Accuracy Atmospheric Correction for Hy-perspectral Data (HATCH) model / Z. Qu, B. Kindel, A.F.H. Goetz // IEEE Trans. Geosci. Remote Sens. - 2003. -V. 41. P. 1223-1231.
11. Marion, R. Atmospheric correction of hyperspectral data over dark surfaces via simulated annealing / R. Marion , R. Michel, Ch. Faye // Proc. of SPIE. - 2005. - V. 5979. 59791T.
12. Деркачева, А. А. Эффективность атмосферных коррекций гиперспектральных снимков Hyperion в регионах с развитым растительным покровом / А.А. Деркачева, О.В. Туту-балина // Современные проблемы дистанционного зондирования Земли из космоса. - 2014. - Т. 11, № 4. - С. 360-368
13. Денисова, А.Ю. Анализ линейной спектральной смеси, инвариантный к атмосферным искажениям гиперспектральных изображений / А.Ю. Денисова, Ю.Н. Журавель, В.В. Мясников // Компьютерная оптика. -2016. -Т. 40, № 3. - С. 380-387.
14. Буканова, Т.В. Разработка региональных алгоритмов коррекции данных спутниковых сканеров цвета для вод юго-восточной Балтики / Т.В. Буканова, С.В. Вазюля, О.В. Ко-пелевич, В.И. Буренков, С.В. Шеберстов, С.В. Александров
// Современные проблемы дистанционного зондирования земли из космоса. - 2012. - Т. 9, № 4. - С. 70-79.
15. Средства отображения информации коллективного пользования. Требования к визуальному отображению информации и способы измерения. Москва. Стандар-тинформ. - 2008. - 27 с.
16. Rozanov, V.V. Differential optical absorption spectroscopy (DOAS) and air mass factor concept for a multiply scattering vertically inhomogeneous medium: theoretical consideration / V.V. Rozanov, A.V. Rozanov // Atmos. Meas. Tech. - 2010. - V. 3. - P. 751-780.
17. Sayer, A.M. A pure marine aerosol model, for use in remote sensing applications / A.M.Sayer, A. Smirnov, N.C. Hsu, B.N.Holben // Journal of geophysical research. -2012. - V. 117. - D05213. doi:10.1029/2011JD016689.
18. Атмосфера стандартная. Параметры. ГОСТ 4401-81. -2004. 165 с.
19. Герцев, М.Н. Восстановление сечений молекулярного поглощения излучения из базы данных HITRAN / М.Н.Герцев. - 2016. - URL: http://library.keldysh.ru/preprint.asp?id=2016-19 (дата обращения 25.10.2017)
20. Баула, Г.Г. Формирование базы данных гиперспектральных оптических характеристик сельскохозяйственных культур в ультрафиолетовой, видимой и ближней инфракрасной областях спектра / Г.Г. Баула, М.Н. Брычихин, М.И. Истомина, А.Ю. Кротков, Е.Ю. Сжёнов, А.А. Ризванов, В.Н. Третьяков // Космонавтика и ракетостроение. - 2013. - № 4(73). - С. 178-184
21. Ризванов, А.А. Гиперспектральные наблюдения системы атмосфера - земля в ультрафиолетовой, видимой и ближней инфракрасной областях спектра с борта международной космической станции // Космонавтика и ракетостроение. - 2015. - № 6(85). - С. 39-44
22. Локальные свойства решений уравнения переноса. / Т.А. Гермогенова. - М.: Наука. 1986. -272 с.
Сведения об авторе
Николаева Ольга Васильевна, 1971 года рождения, в 1993 году окончила факультет вычислительной математики и кибернетики Московского государственного университета им. М.В. Ломоносова. Старший научный сотрудник Института прикладной математики РАН. Область научных интересов: прямые и обратные задачи для уравнения переноса излучения.
ГРНТИ: 27.35.47
Поступила в редакцию 25 сентября 2017 г. Окончательный вариант - 29 декабря 2017 г.
ALGORITHM FOR ELIMINATING GAS ABSORPTION EFFECTS ON HYPERTISPECTRAL REMOTE SENSING DATA
O. V. Nikolaeva1
1Keldysh Institute of Applied Mathematics RAS, Moscow, Russia
Abstract
An efficient algorithm for the elimination of gas absorption effects on the reflectance of sunlight in an atmosphere-ground system is proposed. The algorithm does not require aerosol, surface and gas concentration information. The corrected reflectance is obtained via the correction factor, which is found via analysis of the spectral dependence of reflectance. The algorithm is applicable only to hyperspectral data. Results of testing on model problems are presented.
Keywords: atmospheric correction, reflectance, gas absorption effects.
Citation: Nikolaeva, OV Algorithm for eliminating gas absorption effects on hyperspectral remote sensing data / O.V. Nikolaeva // Computer optics. - 2018. - T. 42, № 2. - C. 328-337. - DOI: 10.18287/2412-6179-2018-42-2-328-337.
References
[1] Tarasenkov, M.V. Software package for reconstruction of reflective properties of the Earth surface in visible and UV ranges / M.V.Tarasenkov, V.V.Belov // Atmospheric and Oceanic Optic. - 2014.- V. 27. № 7. - P. 622-626.
[2] Shanmugam, P. CAAS: an atmospheric correction algorithm for the remote sensing of complex water / P. Shanmugam // Ann. Geophys. - 2012. - V. 30. - P. 2-3-220.
[3] Gao, B.C. Atmospheric correction algorithms for hyperspec-tral remote sensing data of land and ocean / B.G. Gao, M.J. Montes, C.O. Davis, A.F.H. Goetz // Remote Sensing of Environment. - 2009. - V. 113. - P. S17-S24.
[4] Minu, S. Atmospheric Correction Algorithms for Hyper-spectral Imageries: A Review / S. Minu, A. Shetty // International Research Journal of Earth Sciences. - 2015. V. 3, №. 5. - P. 14-18.
[5] Lenit, X. SIERRA: A new approach to atmospheric and topographic corrections for hyperspectral imagery / X. Lenot, V. Achard, L. Poutier // Remote Sensing of Environment. - 2009. - V. 113. - P. 1664-1677.
[6] Hadjit, H. Atmospheric correction of Earth-observation remote sensing images by Monte Carlo method / H. Hadjit, A. Oukebdane, A.H. Belbachir // J. Earth Syst. Sci. - 2013.
- V. 122, №. 5. - P. 1219-1235.
[7] Katkovsky, L.V. Parameterization of outgoing radiation for quick atmospheric correction of hyperspectral images / L.V. Katkovsky // Atmospheric and Oceanic Optic. - 2016. V. 29, № 9. - P. 778-784.
[8] Griffin, M.K. Compensation of Hyperspectral Data for Atmospheric Effects / M.K. Griffin, H.K. Burke // Lincoln Laboratory Journal. - 2003. V. 1, №. 1. - P. 29-54.
[9] Thompson, D.R. Atmospheric Correction for Global Mapping Spectroscopy: ATREM Advances for the HyspIRI Preparatory Campaign / D.R. Thompson, B.C Gao, R.O. Green, D.A. Roberts, P.E. Dennison, S.R.Lundeen // Remote Sens. Environ. - 2015. - V. 167. - P. 64-77/ doi: 10.1016/j.rse.2015.02.010
[10] Qu, Z. The High Accuracy Atmospheric Correction for Hyperspectral Data (HATCH) model / Z. Qu, B. Kindel, A.F.H. Goetz // IEEE Trans. Geosci. Remote Sens. - 2003.
- V. 41. P. 1223-1231.
[11] Marion, R. Atmospheric correction of hyperspectral data over dark surfaces via simulated annealing / R. Marion , R. Michel, Ch. Faye // Proc. of SPIE. - 2005. - V. 5979. 59791T.
[12] Derkacheva, A.A. The effectiveness of atmospheric correction for Hyperion hyperspectral images in regions with
developed vegetation cover / A.A. Derkacheva, O.V. Tutubalina // Current problems in remote sensing of the earth from space. - 2014. - V. 11, № 4. - P. 360-368
[13] Denisova, A.Y. Estimation of parameters of a linear spectral mixture for hyperspectral images with atmospheric distortions / A.Y. Denisova, Y.N. Juravel, V.V. Myasnikov // Computer Optics. -2016. - V. 40, № 3. - P. 380-387.
[14] Bukanova, T. Development of regional algorithms for the atmospheric correction of satellite oceal color data in the South-Eastern Baltic / T. Bukanova, S. Vazyulya, O. Kopelevich, V. Burenkov, S. Sheberstov, S. Aleksandrov // Current problems in remote sensing of the earth from space. - 2012. - V. 9, № 4. - P. 70-79.
[15] Joint use information displaying means. Requirement to visual representation of information and measurement means. Moscow. Standartinform. - 2008. - 27 p.
[16] Rozanov, V.V. Differential optical absorption spectrosco-py (DOAS) and air mass factor concept for a multiply scattering vertically inhomogeneous medium: theoretical consideration / V.V. Rozanov, A.V. Rozanov // Atmos. Meas. Tech. - 2010. - V. 3. - P. 751-780.
[17] Sayer, A.M. A pure marine aerosol model, for use in remote sensing applications / A.M. Sayer, A. Smirnov, N.C. Hsu, B.N. Holben // Journal of geophysical research.
- 2012. - V. 117. - D05213. doi:10.1029/2011JD016689.
[18] Standard atmosphere. Parameters. GOST 4401-81. - 2004. 165 p.
[19] Gertsev, M.N. Reconstruction of molecular absorption cross-sections of radiation from database HITRAN / M.N.Gertsev. - 2016. - URL: http://library.keldysh.ru/preprint.asp?id=2016-19
[20] Baula, G.G. Development of a Database of Hyperspectral Optical Characteristics of Agricultural Crops in the Ultraviolet, Visible and Near Infrared Regions' Spectrum / G.G. Baula, M.N. Brychikhin, M.I. Istomina, A.Yu. Krot-kov, E.Yu. Szhyonov, A.A. Rizvanov, V.N. Tret'yakov. // Cosmonautics and Rocket Engineering. - 2013. - №4(73).
- P.178-184
[21] Rizvanov, A.A. Hyperspectral Observations of the Atmosphere - Earth's Surface System in the Ultraviolet, Visible and Infrared Spectral Region from the International Space Station / A.A. Rizvanov // Cosmonautics and Rocket Engineering. - 2015. - № 6(85). - P. 39-44
[22] Local properties of solutions on thransport equation / T.A. Germogenova. - M.: Nauka. 1986. -272 p.
Author's information
Olga Vasilievna Nikolaeva, (b. 1971) graduated from Computational Mathematic and Cybernetics department of Moscow State University. She works senior researcher in the Keldish Institute of Applied Mathematics RAS. Her research interests are direct and inverse problems for radiation transport equation.
Received September 25, 2017. The final version - December 29 2017.
