ТЕХНОЛОГИИ ОБРАБОТКИ СПУТНИКОВЫХ ДАННЫХ В ПРОГРАММНОМ КОМПЛЕКСЕ PLANETAMONITORING
Василий Валентинович Асмус
Научно-исследовательский центр космической гидрометеорологии «Планета», 123242, Россия, г. Москва, Большой Предтеченский пер., 7, доктор физико-математических наук, профессор, директор, тел. (8499)255-69-14, e-mail: [email protected]
Алексей Александрович Бучнев
Институт вычислительной математики и математической геофизики (ИВМиМГ) СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Лаврентьева, 6, кандидат технических наук, старший научный сотрудник лаборатории обработки изображений, тел. (383)333-73-32,
e-mail: [email protected]
Владимир Анатольевич Кровотынцев
Научно-исследовательский центр космической гидрометеорологии «Планета», 123242, Россия, г. Москва, Большой Предтеченский пер., 7, кандидат физико-математических наук, заведующий отделом обработки спутниковых данных, тел. (8499)252-03-56, e-mail: [email protected] Валерий Павлович Пяткин
Институт вычислительной математики и математической геофизики (ИВМиМГ) СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Лаврентьева, 6, доктор технических наук, профессор, заведующий лабораторией обработки изображений, тел. (383)333-73-32, e-mail: [email protected]
Геннадий Иосифович Салов
Институт вычислительной математики и математической геофизики (ИВМиМГ) СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Лаврентьева, 6, кандидат технических наук, старший научный сотрудник лаборатории обработки изображений, тел. (383)333-73-
32,
e-mail: [email protected]
Представлены некоторые технологии обработки спутниковых данных в программном комплексе PlanetaMonitoring, разрабатываемом совместно ФГБУ «НИЦ "Планета"» и ИВМиМГ СО РАН. Эти технологии были использованы при решении ряда прикладных задач дистанционного зондирования Земли.
Ключевые слова: дистанционное зондирование, спутниковые данные, программные технологии, предварительная обработка, классификация с обучением, кластерный анализ, статистическое выделение линейных и кольцевых структур, пространственные перемещения природных объектов.
SATELLITE DATA PROCESSING TECHNOLOGIES OF SOFTWARE COMPLEX PLANETAMONITORING
Vasiliy V. Asmus
State Research Center «Planeta», 123242, Russia, Moscow, Bolshoy Predtechensky st. 7, director, tel. (8499)255-69-14, e-mail: [email protected]
Aleksey A. Buchnev
Institute of the Computational mathematics and mathematical geophysics (ICM&MG) SB RAS, 630090, Russia, Novosibirsk, acad. Lavrent’ev av., 6, lab of the images processing senior researcher, tel. (383)333-73-32, e-mail: [email protected]
Vladimir A. Krovotyntsev
State Research Center «Planeta», 123242, Russia, Moscow, Bolshoy Predtechensky st. 7, head of satellite data processing department, tel. (8499)252-03-56, e-mail: [email protected]
Valeriy P. Pyatkin
Institute of the Computational mathematics and mathematical geophysics (ICM&MG) SB RAS, 630090, Russia, Novosibirsk, acad. Lavrent’ev av., 6, lab of the images processing head, tel. (383)333-73-32, e-mail: [email protected]
Gennadiy I. Salov
Institute of the Computational mathematics and mathematical geophysics (ICM&MG) SB RAS, 630090, Russia, Novosibirsk, acad. Lavrent’ev av., 6, lab of the images processing senior researcher, tel. (383)333-73-32, e-mail: [email protected]
Presented some of the technology to process satellite data in the software package PlanetaMonitoring, develop joint SRC “Planeta” and ICM&MG SB RAS. Use of these technologies have been applied to solve a number of remote sensing tasks.
Key words: remote sensing, satellite data, software technologies, preprocessing, supervised classification, clustering, statistical detection of linear and circular structures, spatial drifts of nature objects.
Программный комплекс PlanetaMonitoring, в течение длительного времени разрабатываемый совместными усилиями ФГБУ “НИЦ “Планета” и ИВМиМГ СО РАН, является представительным набором программных технологий, позволяющих решать различные задачи обработки данных дистанционного зондирования Земли (ДЗЗ) [1,2]. В программном комплексе PlanetaMonitoring реализованы следующие программные технологии обработки многоспектральной спутниковой информации оптического, инфракрасного и микроволнового диапазонов:
1. Предварительная обработка и геометрические преобразования.
2. Классификация с обучением.
3. Кластеризация.
4. Обнаружение объектов заданной формы (линеаментов и кольцевых структур).
5. Определение пространственных перемещений природных объектов по многоспектральным разновременным спутниковым изображениям.
Представим эти технологии подробнее.
1. Предварительная обработка спутниковых данных. Одна из программных технологий традиционна - технология предварительной обработки данных ДЗЗ, которая включает в себя набор общепринятых операций по яркостным и контрастным преобразованиям изображений. В группу геометрических преобразований включено масштабирование изображений и отображение космических снимков на растровые карты
(геокодирование). Трансформирование спутниковых изображений в картографическую основу является геометрическим преобразованием, отображающим весь снимок либо его часть на заранее подготовленную растровую географическую карту соответствующей территории. Карта строится в выбранной географической проекции с использованием различных баз данных (в том числе шейп-файлов ГИС ArcView) для нанесения на лист карты контурных элементов пространственных объектов. Отображение снимка на карту основано на использовании идентичных опорных точек снимка и карты. Возможно использование одного из двух типов отображений: на основе применения семейства кусочно-аффинных преобразований, которые строятся на множестве треугольников снимка и карты, получаемых в результате триангуляции выпуклой оболочки множества опорных точек снимка, и на основе отображающего полинома заданной степени (первой, второй или третьей). Мозаики спутниковых изображений формируются путем объединения трансформированных на единую картографическую основу изображений. Трансформированные изображения, имеющие область пересечения, объединяются с использованием интерполяции значений пикселов изображений. Интерполяция проводится с использованием весовых коэффициентов, значения которых зависят от степени удаленности пикселов от краев снимков. В результате такого объединения граница пересечения снимков становится незаметной. Программный комплекс включает ряд основанных на видоизменении гистограмм контрастных преобразований. В частности, реализован механизм референтного выравнивания [2], при котором яркостные и контрастные характеристики изображений приводятся либо к характеристикам некоторого изображения, выбранного в качестве референтного (базового), либо к характеристикам гауссово распределения с заданными параметрами. Этот этап является необходимым в тех случаях, когда интересующий исследователя участок земной территории покрывается несколькими разновременными космическими снимками, полученными с разными условиями съемки. Такая ситуация имеет место в задаче обнаружения кольцевых структур на космических снимках, представляющих импактные кратеры на поверхности Земли. На рис. 1 представлен результат выделения импактного кратера Курай, иллюстрирующий выполнение операций программной технологии предварительной обработки данных ДЗЗ. Были использованы разновременные снимки импактного кратера Курай, полученные со спутника SPOT-4 (разрешение 20 м) 08.10.2007 и 14.05.2012 соответственно. Была синтезирована мозаика выровненных по яркости изображений (снимок от 14.05.2012 принят за базовый) и применен непараметрический статистический алгоритм выделения кольцевых структур на космических изображениях (см. ниже).
2. Классификация с обучением. Центральные вопросы тематической обработки (интерпретации) данных ДЗЗ - вопросы повышения качества дешифрирования - непосредственно связаны с проблемой выбора адекватных алгоритмов распознавания [3]. Современный опыт автоматизированного распознавания данных ДЗЗ показывает, что заранее практически невозможно установить, какой алгоритм будет лучше с точки зрения точности классификации. Поэтому в распознающую систему целесообразно закладывать несколько алгоритмов и выбор оптимального алгоритма проводить эмпирически. В классификации данных ДЗЗ чаще других используются методы, которые можно разбить на две группы: контролируемая классификация (классификация с обучением) и кластеризация (классификация без обучения). В классификации с обучением для распределения векторов признаков по представляющим интерес классам используется процесс обучения классификатора способности различать эти классы на основе репрезентативных выборок представителей классов. Эти выборки называются обучающими. Процесс обучения фактически является процессом построения функций правдоподобия (или решающих функций^
для классов. Эти функции будут затем использоваться для классификации каждого вектора признаков как принадлежащему только одному классу (жесткая классификация) либо нескольким классам (нечеткая классификация). Обучение классификатора может быть проконтролировано путем определения вероятности правильной классификации (для этого часть из обучающих выборок переводится в разряд контрольных). Контролируемая классификация в программном комплексе Р1апе1аМош1:огт§ основана на использовании байесовской стратегии максимального правдоподобия для нормально распределенных векторов признаков. Пусть х - N -мерный вектор признаков х = (х1,...,хЛГ)г , где Ы— число спектральных диапазонов. Предполагается, что векторы х имеют в классе о)г нормальное распределение N(т1, Д) со средним щ и ковариационной матрицей Д . В этом случае байесовская стратегия максимального правдоподобия для поэлементного классификатора формулируется следующим образом [3].
Пусть €1 = (со1,...,сот) - конечное множество классов, р(щ) - априорная
вероятность класса а>1. Тогда решающая функция класса щ имеет вид
ёг(х)= ))-0.51п(| Д |) —0.5{х-тг)Т B:\x-тг).
(1)
Классическое решающее правило для жесткой классификации принимает следующий вид: вектор х заносится в класс со1, если #((х) > £\(х)
для всех / ^ /.
Поскольку физические размеры реально сканируемых пространственных объектов, как правило, больше разрешения съемочных систем, между векторами признаков существуют взаимосвязи [3]. Использование информации подобного рода дает возможность повысить точность классификации, если пытаться распознавать одновременно блок смежных векторов квадратной или крестообразной формы. Будем называть такой блок векторов объектом. Рассмотрим объект х = (х1,...,хь)Т, состоящий из смежных 7У-мерных векторов = (например, в окрестности 3*3,
5*5,... элементов). Решение об отнесении центрального элемента объекта тому или иному классу принимается на основе результата классификации всего объекта. Такой подход порождает целое семейство решающих правил. Во-первых, это использование принципа голосования, т.е. независимая классификация элементов объекта и отнесение центрального элемента к тому классу, которому было отнесено большинство элементов объекта. Во-вторых, это применение текстурных операторов (простейший пример - описание объекта Х через вектор средних составляющих его элементов) с последующим отнесением центрального элемента классу, к которому был отнесен параметр, характеризующий Х. В-третьих, описание объекта Х случайным марковским полем. Приведенный ниже рисунок 2 иллюстрируют работу системы контролируемой классификации при решении задачи оценки состояния водной среды Азовского моря.
Рис. 2. Состояние водной среды Азовского моря
3. Кластеризация. В состав программного комплекса Р1апе1аМош1:огт§ входит реализация классического алгоритма жесткой кластеризации -алгоритма ^-средних, широко используемого для разбиения на кластеры больших объемов многомерных данных [3]. Алгоритм ^-средних может быть отнесен к классу параметрических, т.к. он неявным образом предполагает природу плотности вероятности: кластеры стремятся иметь конкретную геометрическую форму, зависящую от выбранной метрики. Используются следующие метрики: Евклидова, Махаланобиса, Чебышева, с11у-Ь1оск
расстояние. Известно также, что результат кластеризации методом ^-средних зависит от задания начальных центров кластеров. Предоставляется выбор одного из трех вариантов, два из которых определяются на основе статистических характеристик набора данных и один основан на случайной выборке. Другой подход, позволяющий получать разбиение векторов измерений на кластеры произвольной формы, основан на предположении, что исходные данные являются выборкой из многомодового закона распределения, причем векторы, отвечающие отдельной моде, образуют кластер [3]. Таким образом, задача сводится к анализу мод многомерных гистограмм. Альтернативой жесткой разделяющей кластеризации является мягкая или нечеткая кластеризация, разрешающая векторам принадлежать всем кластерам с коэффициентом членства ии е [0,1], определяющим степень
принадлежности у-го вектора /-му кластеру:
с ь
2Х- = 1, V/ и 5Х <1 > V/, (2)
'■=1 м
определяя этими соотношениями нечеткую кластеризацию. Здесь С - число кластеров, Ь - количество векторов измерений. В недавнее время нами в состав системы кластеризации программного комплекса была включена
реализация широко используемого алгоритма нечеткой кластеризации, известного как метод С-средних [4]. Это итерационный алгоритм, который используется для разделения смешанных векторов измерений в данных ДЗЗ. Идея метода заключается в описании сходства вектора с каждым кластером с помощью функции уровней принадлежности, принимающей значения от нуля до единицы. Значения функции, близкие к единице, означают высокую степень сходства вектора с кластером. Очевидно, что сумма значений функции уровней принадлежности для каждого пиксела должна равняться единице. Также, как и в алгоритме ^-средних, параметрами соответствующей процедуры (кроме количества кластеров) являются тип метрики и вариант выбора начальных центров кластеров. Дополнительным параметром является показатель нечеткости, значения которого для данных ДЗЗ предлагается брать близкими к двум. Основная часть работы алгоритма С-средних состоит в итеративном перестроении матрицы уровней принадлежности векторов признаков кластерам и пересчете центров кластеров. Алгоритм заканчивает работу при выполнении заданного числа итераций либо при достижении матрицы уровней принадлежности состояния стабильности, т.е. состояния, при котором норма разности матриц в двух последовательных итерациях не превосходит заданного порога.
4. Обнаружение объектов заданной формы (линеаментов и кольцевых структур). В анализе космических изображений при решении прикладных задач ДЗЗ, в мониторинге природной среды возникает проблема обнаружения протяженных объектов заданной формы, например, линейных и кольцевых структур на случайном фоне. При обработке космических снимков с целью обнаружения на них объектов, представляющих интерес, в силу целого ряда причин предпочтение отдается статистическому подходу [1-3]. Основная причина состоит в том, что вследствие случайного характера природных процессов данные дистанционных измерений (спектросовмещенные изображения) содержат много случайных вариаций, маскирующих различия значений яркости изображения в точках области объекта и в точках области фона. Случайные величины, значения которых получаются в результате измерений (наблюдений) в точках изображения, будем называть наблюдаемыми случайными величинами. В подобной ситуации надежные алгоритмы обнаружения могут быть построены только с помощью вероятностного (статистического) подхода. При этом могут быть получены даже оптимальные в том или ином смысле алгоритмы, если распределения вероятностей значений наблюдаемых величин в точках области объекта и в точках области фона известны заранее. Однако эти распределения на каждом изображении могут быть своими и могут изменяться даже в пределах одного изображения, поэтому на практике, как правило, не известны для наблюдателя. Известно лишь, что они непрерывные. В такой ситуации для обнаружения объектов эффективны так называемые непараметрические тесты (критерии), поскольку распределения статистик (функций от наблюдаемых величин), используемых в этих критериях, не зависят от распределений наблюдаемых величин на
изображениях, когда в поле зрения объекты отсутствуют. Предлагается следующая схема обнаружения объектов. Последовательно (или параллельно) анализируются (почти) все возможные положения объектов, интересующих исследователя. Будем считать, что при исследовании каждого возможного положения объекта все наблюдаемые величины берутся в достаточно отдаленных друг от друга точках на изображении так, чтобы они могли рассматриваться как статистически независимые в совокупности, когда в поле зрения объект отсутствует. В каждом положении, чтобы свести к минимуму риск сделать неверный вывод, следует проверить статистическую гипотезу о том, что случайные величины, наблюдаемые в точках проверяемой области объекта и в близлежащих к ним точках фона, одинаково распределены (однородны), означающую отсутствие объекта. Обнаружение плохо видимых объектов известной наблюдателю формы, в случаях, когда однородность фона по диаметру (или ширине) объекта не имеет места, может быть реализовано путем обнаружения их контуров. Пусть для / = 1 С, и - величины, наблюдаемые в точках
пересечений линии I -й нормали к проверяемому положению контура объекта, т.е. в точке проверяемой области объекта и в близкой к ней точке окружающего объект фона соответственно, и пусть наблюдателю из прошлого опыта известно, что яркость в точке области объекта
стохастически больше или меньше яркости £ в точке фона. Тогда для
проверки гипотезы об однородности можно воспользоваться хорошо известным в математической статистике критерием знаков (в данном случае) разностей дг - с, . Для возможности применения более эффективного критерия следует увеличить число наблюдаемых величин на каждой нормали. Обозначим дг , сп , с,2 величины, наблюдаемые на линии \ -й нормали в этом случае. Согласно результатам работ [5,6] введем в
рассмотрение следующие статистики: = Е^>тах(^4)} ,
¿Г = <пш( с,^с12)\ , У = £ - „V - „V , где / - индикаторная функция
события , равная 1 или 0 в зависимости от того, произошло или не произошло событие {. Тогда для проверки гипотезы однородности (объект отсутствует) наиболее эффективным будет статистический критерий вида | |> 21 (£° где Л(г) - наименьшее целое число Л такое, что
где а - гарантированный уровень значимости критерия.
Так как левая часть здесь возрастает с убыванием Л дискретно, то истинный уровень значимости введенного двухстороннего критерия может получиться значительно меньше а. Кольцевой объект представляет собой нечеткое кольцо, близкое к круглому, т.е. замкнутую нечеткую полосу. Эти объекты можно обнаружить с помощью критерия, подобного вышеприведенному. Наблюдаемые величины следует брать аналогично на
(3)
последовательности нормалей, но теперь уже к средней линии проверяемого положения полосы, причем величины и с,2 берутся в точках по разные стороны от проверяемого положения полосы, а величины дг - в точках, принадлежащих самой области полосы. Статистики же £+ , £ , Я0 и критерий определяются с помощью тех же формул. Пример выделения кольцевой структуры по космическим снимкам представлен на рис.1.
5. Определение пространственных перемещений природных объектов по многоспектральным разновременным спутниковым изображениям. Одним из способов определения пространственных перемещений природных объектов объектов по разновременным спутниковым изображениям является метод, основанный на нахождении максимумов коэффициента взаимной корреляции. Корреляция используется как средство поиска эквивалентов объекта-эталона, представленного в виде изображения м>(х,у) размерами JxK , на изображении /(х,у) размерами МхЫ ; предполагается, что J <М и К < N . Коэффициент взаимной корреляции
Е X +я>у+*)-/„ ]
Г(х,у) = -^----
(Ж-\)о-кС7/
(4)
Здесь wm - среднее значение пикселов в эталоне ^ , /т - среднее значение элементов изображения / в области, покрываемой эталоном. Знаменатель в (4) является произведением стандартного отклонения <тк
пикселов эталона м/ на стандартное отклонение <т/ пикселов изображения / в
области, покрываемой эталоном. Согласно [7] поиск позиций найденных эталонов (определение смещений) на следующем изображении серии может быть реализован одним из трех методов: определением максимума
коэффициента взаимной корреляции в пространственной области, определением максимума коэффициента взаимной корреляции в частотной области на основе быстрого преобразования Фурье и нахождением минимума суммы квадратов расстояний. При этом не предполагается при поиске смещений каких-либо преобразований эталона за исключением преобразования переноса. В представляемой работе определение смещений эталонов производится на основе определения максимума коэффициента взаимной корреляции в пространственной области в соответствии с формулой (4). При этом эталон может подвергаться преобразованию, состоящему из масштабирования, поворота и переноса. Реализован соответствующий алгоритм и получены результаты вычислительных экспериментов на изображениях, полученных с КА «METEOSAT-8» [8]. Эти результаты свидетельствуют как о необходимости учета масштабирования и поворота эталона, так и о приемлемом времени соответствующих вычислений. На рис.3 представлено поле перемещений ледяных полей в Азовском море в период с 5 по 6 февраля 2012 года с использованием
предложенной технологии определения перемещений природных объектов по разновременным многоспектральным спутниковым изображениям.
Рис. 3. Дрейф льда в Азовском море (по данным ИСЗ «Метеор-М» №1,
КМСС)
Работа выполнена частично при финансовой поддержке Российского фонда фундаментальных исследований (проект № 13-07-00068.).
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Дистанционное зондирование: количественный подход /пер. с англ. Пяткина В. П., Юдиной О. А.; под ред. Алексеева А. С. - М.: Недра, 1983. - 415 с.
2. Шовенгердт Р.А. Дистанционное зондирование. Модели и методы обработки изображений. - М.: Техносфера, 2010. - 560 с.
3. Асмус В. В. Программно-аппаратный комплекс обработки спутниковых данных и его применение для задач гидрометеорологии и мониторинга природной среды. Диссертация в виде научного доклада на соискание ученой степени доктора физикоматематических наук. На правах рукописи. М. - 2002. - 75 с.
4. Bezdek J. C. Pattern recognition with fuzzy objective function algorithms. Plenum
Press, New York, 1981.
5. Салов Г. И. О мощности непараметрических критериев для обнаружения протяженных объектов на случайном фоне //Автометрия. - 1997. № 3. - С. 60-75.
6. Салов Г. И. Новый статистический критерий для задач с двумя и тремя выборками, более мощный, чем критерии Вилкоксона и Уитни //Автометрия. - 2011. № 4. - С. 58-70.
7. MSG Meteorological Products Extraction Facility. Algorithm Specification
Document. - Doc. No. EUM/MSG/SPE/022. Issue 2.6. 1 June 2004.
8. Бучнев А. А., Пяткин В. П. Мониторинг облачных образований по данным геостационарных спутников Земли // Автометрия. - 2009. - № 4. - С. 40-47.
© В. В. Асмус, А. А. Бучнев, В. А. Кровотынцев, В. П. Пяткин, Г. И. Салов, 2014