CC BY
96
Вычислительные технологии
Том 14, № 1, 2009
Методология обработки и интерпретации радарных изображений при помощи мультифрактального анализа*
Т.Н. ЧимитдоржиЕв, А. В. Дмитриев Отдел физических проблем при Президиуме БНЦ СО РАН, Улан-Удэ, Россия
e-mail: tchimit@pres.bscnet.ru
В работе обобщаются результаты мультифрактального анализа поляриметрических радиолокационных изображений. Показано, что данный подход позволяет классифицировать земные покровы по степени неоднородности. Продемонстрированы результаты, описывающие характерную геометрическую структуру природных неоднородностей.
Ключевые слова: дистанционное зондирование, обработка изображений, фрактальный анализ, поляризационная сигнатура.
Введение
В последние годы отмечается тенденция увеличения количества радиолокационных систем дистанционного зондирования земной поверхности с высоким пространственным разрешением (1-8 м). Это связано с рядом преимуществ: возможность выполнения съемки в любое время суток и при любом покрытии облачностью, высокая проникающая способность радарного излучения, использование различных состояний эллипса поляризации и т. д. Вместе с тем более сложная технология бокового обзора и когерентное излучение требуют принципиально отличных от оптического диапазона методов и алгоритмов обработки данных. Одно из перспективных направлений исследований в данной области — текстурный анализ снимков. Наряду с традиционным подходом к классификации земных покровов по величине обратного радарного рассеяния, анализ текстуры позволяет получить дополнительную информацию о пространственных неоднородностях.
Развитием текстурного анализа стал фрактальный подход [1-8], основное отличительное свойство которого — самоподобие пространственных флуктуаций на различных масштабах. Данное направление в последние годы значительно расширилось применительно к снимкам оптического [1, 2] и микроволнового [3-8] диапазонов. В работах [2, 3] рассмотрены различные методы расчета фрактальной размерности и области применения теории фракталов. Введено понятие мультифрактальности объектов исследования, получены результаты по обработке снимков различных диапазонов [4-8]. Однако не изучены вопросы поляриметрических особенностей фрактальной размерности, значение
* Работа выполнена при частичной финансовой поддержке Российского фонда фундаментальных исследований (грант № 08-02-98010-р_сибирь_а).
© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2009.
которой может трактоваться как количественная оценка разномасштабных пространственных флуктуаций.
Известно, что большинство неоднородных природных сред обладает несколькими фрактальными размерностями, т. е. так называемой мультифрактальной размерностью. Например, структура деревьев (лесного массива) различна в горизонтальной, вертикальной плоскостях и при промежуточных углах наклона, что было показано в [6, 7] при помощи предлагаемого подхода.
Таким образом, настоящая работа является дальнейшим развитием исследований [6-8] и предполагает создание общей методологии обработки и интерпретации полностью поляриметрических радарных данных при помощи вычислительных технологий фрактального подхода.
1. Описание тестового полигона и экспериментальных данных
В качестве тестового был выбран участок побережья озера Байкал, представительный с точки зрения разнообразия земных покровов и охватывающий лесные массивы, гористую местность, сельскохозяйственные поля, луга, дельту реки Селенга, болота и часть озера.
Для исследования были использованы данные поляриметрического радара с синтезированной апертурой (РСА) SIR-C, выполнявшего съемку в октябре 1994 года. Изображения побережья озера Байкал были получены одновременно в двух частотных диапазонах (L- и C-диапазоны с длинами волн 24 и 5.6 см соответственно). Угол падения в ближней зоне составлял 21.747°, в дальней зоне — 25.869°, в центре сцены составлял 23.923°. Режимом работы РСА была предусмотрена однопроходная съемка Земли для всех четырех комбинаций поляризации радиоволны на излучении/приеме.
Введем следующие общепринятые обозначения: LHH — изображение, полученное в L-диапазоне на горизонтальной согласованной поляризации; CHV — изображение, полученное в C-диапазоне на кроссполяризации и т.д. Отметим, что относительная калибровка данных SIR-C была достаточно корректной, и погрешность составляла ±1 дБ. Данные представлены в формате MLC — multi look complex, для которого характерно когерентное усреднение. Данные дистанционного зондирования также дополнены картографическими материалами различных масштабов. Кроме того, для классификации лесных массивов использовались материалы лесоустройства с детальным описанием таксационных характеристик каждого выдела.
2. Постановка задачи
В работах [3-6] установлено, что по величине локальной фрактальной размерности возможно выделить основные типы земных покровов. Результаты исследования [6] приведены в таблице.
Как можно увидеть по табличным данным, шероховатый неоднородный безлесный участок в C-диапазоне на HH-поляризации имеет размерность, близкую к фрактальной размерности ровных однородных участков, а на HV- и VV-поляризациях размерность составляет 2.62... 2.63, что ближе к объемным неоднородностям. На основе данного факта сделан вывод, что неоднородности, на которых происходит обратное рассеяние,
Вид земного покрова Диапазон и поляризация
ьии ьиу ЬУУ сии сиу суу
Ровные однородные участки 2.02 2.03 2.03 2.03 2.08 2.09
Шероховатые участки с мелким кустарником 2.21 2.25 2.24 2.30 2.62 2.63
Леса с полнотой более 0.7 ... 0.8 и густой крупный кустарник 2.76 2.75 2.75 2.78 2.78 2.79
Леса с полнотой менее 0.7 ... 0.8 с подлеском 2.92 2.90 2.89 2.93 2.91 2.91
имеют в большей степени вертикальную структуру, а горизонтальные размеры этих неоднородностей существенно меньше вертикальных размеров.
Известно, что поляризационные сигнатуры [9] позволяют воссоздать поляриметрические особенности обратного радарного рассеяния при различных углах эллиптичности и наклона плоскости поляризации (рис. 1). По аналогии с поляризационными сигнатурами можно предположить (и частично это описано в [6, 7]), что существует некоторая зависимость локальной фрактальной размерности от поляриметрической комбинации на излучении/приеме. Поскольку локальная фрактальная размерность, полученная по радиолокационным изображениям, характеризует пространственные флуктуации обратного рассеяния, то фактически поляриметрические особенности фрактальной размерности можно трактовать как поляризационную сигнатуру пространственных флуктуаций.
Рассмотрим профиль поляризационной сигнатуры (рис. 1) с углом эллиптичности равным 0о, т. е. линейную поляризацию при различных углах наклона плоскости поляризации.
С этой целью на основе полностью поляриметрических радиолокационных изображений ЯШ-С (НН, НУ, УУ) синтезированы промежуточные поляризации на приеме. При неизменной горизонтальной поляризации на излучении угол наклона линейной поляризации на приеме составил 0 (НН), 20, 40, 60, 80, 90 (НУ), 100, 120, 140, 160о. На рис. 2
Рис. 1. Поляризационная сигнатура лесного массива
Рис. 2. Обратное рассеяние для лесного массива с полнотой менее 0.7... 0.8
показаны зависимости коэффициента ослабления обратного радиолокационного рассеяния от угла наклона плоскости поляризации на приеме при неизменной горизонтальной поляризации на излучении.
Коэффициент ослабления в Ь- и С-диапазоне определялся в окне размером 25 х 25 пикселов для 30 различных участков лесного массива с полнотой древостоя менее 0.7... 0.8. Полнота — таксационный параметр, характеризующей степень плотности размещения деревьев в древостое, т. е. долю использования ими занимаемого пространства. Значение полноты древостоя 0.7... 0.8 выбрано в качестве порогового, поскольку большая часть лесов рассматриваемого полигона имеет полноту менее указанной величины, а практический интерес представляют леса с высокими значениями полноты. На рис. 2 по оси абсцисс отложен угол наклона линейной поляризации Р на приеме от 0о до 180о. По оси ординат отложен коэффициент ослабления К в децибелах. Сплошной линией обозначена зависимость в Ь-диапазоне, штриховой линией — зависимость в С-диапазоне.
По такой же методике были определены аналогичные зависимости для 30 участков леса с полнотой более 0.7.. .0.8. Значения коэффициента ослабления были выше, чем для более редкого леса, в среднем на 3 дБ при всех углах наклона поляризации. Однако характер изменения полученных кривых был идентичен кривым на рис. 2, поэтому данные зависимости не приводятся.
Представленные на рис. 1 и 2 поляризационная сигнатура обратного рассеяния и формы профилей линейной поляризации для леса с полнотой менее 0.7. ..0.8 полностью согласуются с теоретическими и экспериментальными сополяризованными сигнатурами шероховатой поверхности [9]. Однако по данным профилям не представляется возможным оценить пространственную структуру рассматриваемого объекта. Следовательно, необходимо исследовать более сложные явления, сопровождающие обратное рассеяние, — поляриметрические особенности пространственных флуктуаций отраженной волны.
3. Алгоритм расчета
По радиолокационным поляриметрическим изображениям, представленным в градациях серого 0... 255, на первом этапе строится поле локальных фрактальных размерностей (далее — фрактальное изображение). Для создания фрактальных изображений использован метод вариограмм [10, 11], основанный на статистическом моделировании изображений. Предполагается, что изображение может быть аппроксимировано случай-
ным стохастическим процессом — фрактальным броуновским движением, когда среднеквадратичное смещение (дисперсия) диффундирующей частицы зависит от времени не как обычно линейным образом, а аномальным степенным образом [12, 13]:
((X(t) - X(0))2>« t2H,
здесь H — критический индекс аномальной диффузии.
Таким образом, ставится задача — определить критический индекс, он и будет характеристикой фрактальной поверхности. В соответствии с этим методом статистическое соотношение, которое существует для расстояния между двумя пикселами и дисперсией их значений (яркостей), может быть представлено в логарифмическом виде следующим соотношением:
bg([g(x + e) - g(x)]2) = 2Hlog е + log C,
где g — значение яркости пиксела в точке х, е — расстояние между точками, H — критический индекс фрактального броуновского движения, описывающий шероховатость поверхности, C — постоянный коэффициент. Он несуществен для дальнейших вычислений, поскольку фрактальная размерность (критический индекс) вычисляется по тангенсу угла наклона. Тогда фрактальная размерность может быть вычислена следующим образом: D = 3 — H = 3 — B/2, где B — тангенс угла наклона линии регрессии дисперсии g в зависимости от расстояния е.
В данной работе расчеты производились так: на исходном изображении, начиная с его верхнего левого угла, выделяется “окно” размером 25x25 пикселов, по которому определяются минимальное и максимальное расстояния между пикселами на изображении (Rmin, Rmax). Затем вычисляется величина интервала А:
R — R ■
max min
= N ,
где N — число интервалов, в данном случае N = 5. Для всех возможных пар пикселов вычисляются расстояния между ними, и определяется интервал, к которому они принадлежат. Для каждого интервала рассчитывается дисперсия яркостей пикселов, после чего строится зависимость логарифма дисперсий каждого кластера от логарифма верхней границы интервала. Строится линейная регрессия между двумя этими переменными, и значение угла наклона данной линии регрессии затем используется для расчета фрактальной размерности. Вычисленная размерность (как правило, в пределах от D = 2 до D = 3) пересчитывается в значения яркости 0... 255. Полученное значение яркости присваивается первому пикселу первой строки нового изображения. Далее “окно” сдвигается на 1 пиксел вправо по текущей строке, описанная выше процедура повторяется до конца строки и переходит на следующую. Отметим, что приведенные выше размеры “окна” и число интервалов N были выбраны эмпирически путем сопоставления с наземными данными.
4. Анализ поляриметрических сигнатур фрактальных размерностей
Как правило, поляризационные сигнатуры [9] анализируются для всех возможных случаев угла эллиптичности, в данном случае, как и выше, ограничимся одной составля-
ющей — линейной поляризацией волны. По синтезированным поляриметрическим изображениям с углами поляризации на приеме, указанными выше, построены поля локальных фрактальных размерностей.
Изображение Ь-диапазона, синтезированное для линейной поляризации с углом 120° на приеме, представлено на рис. 3, а. Поле локальных фрактальных размерностей, построенное по данному изображению, представлено на рис. 3, б. В градациях серого представлены полученные локальные фрактальные размерности от 2 до 3, темными тонами изображены участки с размерностью, более близкой к 2, увеличение яркости соответствует увеличению размерности.
По той же методике, которая была использована выше для анализа величины обратного рассеяния, получены зависимости, изображенные на рис. 4. По оси абсцисс также отложен угол наклона линейной поляризации на приеме Р в градусах. По оси ординат
б
Рис. 3. Изображение тестового полигона: а — исходное, б — фрактальное
а
Рис. 4. Зависимость фрактальной размерности от угла наклона линейной поляризации на приеме для леса: а — с полнотой древостоя более 0.7... 0.8; б — с полнотой древостоя менее 0.7. ..0.8
отложено значение фрактальной размерности Б в промежутке 2... 3. Сплошной линией обозначена усредненная по 30 реализациям зависимость в Ь-диапазоне, прерывистой линией — усредненная по 30 реализациям зависимость в С-диапазоне.
В С-диапазоне (рис. 4, а) поляриметрическая зависимость фрактальной размерности от угла наклона поляризации на приеме указывает на значительные пространственные флуктуации при углах 50... 130°. Это вызвано рассеянием от сучьев, веток, хвои, листвы, т. е. структурами, которые расположены под подобными углами. Данный факт подтверждается исследованиями [12], где показано, что для С-диапазона характерно обратное рассеяние от густого древостоя, и следует полагать, что флуктуации определяются кроной лесного массива. В то же время обратное рассеяние [12] и полученные пространственные флуктуации обратного рассеяния при иных углах меньше (рис. 4, а), т. е. объемная доля горизонтальных рассеивающих элементов лесной кроны мала.
Фрактальная размерность лесного массива с полнотой более 0.7... 0.8 в Ь-диапазоне близка к 3, т. е. это объемная неоднородность, и, как видно по рис. 4, а, размерность изменяется незначительно. Соответственно, пространственные флуктуации величины обратного рассеяния не зависят от угла наклона поляризации. Данное обстоятельство связано с тем, что рассеяние в дециметровом (Ь) диапазоне определяется всей толщей леса, неоднородности последнего в свою очередь ориентированы равномерно под углами от 0 до 180°.
В случае редкого леса с полнотой менее 0.7... 0.8 (рис. 4, б) наблюдается обратная ситуация: в С-диапазоне размерность достаточно постоянна; для дециметровых волн значения фрактальной размерности на углах 60... 120° ниже, чем на других. Данные результаты связаны с уменьшением количества вертикально ориентированных стволов, что важно для дециметрового диапазона.
Отсутствие существенной угловой зависимости для сантиметрового диапазона указывает только на тот факт, что лес достаточно редкий и рассеяние определяется неоднородностями леса и подстилающей поверхности, размеры которых соизмеримы с длиной волны 5.6 см. Общий (рис. 4) низкий уровень фрактальной размерности в С-диапазоне (небольшие пространственные флуктуации) порядка 2.2 позволяет сделать предположение о диффузном характере обратного радиолокационного рассеяния в сантиметровом диапазоне, соответственно высокий уровень размерности в Ь-диапазоне (значительные пространственные флуктуации) — об уголковом механизме рассеяния (подстилающая
3.0т о
3.0т о
2.8
2.4-
2.6
2.6
2.8
2.4
2.2
2.2
Л град
Л град
2.0
2.0
0 30 60 90 120 150 180
0 30 60 90 120 150 180
а
б
Рис. 5. Зависимость фрактальной размерности от угла наклона плоскости поляризации на приеме для: а — редкого кустарника, б — для лиственного леса и крупного кустарника
поверхность—стволы деревьев) для дециметровых волн.
Полученные результаты согласуются со статистикой наземных наблюдений на тестовом полигоне и показывают:
1) реальные неоднородности верхнего яруса густого хвойного леса ориентированы именно под углами от 50 до 130° от горизонтали;
2) неоднородности всей толщи густого хвойного леса, соизмеримые с длиной волны Ь-диапазона, располагаются равномерно в промежутке углов 0 ... 180° от горизонтали.
На рис. 5 представлены аналогичные зависимости для кустарника и лиственного леса. Высокий уровень пространственных флуктуаций (фрактальная размерность порядка 2.8) в сантиметровом диапазоне указывает на тот факт, что обратное рассеяние возникает большей частью от природных неоднородностей, соизмеримых с длиной волны (листва, мелкие и крупные ветви кустарника и лиственного леса). При этом неоднородности ориентированы равномерно по всем указанным на графике углам наклона.
Что касается вариаций обратного рассеяния от мелкого кустарника в дециметровом диапазоне, то увеличение размерности при углах 40... 140° с локальным минимумом на кроссполяризации объясняется характерной для кустарника ориентацией ветвей под этими углами. Небольшое уменьшение на поперечной поляризации связано с тем, что отсутствуют ярко выраженные вертикально ориентированные стволы. Аналогично трактуется рис. 5, б, с той лишь разницей, что локальный минимум менее выражен вследствие наличия в нижнем и среднем ярусах вертикальных стволов лиственных деревьев, в то время как верхний ярус не имеет ярко выраженной вертикальной ориентации. Большая размерность при других углах свидетельствует о наличии неоднородностей (ветвей) в нижнем ярусе.
Заключение
Полученные при помощи математического аппарата фрактального анализа поляриметрические зависимости пространственных флуктуаций обратного радарного рассеяния на линейной поляризации позволяют лучше разделять леса по полноте древостоя по сравнению с подобной поляриметрической зависимостью усредненной величины обратного рассеяния. Отличия в значениях фрактальных размерностей на различных поляризациях позволяют сделать вывод о преобладающей геометрии рассеяния и соответственно оценить пространственную структуру объекта исследования. Значительными полученными результатами являются:
1) понижение фрактальной размерности в дециметровом диапазоне для редкого леса, что свидетельствует об уменьшении общего количества вертикально ориентированных стволов на единицу площади;
2) увеличение локальной фрактальной размерности в С-диапазоне для леса с полнотой древостоя более 0.7. . . 0.8 при углах наклона поляризации на приеме в пределах 50... 130°, что позволило утверждать о подобном угловом распределении ветвей деревьев верхнего яруса густого леса;
3) локальный минимум на кроссполяризации, при общем увеличении размерности в промежутке углов 40 . . . 140° в дециметровом диапазоне длин волн для кустарника и лиственного леса, указывает на отсутствие ярко выраженных вертикально ориентированных стволов для верхнего яруса лиственного леса и для всего объема кустарника.
Таким образом, составляющая предлагаемой поляризационной сигнатуры локальных фрактальных размерностей, в данном случае — линейная поляризация, позволяет
существенно увеличить возможности классификации и количественной оценки геометрической структуры лесных сообществ. Дальнейшие исследования могут быть направлены на анализ трехмерной поляризационной сигнатуры локальной фрактальной размерности неоднородных природных сред.
Список литературы
[1] Myint S.W. Fractal approaches in texture analysis and classification of remotely sensed data: comparisons with spatial autocorrelation techniques and simple descriptive statistics // Intern. J. of Remote Sensing. 2003. Vol. 24, N 9. P. 1925-1947.
[2] Sun W., Xu G., Gong P., Liang S. Fractal analysis of remotely sensed images: A review of methods and applications // Intern. J. of Remote Sensing. 2006. Vol. 27, N 22. P. 4963-4990.
[3] ПотАпов А.А. Фракталы в дистанционном зондировании // Зарубежная радиоэлектроника. Успехи современной радиоэлектроники. 2000. № 6. С. 3-65.
[4] Иванов В.К., Кучук Г.А., Стадник А.М., Яцевич С.Е. Методы многочастотного радиолокационного зондирования лесов // Зарубежная радиоэлектроника. Успехи современной радиоэлектроники. 2005. № 7. С. 57-72.
[5] Иванов В.К., Пащенко Р.Э., СтАдник А.М., Яцевич С.Е. Применение фрактального анализа при обработке радиолокационных изображений сельскохозяйственных полей // Электромагнитные волны и электронные системы. 2006. Т. 11, № 7. С. 29-36.
[6] ЧимитдоржиЕв Т.Н., Архинчеев В.Е., Дмитриев А.В., Цыдыпов Б.З. Фрактальный анализ радиолокационных поляриметрических данных для классификации земных покровов // Исследование Земли из космоса. 2007. № 4. C. 27-33.
[7] ЧимитдоржиЕв Т.Н., Архинчеев В.Е., Дмитриев А.В. Поляриметрическая оценка пространственных флуктуаций радарных изображений для восстановления структуры лесного полога // Исследование Земли из космоса. 2007. № 5. C. 80-82.
[8] ЧимитдоржиЕв Т.Н., Архинчеев В.Е., Дмитриев А.В., Цыдыпов Б.З. Поляриметрическая оценка пространственных флуктуаций радиолокационной фазы для классификации земных покровов // Исследование Земли из космоса. 2008. № 1. С. 24-30.
[9] Van Zyl J.J., Zebker H.A., Elachi C. Imaging radar polarization signatures: theory and observation // Radio Science. 1987. Vol. 22, N 4. P. 529-543.
[10] Mark D.M., Aronson P.B. Scale-Dependent fractal dimensions of topographic surfaces: An empirical investigation with applications in geomorphology and computer mapping // Math. Geology. 1984. Vol. 16, N 7. P. 671-683.
[11] Lam N.S.-N., Qiu H.-L, Quattrochi D. An evaluation of fractal surface measurement methods using ICAMS (Image Characterization and Modeling System) // ACSM/ASPRS Annual Convention. Seattle. 1997. Vol. 5. P. 377-386.
[12] Архинчеев В.Е. Случайные блуждания на иерархических фрактальных (гребешковых) структурах // ЖЭТФ. 1999. Т. 115, № 4. С. 1285-1294.
[13] Arkhincheev V.E. Anomalous diffusion and charge relaxation on comb model: exact solutions // Physica A. 2000. Vol. 273. P. 204-215.
Поступила в редакцию 6 марта 2008 г.
