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

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

ГУСЕВ Владимир Юрьевич

МЕТОДЫ И СРЕДСТВА РАДИОМЕТРИЧЕСКОЙ И ГЕОМЕТРИЧЕСКОЙ ОБРАБОТКИ СКАНОВЫХ ИЗОБРАЖЕНИЙ ЗЕМНОЙ ПОВЕРХНОСТИ

05.13.11 - математическое и программное обеспечение вычислительных машин,

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

АВТОРЕФЕРАТ диссертации на соискание учёной степени кандидата физико-математических наук

2/НОЯ 2014

Москва - 2014

005555751

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

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

доцент Юров H.H.

Научный консультант: кандидат физико-математических наук,

доцент Зайцев В.Е.

Официальные оппоненты: доктор физико-математических наук, профессор, зав. сектором ИПМ РАН им. Келдыша М.В. Соколов С.М., доктор физико-математических наук, профессор, профессор МИЭТ (НИУ) Умняшкин C.B.

Государственный научно-исследовательский институт авиационных систем Защита состоится «19» декабря 2014 г. в 12 ч. 00 мин.

на заседании диссертационного совета Д 212.125.04 в Московском авиационном институте по адресу: 125993, г. Москва, Волоколамское шоссе, д. 4. Предварительный заказ пропусков по телефону: 8-499-158-40-90. С диссертацией можно ознакомиться в библиотеке МАИ или по ссылке http://goo.gl/j5yxo4.

Отзывы на автореферат, заверенные гербовой печатью организации, просьба направлять по указанному адресу. Автореферат разослан «// » 2014 г.

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

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

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

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

Задачи обработки изображений, в том числе полученных при космической съемке, решались в таких научных центрах как: Институт систем обработки изображений РАН, Самарский аэрокосмический университет, ГосНИИАС, НИИ ТП, РГРТУ, ИПМ РАН, МГУ. Свой вклад в развитие методов геометрической и радиометрической коррекции спутниковых изображений внесли Желтов С.Ю., Лобанов А.Н., Журкин И.Г., Злобин В.К., Еремеев В.В, Бугаевский JI.M., Цветков В .Я., Rosenfeld A., Pratt W., Kronberg Р., Heiava U. Большой вклад в обработку изображений внесли: Баяковский Ю.М., Галактионов В.А., Сергеев В.В., Сойфер В.А., Соколов С.М., Умняшкин C.B., Конушин A.C., Лемпицкий B.C. и др. Методы обработки спутниковых изображений реализованы в отечественных и зарубежных системах: OrthoNormScan, Фотоплан, Талка, Erdas, ArcGIS, ENVI и др. Однако, некоторые актуальные классы изображений и востребованные практикой виды обработки не имеют удовлетворительных решений.

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

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

3

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

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

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

• разработка методов удаления структурных искажений - помех в виде вертикальных полос, шумовых горизонтальных полосок;

• разработка высокоточных методов совмещения изображений, полученных соседними ПЗС-датчиками;

• оценка качества разработанных алгоритмов с помощью различных мер сравнения изображений и специально разработанных программных средств;

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

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

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

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

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

Автором реализовано и внедрено системное и прикладное программное обеспечение для предложенных методов.

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

• ОКР "Прибор-ОЭК-ИСЗ 5" по заказу ОАО "ЛОМО".

• ОКР "НПО-АСР" по заказу ОАО «Красногорский завод им. СЛ. Зверева». Соответствие диссертации паспорту специальности. Работа выполнена в

областях, предусмотренных паспортами специальностей: 05.13.11 область 7 (модели, методы, алгоритмы и программные средства машинной графики, визуализации, обработки изображений) и 05.13.01 область 12 (визуализация, трансформация и анализ информации на основе компьютерных методов).

Апробация работы. Основные результаты диссертации докладывались и обсуждались на следующих научных конференциях и семинарах:

• "Отечественные разработки в области геодезии и картографии и их применение в хозяйственной и оборонной деятельности страны", доклад "К вопросу фильтрации структурных искажений на космических оптико-электронных изображениях", Москва, ЦНИИГАиК, 1 марта 2012 г.;

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

• "Ломоносов-2013", доклад "Геометрическая коррекция спутниковых изображений", Москва, МГУ, 11 апреля 2013 г.;

• НТС кафедры вычислительной математики и программирования МАИ, 13 декабря 2013 г.;

• Семинар лаборатории компьютерного машинного зрения ГосНИИАС под руководством д.ф.-м.н. Визильтера Ю.В., 11 июня 2014 г.;

• Семинар ИПМ РАН под руководством д.ф.-м.н. Платонова А.К., д.ф.-м.н. Павловского В.Е., д.т.н. Пряничникова В.Е., 18 сентября 2014 г.;

• "ГрафиКон'2014", статья "Комбинированный метод коррекции последствий неравномерности чувствительности датчиков спутниковой съемочной системы сканового типа", Ростов-на-Дону, 30 сентября - 3 октября 2014 г. Представлено в секцию обработки изображений и видео.

Публикации. Материалы диссертации опубликованы в 10 научных работах: шесть статей в пяти разных рецензируемых изданиях, рекомендованных ВАК [1-6], и четыре — в материалах конференций [7-10].

Структура и объем диссертации Диссертация состоит из введения, четырёх глав, заключения, библиографии и приложений. Объем диссертации составляет 157 страниц. Объем приложений - 26 страниц. Работа содержит 94 рисунка и 19 таблиц. Библиография включает 87 наименований.

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

В первой главе рассматриваются методы удаления вертикальных полос. Изучаются особенности съемочной аппаратуры сканового типа и процесса съемки, анализируются возникающие специфические помехи. Так, вертикальные полосы возникают из-за деградации характеристик пикселей регистрирующей аппаратуры. Производится анализ существующих методов, в том числе известного метода на основе анализа модов гистограмм и фрагментарного анализа изображения*', выявляются их недостатки. Сложность задачи в том, что корректное восстановление яркости по всему снимку требует учета характера самого изображения, на котором могут встречаться значительные перепады яркостей, например, облака и контуры природных и искусственных объектов. Многие известные методы не универсальны, либо некорректно обрабатывают изображения с существенными неоднородностями по текстуре и яркости.

•) Злобин В.К., Еремеев В.В. Обработка азрокосмических изображений. - М.: ФИЗМАТЛИТ, 2006. - 288 с.

6

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

g(x,y) =а(хУЛху)Щх), (1)

где fix у) - изображение без помех, g(xy) - изображение с искажениями, а(х), Ь(х) -коэффициенты, зависящие от номера столбца. Задача заключается в получении приближенного неискаженного изображения /(х, у).

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

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

РтЬ>\0) = Api(v\s) + ( 1 -А)рс(у;Л (2)

где вектор параметров в = [s; сг2; A], v - разность яркостей, a pL и рс следующие функции распределения вероятностей Лапласа и Гаусса:

1 _М 1 - "2

Р,=~е *,Рс=-£^Гае (3)

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

В каждый момент обработки рассматриваются два очередных соседних столбца р2 и р,. Осуществляется выравнивание р2 к яркостям столбца р,, то есть предполагается, что пиксели р, имеют эталонный уровень яркости, а пиксели из р2 линейно искажены. Далее строится двумерная гистограмма Нп яркостей для р2 и р,. Значение IIп в некоторой позиции (v,, v2) соответствует числу соседних пар пикселей, для которых яркость 1-ого пикселя Vj, а 2-ого v2.

Для уровней яркости q = 0, Vmax — 1 пикселей первого столбца, где Vmax - число градаций яркости, находится соответствующий уровень яркости пикселей второго столбца m{q) = argmax[W12((/,))]. Задача состоит в определении коэффициентов а2 и

b2 линейного преобразования яркостей. Имея массив соответствий яркостей q и m(q),

*) David Odoma and Peyman Milanfar. Modeling Multiscale Differential Pixel Statistics. Electrical Engineering Department, University of California, Santa Cruz CA. 95064 USA. 7

применяя модель линейной регрессии и метод наименьших квадратов (МНК), находим данные коэффициенты.

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

аобщ'-~ ai • о общ, Ьобиа2 • Ьобщ + Ь2, (4)

для предотвращения распространения и усиления ошибок важно как можно более точное вычисление коэффициентов.

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

В случае метода с сегментацией, по мере продвижения по изображению, группы столбцов пикселей шириной в размер блока подвергаются блочной сегментации. Сначала вертикальный фрагмент разбивается на блоки и производится вычисление характеристик блоков, таких как МО и среднеквадратическое отклонение (СКО). МО и СКО считается для каждого столбца блока. Общее СКО блока вычисляется как среднее СКО по ширине блока, чтобы избежать влияния вертикальных полос. Далее производится слияние соседних блоков в регионы по критериям подобия МО и малого значения СКО. Для корректировки яркостей регионы сортируются по размеру. Коэффициент ¿>2 вычисляется как средняя разность по пикселям, которые принадлежат наибольшим регионам, и площадь которых больше порогового значения. Если таких областей не набирается, то Ь2 вычисляется по всем пикселям.

Второй предлагаемый метод учитывает локальные характеристики яркости по столбцам. Для р2(у) — PiOO ПРИ каждом значении у вычисляются МО М(у) и СКО D(y) по окну Д пикселей. Далее производится сортировка массива M с соответствующим переупорядочиванием D. Крайние значения из M и D удаляются, пока диапазон различий в M не будет меньше порогового значения. Далее выбираются элементы M(i), для которых D(i) меньше установленного порога. Коэффициент Ь2 вычисляется как

8

среднее отобранных значений из М, если их число превышает пороговое значение, иначе Ь2 вычисляется по всем пикселям.

Предотвратить накопление ошибок из-за последовательной обработки по столбцам можно, выполняя корректировку только в случае превышения определенного порога разности между столбцами, тем самым снижая частоту перевычисления аобщ, Ъ0вщ. Кроме этого, предложен этап постобработки. Производится повторная более грубая корректировка g(x,y) с удалением только явно выделяющихся полос (к примеру, порог перепада средних значений яркости больше 10), получаем изображение fix,у). По значениям £f=1/(z,i) для всех х строится массив значений А ( скользящего среднего с довольно большим окном. Задание большого окна позволяет учесть на изображении общую информацию о яркости и сгладить влияние оставшихся полос. Далее подобное скользящее среднее Л2 строится по fix,у), и затем значения fix,у) обновляются по формуле fix,у) ■= fix,у) + At (х) - Аг О).

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

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

Рис. 1. Исходное изображение с полосами и результат коррекции.

Также предлагается глобальный метод коррекции. Идея метода - в задании специальной функции энергии изображения, минимизация которой приведет к удалению полос. Энергетическая функция определяется следующим образом:

Е(Ь) = S^xZS (fix + 1,у)~ fix,y)f ■ Sxy =

£?=1 irx:Kaix +1.У) - 9ix,у)+ £(* + i)- S(*))2 • (5)

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

Для каждых двух столбцов х и х + 1 корректируемого изображения производится

деление гистограмм на классы. Разбиение гистограммы H для каждого столбца

начинается с определения порога ts деления на классы методом Оцу. Данное разбиение

принимается, если выполняются все следующие условия:

'тп! > m • т2 > m-k-i

max(N1) -vt> maх^) • кг (6)

тах(Н2) - vt > шах(Я2) ■ кг где Н\ и #2 - полученные части гистограммы, т1 = (t), m2 = Si W2 (i), m = Si H СО» vt — усредненное по окрестности значение гистограммы в точке разбиения ¿2 - установленные коэффициенты (0 - 1). Рекомендуется к\ = 0.1, к2 = 0.5. Данные условия говорят о том, что оба класса не должны быть слишком малы, и максимумы в данных классах должны возвышаться над усредненным значением между классами. Усреднение необходимо, так как гистограмма может иметь форму гребенки. Если условия не выполняются, учитывается вся гистограмма, а если выполняются, берется часть, соответствующая большему значению из тг и т2. Границы класса могут быть дополнительно сужены путем отбрасывания тех участков, отклонение которых от МО класса больше к3 СКО класса (к3 = 2).

Описанным выше способом получаем границы учитываемых диапазонов и Dj значений пикселей для первого и второго столбца. Отбираются пары пикселей, значения которых принадлежат соответствующим диапазонам. Если число отобранных пар пикселей Nx меньше kt-h, где - установленный коэффициент (кА = 0.5), то учитываются все пары пикселей, и D.J, Dj в этом случае равны всему множеству значений пикселей. Медиана medx вычисляется по отобранным парам пикселей. Выражение для Sxy имеет вид:

1 -W*.y)~ 3(x-4,y1-tnedx)2

. Ъу = —е , 9(х,у) 6 Dl Л д(х + 1,у) е D| ^ (?)

k = 0, д(х,у) ÎD] V д(х + 1,у) g £)|

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

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

В работе даются временные и пространственные оценки предлагаемых алгоритмов. Производится экспериментальный анализ методов при работе на реальных и специально созданных изображениях с полосами. Для получения количественных оценок тестирование проводилось на данных с известными правильными результатами обработки J[xj>). В качестве одного из способов оценки, применялось вычисление СКО разности двух изображений. Для учета локальной гладкости выравнивания предложен специальный способ оценки. Производится вычисление разностей средних по столбцам fix,у) и jixy) и суммирование приращений получившихся разностей от одного столбца к другому. Таким образом, будет учитываться локальная гладкость выравнивания, а сильные перепады будут учтены только единожды. Формула имеет вид:

_ ly.tipg+p-pcni п,л _ si1,^о.о si1,trem m

5 =-—-,D О)--1-;---— (8)

где w, h - ширина и высота изображений fix, у) и f(x, у).

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

Выравнивание средних по столбцам Выравнивание средних полиномом Комбинир. метод Глобальный метод

СКО 12.24 12.45 3.93 9.26

S 0.48 0.43 0.11 0.19

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

Предлагается частотный метод устранения шумов в виде горизонтальных полосок на базе дискретного преобразования Фурье (ДПФ). Скорректированное изображение получается как /(х,у) = £(х,у)-м!(х,у)г1(х,у), где г\{х,у) - выделенное фильтром изображение шума, х,у) - весовые коэффициенты.

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

В данной задаче фильтрации частоты шумов достаточно размыты по спектру, поэтому подходящим способом является массовое подавление частот вдоль центральной вертикали путем умножения на маску. Были созданы маски фильтров с использованием зависимости из фильтра Баттерворта:

(9)

I О0 \

где О0 - заданная величина, а ¿>(м,у) обозначает расстояние от точки (и, у) до центра маски фильтра, п - параметр, регулирующий плавность кривой фильтра. Данный вид фильтра выбран ввиду его близости к фильтру Гаусса, но при этом он может иметь более крутые края.

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

Рис. 3. Изображение прямоугольного и треугольного фильтров.

Рис. 4. Шум, выделяемый прямоугольным и треугольным фильтрами.

Так как фильтр из-за своей ширины вызывает артефакты в определенных ситуациях, в дополнение ко всему предлагаются варианты взвешенного вычета шума, путем задания w(x,y). Первый предлагаемый способ вычисления - на основе минимизации локального СКО fix, у) при дополнительных ограничениях по значениям локальных СКО g(x,y) и ¡w(x,y)rj(x,y)\. Второй - w(xy) зависит только от локальных СКО g(x,y) и локальных средних ц(х,у). Различные пороговые значения для ограничений обосновываются. Меры направлены на вычет небольших значений шума на относительно гладких участках, что позволяет одновременно сохранить полезную часть информации и повысить визуальное качество изображения. На Рис. 5 представлен результат обработки фрагмента изображения.

Рис. 5. Фрагмент изображения и результат использования фильтра.

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

13

обрабатываемого блока. Благодаря этому сохраняется больше полезной информации. Картина шума ц{х,у) составляется из блоков, а вместо задания w(x,y) производится фильтрация г]{ху) по столбцам на основе среднего МО и среднего СКО для столбцов блоков в общем столбце путем ослабления и отсечения выбросов.

В работе даются временные и пространственные оценки сложности алгоритмов, демонстрируется работа предложенных и известных методов. Для тестирования используются реальные изображения с помехами для визуальной оценки, а также изображения с искусственно смоделированными шумами для числовой оценки. Для сравнения изображений используются оценки качества СКО и PSNR (peak signal-to-noise ratio). При сравнении с известным билатеральным фильтром и с аналогом клиновидной маски, используемой в зарубежной системе Erdas, получены следующие данные:

Клиновидная маска Билатеральный фильтр Треугольный фильтр с весовыми коэф. Блочный метод

СКО 14.39 4.74 2.41 2.25

PSNR 24.97 32.2875 41.2 41.07

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

Анализируются известные фотограмметрические и нефотограмметрические методы сшивки полос спутниковых изображений, в том числе используемые в известных системах ОгЙюЫогтЗсап, ЫогтЗсап, Фотоплан, выявляются их недостатки с учетом актуальных технических заданий на ОКР "Прибор-ОЭК-ИС35" и "НПО-АСР". Существующие методы сшивки с использованием фотограмметрии слишком медленны, сложны в реализации и требуют дополнительных данных. Необходимы универсальные, быстрые и точные методы, которые будут работать только с использованием самих изображений и съёмочной паспортной информации. Также требуется более совершенный фотограмметрический метод, который обеспечит в необходимых случаях точное следование модели съемки и учет всех имеющихся данных о КА.

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

14

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

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

ЯЛ) = ||В(5)||2 = ^^[Ши.иУ?,)-^)2 + + (Ю)

ЖЕ^ [{ГнЧ^'Укд-Г^коУкд)2 + (/Д*Й-УЙ) -

где В(5) - вектор-функция невязок точек, ы - вектор из всех параметров &>к, М- число изображений, /V* - число фиксирующих точек (ФТ) на изображении к, - число СТ между изображениями к и А+1, у^(), (х* у*() - координаты фиксирующей точки на изображении к и ее координаты на сшитом изображении, у™(), у™() -координаты сопоставленных СТ на изображениях к и А+1.

Были исследованы различные способы преобразований, в результате было выбрано проективное преобразование координат:

ГХ(Г .л _ аЦгСО+вЦР)х+ык(3)у ,у, . _ ы1с(6)+й)ц(7)д^+йц.(8)у ....

,к К ,У) 1+Шк(4>+Шк(5)у >1к*<Х'У)- 1+й,к(4)х+а1с(Б)у

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

Задача минимизации решается методом Гаусса-Ньютона с добавлением специальной весовой диагональной матрицы Р:

ш/+1 - Щ + (12)

где у - номер итерации, ЛДЙу) - матрица Якоби вектор-функции которая

вычисляется аналитически, /?у - параметр, регулирующий шаг. Благодаря диагональной матрице Р учитывается неточность данных о положении ФТ. Реальное положение ФТ может отличаться в пределах от нескольких пикселей до десятка. Веса балансируют записи в системе: ФТ присваиваются малые веса, а СТ - большие. В работе проводится обоснование выбора значений весов исходя из оценок точности совмещения по тестовым данным.

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

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

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

Со спутника поступают параметры полета и ориентации Я, Р, Ь, V, в, А (геоцентрический радиус, широта, долгота центра проектирования, углы тангажа, крена, рысканья). Уточняются положения матриц в фокальной плоскости х„ у, (/ - номер матрицы) и коэффициенты полиномов, аппроксимирующих Я, Г, Ь, V, в, А. Для уточнения коэффициентов используется информация об опорных и связующих точках. Опорные точки помогают уточнить ЭВО, а связующие точки — совместить изображения. Отдельно для каждой ПЗС матрицы уточняются положения матриц х„ у, и угол А, так как он находится в одной плоскости с углами положения матриц в фокальной плоскости, которые требуют уточнения. Остальные параметры являются общими для

всех ПЗС матриц. Все вместе они определяют новый вектор параметров си. Решается задача минимизации вектор-функции В (5) невязок точек в координатах на местности:

ЯоЗ) = ||В(3)||2 = ХК^ЕГЛ - *£!)2 + (С^.УГО - К")2] + (13)

Ж [(СМ.УЙ) - ^+1(хй,Яв())2 + (С*У(*Й-УЙ) - - гпт,

где из новых обозначений - число опорных точек на изображении к, (Х/$,Уы), Ук") ~~ координаты опорной точки на изображении к и се координаты на местности,

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

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

= + (А/Р А] + <Э)"М[РВу, (14)

где в данном случае Р — весовая диагональная матрица для опорных и связующих точек, 2 — весовая диагональная матрица для параметров ш. Таким образом, обеспечивается регуляризация задачи оптимизации. Веса вычисляются с использованием априорной информации о точностях определения точек и параметров К А. Матрица Якоби ЛДшу) вычисляется численным образом.

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

получены следующие оценки на тестовых изображениях при моделировании дефектных сдвигов и поворотов ПЗС матриц:

—___Метод Полином. Проективное Мультиквадр. Фотограм.

Оценка '— метод преобразование уравнений метод

По связующим 0.08 0.08 0.04 0.27

точкам (по х, У) 0.06 0.06 0.05 0.26

По ортофотоплану 0.86 0.66 0.92 0.54

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

Рассматривается экспериментальный моделирующий комплекс (ЭМК), разработанный в рамках ОКР "Прибор-ОЭК-ИС35", в котором используются предложенные алгоритмы. ЭМК производит моделирование процесса съемки и полный цикл обработки снимков до создания ортофотопланов. Более детально рассматривается подсистема радиометрической и геометрической коррекции (КП "Инфоплан-С ОЭК"), разработанная в ходе данной работы и включающая предложенные методы коррекции спутниковых изображений.

Приведены результаты экспериментального исследования параллельных методов обработки с применением библиотек Intel Cilk Plus, TPL .NET.

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

Основные результаты работы

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

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

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

4. На базе разработанных методов автором реализованы программные средства точной, качественной и быстрой обработки изображений с использованием 64-битной архитектуры и параллельных вычислений, которые используются в реальных комплексах обработки спутниковых видеоданных для ОАО «ЛОМО» и ОАО «Красногорский завод им. С.А. Зверева».

На защиту выносятся

1. Комбинированный локальный и глобальный методы удаления вертикальных полос на спутниковых изображениях.

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

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

4. Программные средства обработки спутниковых изображений и их внедрение в реальные комплексы в рамках ОКР "Прибор-ОЭК-ИС35" и "НПО-АСР".

Публикации в изданиях, входящих в перечень ВАК

1. Гусев В.Ю., Крапивенко A.B. Методика фильтрации периодических помех цифровых изображений. Электронный журнал Труды МАИ №50, 2012.

2. Гусев В.Ю. К вопросу фильтрации структурных искажений на космических оптико-электронных изображениях. Геодезия и картография, №12, 2012. -М.: ФГУП "Государственный картографический и геодезический центр". С. 111-114.

3. Барабин Г.В., Гусев В.Ю., Зайцев В.Е., Юров H.H. Комбинированные методы радиометрической коррекции характерных артефактов спутниковых изображений. Вестник МАИ, №3, 2013. -М.: МАИ. С. 160-165.

4. Барабин Г.В., Гусев В.Ю. Фотограмметрический метод построения единого изображения при спутниковой съёмке секционированным датчиком изображений. Электронный журнал Труды МАИ №71, 2013.

5. Барабин Г.В., Гусев В.Ю., Зайцев В.Е., Юров H.H. Методы построения единого изображения при синхронной спутниковой съёмке секционированным датчиком изображений. Вестник компьютерных и информационных технологий, №4, 2014. С. 15-21. -М.: ООО "Издательский дом "Спектр".

6. Гусев В.Ю. Глобальный метод радиометрической коррекции характерных артефактов спутниковых изображений. Нейрокомпьютеры: разработка, применение. -М.: Издательтво "Радиотехника", (в печати)

Публикации по теме диссертации в других изданиях

7. Гусев В.Ю. Методы сшивки перекрывающихся полос спутниковых изображений на начальной стадии обработки. 11-ая Международная конференция "Авиация и космонавтика- 2012". Тезисы докладов. -СПб.: Мастерская печати, 2012. С. 253-254.

8. Гусев В.Ю. Геометрическая коррекция спутниковых изображений. В сб.: Ломоносов - 2013: XX Международная научная конференция ... секция «Вычислительная математика и кибернетика»: -М.: Издательский отдел факультета ВМиК МГУ; МАКС Пресс, 2013. С. 102-104.

9. Гусев В.Ю., Хлебников Е.Е. Параллельная обработка спутниковых изображений. Сборник трудов участников III Всероссийской научно-практической конференции «Современное непрерывное образование и инновационное развитие». - Серпухов: МОУ «ИИФ», 2013. С. 116-119.

10. Гусев В.Ю. Комбинированный метод коррекции последствий неравномерности чувствительности датчиков спутниковой съемочной системы сканового типа. Труды конференции ГрафиКон'2014: 24-я Международная конференция по компьютерной графике и зрению. Ростов-на-Дону, Академия архитектуры и искусств ЮФУ. С. 69-71.

20

Подписано в печать: 11.10.2014 Тираж: 100 экз. Заказ № 1250 Отпечатано в типографии «Реглет» г. Москва, Ленинградский проспект д.74 (495)790-47-77 www.reglet.ru