CC BY
45
КЛАСТЕРНАЯ ТЕХНОЛОГИЯ ФОРМИРОВАНИЯ И ПАРАЛЛЛЕЛЬНОЙ ФИЛЬТРАЦИИ БОЛЬШИХ ИЗОБРАЖЕНИЙ
С.Б. Попов, В.А. Сойфер, А.А. Тараканов, В.А. Фурсов Институт систем обработки изображений РАН Самарский государственный аэрокосмический университет
Формулировка проблем
В настоящее время существует большой интерес к изображениям, полученным при дистанционном зондировании Земли, с разрешением около 1м и широкой полосой захвата. На начальных этапах формирования таких изображений по данным аэрокосмического мониторинга обычно приходится выполнять фильтрацию исходных изображений с целью улучшения их качества и решать задачу совмещения больших фрагментов изображений при составлении мозаики из нескольких снимков на протяженные участки и территории (фотокарты местности). Размеры регистрируемых изображений обычно таковы, что их обработка возможна лишь с применением многопроцессорных, в т.ч. кластерных, систем [1]. Эффективность параллельных программ напрямую зависит от соответствия алгоритма архитектуре вычислительной системы. Поиск и реализация алгоритмов, структура которых допускает эффективное отображение решаемой задачи на многопроцессорную вычислительную систему с укрупненным распараллеливанием (MPP-систему) является, безусловно, важной проблемой.
В настоящей работе рассматриваются две задачи, решаемые при обработке изображений больших размеров. Это задача оценки и коррекции искажений изображений, обусловленных влиянием оптической системы, ПЗС-приемников и динамики процесса видеоизмерений. А также задача автоматического совмещения фрагментов изображений путем совместной обработки перекрывающихся зон (краев) соседних исходных изображений.
В настоящей работе обсуждаются возможные технологии формирования и обработки изображений и описана оказавшаяся наиболее эффективной, с точки зрения ускорения, схема обработки изображений больших размеров, реализованная на кластере Самарского научного центра РАН.
Построение алгоритма фильтрации
В работе исследовалась возможность использования для обработки больших изображений линейного фильтра с бесконечной импульсной характеристикой (БИХ-фильтра):
8(пъ П2 ) = - Е Е ат1,т2 8(п1 - П2 - т2 )
(т1,т2
+ Е Е Ьт1,шг /(п1 - ^ п2 - т2 ) , (1)
(т1,т2 )&2/
где / (п1, п2), 8 (п\, п2) - отсчеты входного и выходного изображений соответственно, а а,т, Ь -коэффициенты фильтра. Предполагается, что пара-
метры фильтра определяются по малым фрагмента изображений так, как описано в работе [2]. На этом этапе решается также задача обеспечения устойчиво -сти фильтра (1).
В случае использования БИХ-фильтра при фиксированной расфокусировке можно существенно уменьшить размер опорной области по сравнению со случаем использования фильтра с конечной импульсной характеристикой (КИХ-фильтра). Это позволяет существенно снизить вычислительную сложность алгоритма на одной итерации.
Однако возникают проблемы реализуемости БИХ-фильтра, связанные с тем, что для вычисления выходных отсчетов должны использоваться отсчеты выходного изображения, которые для наиболее часто используемых форм опорной области, как правило, еще не определены. Для преодоления этого затруднения применяют итерационную схему вычисления отсчетов выходного изображения [3]
у(п\, П2) = а(п^ П2)* *х(п, П2) + е(пг, * у(п, П2), (2)
где звездочки означают двумерную свертку, а
c(n\, п2 )= 1 - Ь(n\, п2 ) .
В частотной области итерационная схема по аналогии описывается соотношением
У (в*1, е*2) = А(е*1, е*2) X (е*, е*2) +
+ С(е**, е*2)У (е*, е*2), (3)
где С (в*1, е*2) = 1 - В(е*, е*2) .
Для обеспечения устойчивости итерационной схемы (2), (3) вводится параметр Я, удовлетворяющий неравенствам
0 <Л< 2/тах В(е*1,е*2) (4)
( *1, *2 )
так, что частотный отклик фильтра принимает вид
У (е*1, е*2) = ЯА(е*, е*2 )Х(е*1, е*2) +
+ С(е*1, е*2)У-_1(е*1, е*), (5)
где С (е*1, е*2) = 1 -ЯВ(е*1, е*2).
Имеет место связь точности восстановления с величиной параметра Я и числом итераций. Ниже приводится пример выбора этих параметров, обеспечивающих компромисс вычислительной эффективности и качества восстановления.
Совмещение фрагментов изображений
На этапе формирования выходного изображения производится совместный анализ и обработка полос. Сначала, путем попарного сравнения перекрывающихся зон (краев) соседних полос, осуществ-
ляется субпиксельная оценка их геометрического рассогласования. По результатам этой оценки производится уточняющий пересчет параметров геометрического положения строк входных изображений. После этого формируется выходное изображение. Его отсчеты получаются путем интерполяции значений функции яркости по ближайшим отсчетам входных изображений, с заданными шагами равномерной пространственной дискретизации в двумерной системе координат фокальной плоскости.
Алгоритм оценки геометрического рассогласования полос предназначен для определения относительного локального сдвига строк соседних полос изображений. Для этого используются одноименные фрагменты изображений, имеющиеся на перекрывающихся участках (в зонах перекрытия) полос. Основные этапы работы алгоритма следующие:
- выделение и определение координат информативных фрагментов на перекрывающемся участке одного из изображений;
- определение координат соответственных фрагментов на перекрывающемся участке соседнего изображения;
- определение параметров геометрической трансформации положения строк по найденным координатам соответственных точек в зонах перекрытия соседних полос.
Условием применимости алгоритма является наличие на соседних полосах изображения перекрывающихся участков, на которых зарегистрирован один и тот же участок подстилающей поверхности, и имеются так называемые «информативные фрагменты», функции яркости которых позволяют однозначно определять соответственные координаты в любом из заданных направлений при их совмещении по критерию максимальной близости.
Метод определения координат информативных точек заключается в следующем [3]. Для каждого положения прямоугольного фрагмента формируется матрица X, для которой затем вычисляется мера информативности Ф(Ц). В качестве мер информативности используется так называемый показатель диагонального преобладания соответствующей матрицы сопряженности [4]:
Ф(К) = N 2/(Ж + X ) I, у = 1, N, I Ф у.
(6)
Элементы гуу ,1, у = 1, N,у ф у матрицы сопряженности
1
'21
Г
г
12
г
1т
Чт
' т1 Гт2 ••• 1
определяются по соотношениям
(7)
гу ='
( Xу XI ) ( Xу )
(8)
Матрица X и соответствующая ей матрица Я составляются отдельно для каждого фрагмента, «отвечающего» за свое направление. В частности, из соображений обеспечения высокой надежности алгоритма, предусматривается вычисление показателей информативности Ф(Я) для двух направлений: вдоль и поперек полосы, а также под углом 45° к направлению маршрута. Число строк Ь} и число столбцов Ь2 матрицы X для каждого прямоугольного фрагмента определяется размерами Ь]хЬ2 соответственно и может быть различным для разных направлений ориентации прямоугольников.
Отбор информативных фрагментов осуществляется по критерию максимума минимального значения Ф(Ц) на множестве показателей, подсчитанных для каждого заданного положения фрагмента, при условии, что для всех 5=1,4 величина
Н
'/М =
(хГ х )
I/ М, I = 1, м
(9)
где х1, х у - строки с номерами I и у ^х^-матрицы X.
не менее некоторого заданного порогового значения Н. Если
Н<Ндоп , (10)
то соответствующий фрагмент пропускается. Введение этого дополнительного требования снижает опасность ложного отнесения фрагментов к числу информативных при низком уровне яркости, т. к. высокое значение показателя диагонального преобладания может быть связано с наличием шумов, а не с характером полезного сигнала на фрагменте.
Достоинством предложенного метода определения информативности фрагментов является его согласованность с корреляционным методом определения соответственных точек. Относительно большая величина показателя Ф(Я) свидетельствует о слабой корреляции строк матрицы X с различающимися индексами. Следовательно, функция корреляции отобранных по этому показателю фрагментов будет иметь выраженный максимум именно в точке совпадения одноименных фрагментов.
Фактические выходные координаты одноименных фрагментов изображений в зонах перекрытия определяются на основе величины, измеряющей корреляцию между функциями яркости фрагментов изображений заданных размеров в зоне поиска.
Вычислительная схема распараллеливания по данным
Информационная технология обработки данных аэрокосмического мониторинга допускает различные схемы распараллеливания. В настоящей работе исследуется схема декомпозиции по данным, учитывающая специфику и формат обрабатываемой информации. На этапе фильтрации параллельно обрабатываются неперекрывающиеся полосы изображения. При наличии достаточного количества процессоров в параллельной вычислительной системе изображение распределяется на такое количество полос, которое обеспечивает равномерную загрузку всех процессоров, участвующих в вычислениях. При этом размер фрагментов (полос) будет существенно меньше всего изо-
х
1
бражения, что позволяет целиком размещать данные в оперативной памяти, не затрачивая дополнительное время на буферизацию изображения в памяти постоянного хранения. При такой организации процесса параллельных вычислений необходимо на всех шагах итерационного процесса обеспечить обмен крайними строками между процессами, обрабатывающими соседние фрагменты изображения.
Такая схема слабо чувствительна к скорости коммуникационной среды, однако очень чувствительна к правильной реализации обмена данными. Заметим, что это преимущество связано с использованием БИХ-фильтров. Применение КИХ-фильтров приводит к существенному увеличению объемов данных, которыми должны обмениваться процессоры.
Аналогичная проблема возникает в задаче автоматического совмещения фрагментов изображений при совместной обработке перекрывающихся зон (краев) соседних исходных изображений.
Пример реализации БИХ-фильтра на кластере Был реализован БИХ-фильтр, описываемый передаточной функцией
° ('22 )= В-)
В (21 / 22 )
(11)
где
А (21 / 22 ) = + а1 (21 + 21-1 + 22 + 221) + а2 (2122 + 21-12-1 + 2122 1 + 2-22 ) В ( 21, 22 ) = 1 + Ъ1 ( 21 + 2—1 + 22 + 2-1) +Ъ2 (2122 + 2—12-1 + 2122Х + 21-122 ) .
Параметры фильтра - а0 = 35,2857, а1 = -1,2,
а2 = 0,27 , I) 1 = 1,2, I)2 = 0,27 задавались, исходя из требования обеспечения качества восстановления, полученного моделированием исходного изображения при условии устойчивости БИХ-фильтра.
Передаточной функции (11) соответствует опорная область 3x3. Следовательно, для обеспечения работы описанного выше параллельного алгоритма после каждой итерации процессы, соответствующие соседним полосам изображения, должны обмениваться одной строкой. Параметр Я итерационной схемы (5) принимался равным 0,01. При таком значении качество восстановления изображения получается вполне приемлемым при количестве итераций 5-10.
200 180 160 140 120 100 80
/ ч\ '••л\ \ч\ —
\\
— — —
300
260
220
180
140
100
ч ч \ ч. ч •• \
0 2 4 6 8 10 12 14 16 0 2 4 6 8 10 12 14 16
Рис. 1. Зависимость времени выполнения программы от количества процессоров для 1, 5, 10 итераций
На рис. 1 (а) приведены графики изменения времени обработки изображения размером 4096x4096 для 1-й, 5-й, 10-й итераций при реализации алгоритма фильтрации на разном количестве процессоров (нижняя кривая соответствует 1 итерации, верхняя - 10). Одному процессору соответствует обычный, то есть последовательный алгоритм.
На рис. 1 (б) приведены графики, полученные при тех же условиях для 10-й и 20-й итераций. Из графика видно, что при двадцати итерациях, обработка ускоряется в два с половиной раза. Интересно, что при увеличении числа процессов более восьми время обработки практически не изменяется. Это означает, что увеличение времени, затрачиваемого на передачу, происходит также быстро, как и сокращается время вычислений за счет увеличения числа процессоров.
На рис. 2. приведены исходное и обработанное изображения, полученные при значении параметра Я =0,01 и числе итераций 10.
Заключение Применение БИХ-фильтров с малыми размерами опорной области позволяет строить эффективные алгоритмы параллельной обработки больших изображений с небольшими накладными расходами. Эффективность алгоритмов зависит от размеров обрабатываемых изображений. Алгоритм совмещения изображений по информативным фрагментам также допускает эффективное распараллеливание по данным, что открывает возможность строить эффективные кластерные технологии обработки данных аэрокосмического мониторинга.
Рис. 2. Пример обработки участка изображения (использован демонстрационный снимок с www.sovinformsputnik.com)
Работа выполнена при частичной поддержке
РФФИ (гранты №01-01-00097, 00-01-05001).
Литература
1. Попов С.Б., Фурсов В.А. Кластерная технология конвейерной фильтрации // Труды международной конференции. Математическое моделирование - 2001. Самара, 13-16 июня 2001. С. 137-140.
2. Сойфер В.А., Котляр В.В., Фурсов В.А. Построение алгоритмов оперативной коррекции искажений на изображениях в оптико-электронных
системах наведения и целеуказания // В сб. Современные методы проектирования и отработки ракетно-артиллерийского вооружения. Саров. ВНИИЭФ, 2000. С. 340-345
3. Методы компьютерной обработки изображений / Под ред. Сойфера В.А., М., Физматлит, 2001.
4. Фурсов В.А. Идентификация моделей систем формирования изображений по малому числу измерений // Самара, ИПО СГАУ, 1998.
