CC BY
24
ОБРАБОТКА ИЗОБРАЖЕНИИ
оценивание пространственно-зависимых искажений и корректирующих фильтров кусочно-постоянными линейными моделями
В.В. Сергеев, В.А. Фурсов, М.В. Максимов* Институт систем обработки изображений РАН, *) Самарский государственный аэрокосмический университет
1. Введение.
Для решения задач коррекции (реставрации) изображений необходимо знание модели искажений. Один из широко распространенных видов искажений - дефокусировка. В оптических системах формирования изображений искажения этого типа часто оказывается неизопланатичными. Известно, что задача компенсации указанных искажений может быть решена применением линейных восстанавливающих фильтров с перестраиваемыми по пространственным координатам параметрами. В некоторых случаях модель неизопланатичных искажений, вносимых оптическими системами, может быть представлена в виде параметрического семейства двумерных функций. При этом на этапе восстановления параметры восстанавливающего фильтра могут быть вычислены в каждой точке обрабатываемого изображения. Параметры этой функции могут быть определены в результате решения задачи идентификации по тестовым изображениям. Из соображений упрощения целесообразно задачу идентификации модели пространственно-зависимых искажений свести к совокупности задач оценки пространственно-инвариантных моделей на малых, насколько это возможно, фрагментах изображения. При этом возникает проблема отбора малых фрагментов изображений, на которых данные наблюдений являются информативными, т. е. таких на которых задача идентификации может быть решена с достаточной точностью. Связанному с этой проблемой кругу вопросов и посвящена настоящая работа.
2. Постановка задачи.
Характерным для искажений, вносимых оптическими системами, является неравномерность расфокусировки связанная с аберрациями различных порядков. Предположим, что эти искажения могут быть описаны уравнением фильтра с конечной импульсной характеристикой (КИХ- фильтра) [1]. Тогда для опорной области в виде квадрата со стороной N (К нечетно) связь отсчетов - х(п1,п2) входного (неискаженного) и отсчетов у(п1, п2) выходного (подвергшегося искажениям) изображений описывается следующим образом
У (Пу, щ) = X к( к1,к 2)
X
к1,к2
X (п1 — к1, п2 — к2) + Х(п1, п2),
= ЖЖ^ = NN
к л --, ,к п. --,
1 2 2 2 2 2
(1)
где к(к1,к2) - значения импульсной характеристики, а Х(п1, п2) дискретная последовательность, включающая погрешности измерений и ошибки связанные с ограничением порядка модели.
Предположим, что значения импульсной характеристики задаются параметрическим семейством функций вида [2]
к(к1, к 2) =
1
(
л/2рс
-ехр
г{к 1, к2 ) 2а2
2
(2)
с параметром а зависящим от пространственных координат: а = а(п1,п2). Здесь к(к1,к2), - значение импульсной характеристики в точке с координатами
(к1,к2), а г(к1,к2 ) = ^к12 + к2 - расстояние от
центра опорной области до этой точки.
Фрагмент на изображении, для решения задачи идентификации, выберем малым настолько, чтобы в разностном уравнении (1) параметры к(к1,к2), в
силу (2) зависящие от параметра а, можно было считать постоянными. Тогда, если этот фрагмент содержит М различных отсчетов у( п1, п2) выходного изображения, в соответствии с (1) можно записать матричное соотношение
у = ХЬ+х,
(3)
где у -Ых1 -вектор, составленный из отсчетов искаженного изображения; X - ЫхМ -матрица, каждая строка которой соответствует одной опорной области и составлена из М отсчетов х (п1 — к1, п2 — к2)
неискаженного изображения; Ь - вектор составленный из М соответствующих значений импульсной характеристики - к(к1,к2); а - £ - лх1-вектор, составленный из ошибок Х( п1, п2).
Для модели (3) МНК- оценка к вектора параметров И имеет вид [3]
h = [xTx]"xTy.
(4)
Известно, что эта оценка весьма чувствительна к (особенно выделяющимся) ошибкам в исходных данных. В работах, [4,5] предложено эту оценку использовать лишь как начальное приближение для последовательности МНК-оценок вида (4), которые
строятся на наборах данных Xполученных из
исходных путем преобразований вида
X = СХ,у = Су, где О - диагональная ЫхЫ-
матрица вещественных чисел. Матрица С задается после каждого шага вычисления оценок по соотношению (4) на основе анализа вектора невязок X = у - ХИ , так чтобы придавался меньший вес наблюдениям, для которых соответствующая компонента вектора невязок оказалась большой.
При реализации описанной схемы идентификации основная проблема заключается в том, что ошибка идентификации
Ah = h - h = [xt х]-1 ХТ X
(5)
может быть чрезвычайно большой, даже если норма вектора ошибок в исходных данных - X невелика. Это может происходить как вследствие того что на малом фрагменте изображения отсутствуют какие-либо информативные элементы (фрагмент попал на участок фона), так и из-за задания небольшого числа выделяющихся весовых коэффициентов в матрице С. В любом из указанных случаев задача идентификации может оказаться плохо обусловленной, или даже вырожденной. Следовательно, в данной схеме идентификации актуальной является задача построения процедур автоматического контроля информативности видеоданных.
Для обнаружения факта вырожденности и количественной оценки обусловленности в работах [4,5] предложены решающие правила, основанные на оценке величины диагонального преобладания матрицы А = ХТХ . В частности, показано, что задача хорошо обусловлена и гарантируется решение задачи определения параметров с высокой точностью, если имеет место неравенство
M-1 + £<ф<M ,
(6)
где ф = №2А/¡ГА ,1 <ф<М , 1гА - след матрицы А, а е - заданное малое (не более 0.1) положительное число. Важное достоинство достаточных условий (6) - их вычислительная простота. Однако, при уменьшении размеров фрагментов ситуация часто оказывается такой, что неравенство (6) не выполняется, хотя задача, в принципе, еще может быть решена с удовлетворительной точностью.
В настоящей работе исследуется возможность применения для отбора малых фрагментов изображений более эффективной количественной оценки обусловленности, основанной на непосредственном анализе степени сопряженности векторов-столбцов матрицы Х c нуль-пространством этой матрицы.
3. Меры информативности видеоданных на фрагментах
Причиной плохой обусловленности матрицы A = ХТ Х является "почти" линейная зависимость (мультиколлинеарность) [3] векторов-столбцов матрицы Х. В качестве мер мультиколлинеарности предлагается использовать величины
S . = min|S.|, (7)
min i ' v 7
S i = I St
(8)
f
где St =
I b2
Л/2
V j=
а РJ = (х ■1,1 j) - проекция вектора х. на ^й вектор 1 базиса, образованного собственными векторами,
соответствующими нулевым собственным значениям матрицы ХМ-1ХМ-1. Здесь ХМ_1 - Ых(М-1) матрица составленная из М-1 нормированных (||х.. || = 1,. = 1, п ) векторов-столбцов матрицы Х путем исключения вектора х.. Геометрически Si - "косинус" угла между вектором х.. и нуль- пространством подпространства, натянутого на остальные М-1 векторов. Покажем это.
В соответствии с указанной выше геометрической трактовкой определим характеристику Si в виде
(х.,(х.
S.. =i
x, X x,
(9)
Здесь х. - один из векторов размерности Ых1, из которых составлена матрица Х,
х. = Хм. (10)
- проекция этого вектора на подпространство, натянутое на остальные М-1 векторов, а
z t =[хМ-1Х M-1 Х -1х
(11)
разложение вектора х. по базису, образованному пространством матрицы ХМ -1 .
Подставляя (11) в (10), а затем в (9) после соответствующих преобразований [6] получаем
S. =(х Т T TT х i )/(х Т х i )>
(12)
где Т0 - матрица, размерности Ых-М+1 составленная из Ы-М собственных векторов, соответствующих ну-
левым собственным значениям матрицы ХмЧХМч . Числовые значения величины вычисленные по соотношению (11), совпадают с вычисленными по исходному соотношению (7) при ||х = 1.
4. Сравнительный анализ мер
Можно показать тесную связь меры (7) с максимальной сопряженностью [3], определяемой как:
R = max Rt R = <XM_1z, х,)Д = 1,M
(13)
где z = [XM-1XM-1 ] 1XM-Iхi ,
Геометрически R. - "косинус" угла между вектором х i и подпространством, натянутым на остальные M-1 векторов. Можно показать, что
M
Ri = ^ а2 , где aj = (х ., f) - проекция вектора
j=1
х i на j-й вектор f j базиса, образованного собственными векторами, соответствующими ненулевым собственным значениям. Ясно, что с ростом величины R величина Smin стремится к нулю. Если задача вырождена, R=1, при этом Smin =0. Поэтому меру (7) можно назвать показателем минимальной сопряженности с нуль- пространством. Меру (8) по аналогии назовем показателем суммарной сопряженности с нуль- пространством.
Мера (7) по сравнению с (13) обладает важными преимуществами. Во-первых, для определения матрицы T0 - не требуется решение полной проблемы собственных значений, что при сильной мульти-коллинеарности представляет значительные трудности. В данном случае достаточно определить лишь какие строки матрицы являются линейно-зависимыми. Исследования показали также, что процедура оценивания показателей Smin и s^-обладает более высокой вычислительной устойчивостью именно в ситуациях, когда обусловленность задачи вызывает сомнения.
Объем вычислительной работы для вычисления мер (7), (8) по сравнению с принятием решения на основе неравенства (6), конечно, намного больше. Поэтому использование этих мер может быть оправдано лишь высокой стоимостью риска, связанного с возможными ошибками идентификации. Необходимо иметь ввиду, что объем вычислительной работы при использовании меры (7) при фиксированном N уменьшается с ростом m, а вычислительные преимущества по сравнению с мерой (13) проявляются в полной мере при достаточно малых значениях N/M. Относительная вычислительная эффективность метода возрастает с уменьшением этого отношения (на малых фрагментах изображения). В частности, если объем вычислительных операций в расчете на каждый из векторов t j и f j примерно одинаков, то
метод эффективнее, если размерность нуль- пространства, по крайней мере, не превышает размерности пространства параметров.
5. Связь погрешности идентификации с мерами информативности
Проводились вычислительные эксперименты с целью оценки эффективности процедуры контроля информативности данных на малых фрагментах изображений с использованием показателя минимальной сопряженности с нуль пространством - 8тйп. В качестве тестового для этой цели использовался участок типового изображения типа "текст", показанный на рис. 1.
Для наглядной визуальной оценки связи показателя 8т„ с информативностью исходных данных он был рассчитан по соотношению (11) в каждой точке исходного изображения "текст". Результаты расчетов представлены в виде изображения на рис. 2. Обращает на себя внимание факт, что показатель "безошибочно замечает" строки и слова, т.е. наиболее информативные, собственно представляющие интерес для идентификации, участки изображения.
06 работка изображ связана с решением кик задач, в катары: входные, и выхадны ные являются изобр ниями. Одним из пр рав служат системы редачи изображена оазцаботчики стал кI
Рис. 1. Исходное изображение.
Рис. 2. Поле значений показателя S.
Рис. 3. Расфокусированное изображение.
Наряду с указанной оценкой качественного характера, исследовались также количественные оценки связи погрешности идентификации с величиной показателя минимальной сопряженности с нуль-пространством и с величиной показателя суммарной сопряженности. Для проведения этих исследований наряду с тестовым изображением (рис. 1) использо-
валось изображение, полученное из него "расфокусировкой", показанное на рис. 3. При моделировании этого изображения параметр расфокусировки принимался пока одинаковым по всему изображению и равным <г=1.5. Соответствующий пространственно-инвариантный КИХ-фильтр имел опорную область 9х9.
Далее считая изображение на рис.1 входным, а изображение на рис.3 выходным, проводилась идентификация на малых фрагментах изображений с оценкой их информативности по показателям Smin и SS. Для получения количественных оценок связи погрешности идентификации со значениями показателей Smin и SJJ на искаженное тестовое изображение (рис.3) дополнительно накладывался аддитивный белый гауссовский шум, при отношениях шум-
сигнал 0.005, 0.01, 0.03 и 0.05. Изображение разбивалось на фрагменты 6x6, на каждом фрагменте определялся вектор оценок параметров импульсной характеристики и вычислялись относительные погрешности идентификации.
Полученные на всех фрагментах значения показателей были выстроены в вариационный ряд, и разбиты на группы по 20 значений в каждой группе. Затем показатели Smin, SS и соответствующие им относительные погрешности идентификации усреднялись для каждой группы отдельно. Полученные описанным способом результаты представлены в виде графиков на рис. 4., 5. На графиках видна устойчивая тенденция к снижению погрешности идентификации по мере увеличения показателей Smin и ss.
7.
К,
U \ / \
'Ч " \ 1 V
I ч ■
v
----0 6i 5%
......0 6i 3%
— 0 6i 1%
-0 6i 0.5%
0,8 0,7 0,6 0,5 0,4 0,3 0,2 0,1 0,0
0,0005 0,0295 0,0422 0,0553 0,0670 0,0788 0,0908 0,1166 Рис. 4. Связь погрешности идентификации с показателем Smin.
0,8 0,7 0,6 0,5 0,4 0,3 0,2 0,1 0,0
0,6829 1,0125 1,2447 1,4510 1,7550 2,0711 2,3606 3,0499
Рис.5. Связь погрешности идентификации с показателем SS.
Приведенное семейство зависимостей показывает возможность обоснованно задавать пороговое значение показателя Smin или ^^ при котором можно ожидать заданную точность оценивания модели. Для повышения надежности такого прогноза, указанные
семейства зависимостей следует строить заново для каждого нового класса изображений и характера помех в исходных данных. Заметим, что большие значения показателя SJJ не могут гарантировать невырожденность задачи, поэтому далее для отбора ин-
формативных фрагментов используется показатель с
^min ■
6. Идентификация корректирующих фильтров для восстановления пространственно-зависимых искажений
Решение задачи идентификации параметров модели пространственно-зависимых искажений и построение корректирующих фильтров предлагается осуществлять в виде следующей последовательности этапов.
1. Оценка параметров импульсной характеристики на малых фрагментах изображения.
2. Оценка на каждом малом фрагменте параметра расфокусировки о .
3. Оценка параметров функции, описывающей изме-
нения параметра расфокусировки по полю изображения (модели неизопланатичности).
4. Оценка параметров восстанавливающих фильтров, соответствующих полученным значениям параметра расфокусировки о.
Фрагменты для идентификации отбираются по заданному пороговому значению показателя Smin. Для каждого выбранного фрагмента задача идентификации параметров импульсной характеристики (этап 1) решается, как описано в разделе 2 (более подробное описание алгоритма можно найти в работе [5]).
По полученным на каждом фрагменте оценкам значений импульсной характеристики -
h (k1, к2), k1, к2 = 1,М далее строится оценка о
параметра о гауссиана (2). Для этого решается задача определения
0* -Q(о* ) = mi? Q(о),
оеп
где s - номер фрагмента, а х. - область допустимых значений параметра расфокусировки. В частности, для получения МНК - оценки О s параметра О s,
N 2
N 2
(
Q(оs)= Z Z
1
к =--к? =--V
1 2 2 2
V2Pi
-X
яо„
x exp
ÍK к 2)
2
2 о
"h (k1, к 2 )
У
У
Далее (этап 3) по совокупности а строятся оценки параметров модели неизопланатичности а=а(п1,п2). Для этого необходимо задать параметрическое семейство функций "хорошо" описывающих характер неизопланатичности. Для выбора этого класса функций можно осуществлять предварительный визуальный просмотр полученного множества оценок а 5. Например, если расфокусировка "плавно" увеличивается по мере удаления от центра изо-
бражения, удобно изменение параметра а описать параметрическим семейством функций вида
а(п1, п2 ) = 1 П + 1 уп22 +
+ 1 ху П1П2 +Рхп1 +Р уп2 +ао где п1,п2 - текущие координаты точки, а 1 х ,1 у, 1 ху ,Ь х Ь у ,а0 - константы, зависящие от степени не-
изопланатичности.
В настоящей работе мы ограничились рассмотрением частного случая указанной зависимости, когда неравномерность расфокусировки обладает радиальной симметрией. Параметрическое семейство функций описывающих такую модель неизоплана-тичности имеет вид
(14)
о( П1, n2) = 1 xn\ + 1 n +о 0
где Ях =Яу.
Четвертый заключительный этап идентификации - определение параметров восстанавливающих фильтров. С использованием оцененной на третьем этапе модели неизопланатичности, в принципе, в любой точке изображения можно строить соответствующий восстанавливающий фильтр. Однако при этом объем вычислительной работы оказывается весьма значительным. Для сокращения времени обработки изображений целесообразно заранее построить множество фильтров, соответствующих различным значениям параметра а, в виде некоторой таблицы соответствия, а затем в ходе обработки выбирать "ближайший".
Для составления такой таблицы можно воспользоваться методом непосредственной идентификации характеристик инверсного тракта, описанным в работе [7]. Для этого необходимо осуществлять идентификацию параметров восстанавливающего фильтра при различных искажениях, полагая входным - искаженное, а выходным - неискаженное изображение. Достоинство такого подхода заключается в том, что сразу на этапе оценивания параметров восстанавливающих фильтров можно учесть влияние помех, подвергая искаженное изображение за-шумлению с заданной интенсивностью. Получаемые таким образом фильтры оказываются регуляризо-ванными в смысле критерия, заданного при идентификации инверсного тракта.
7. Результаты экспериментальных исследований.
Для экспериментального исследования эффективности предложенной схемы построения корректирующих фильтров использовалось изображение, полученное из исходного (рис.1) моделированием пространственно- зависимых искажений. Расфокусировка осуществлялась КИХ-фильтром (1), значения импульсной характеристики которого определялись по формуле (2), а параметр а в каждой точке изображения рассчитывался по соотношению (14) для
1 х =1 у =0.000227534;а0 =0.931981.
2
При указанных параметрах модели неизопла-натичности размер опорной области увеличивался от центра изображения к периферии, так что минимальный размер маски был равен 5, а максимальный - 13. Полученное при этих параметрах модели искаженное изображение показано на рис 5.
Идентификация осуществлялась в виде последовательности описанных в предыдущем разделе четырех этапов. В результате реализации первых двух этапов на фрагментах изображения, отобранных по показателю Smin, было получено множество оценок параметра расфокусировки. Далее с использованием полученного множества оценок параметра расфокусировки о реализовывался третий этап - оценка параметров функции, описывающей изменения параметра о по полю изображения (модели неизопланатично-сти). Результаты оценивания приведены в таблице 1.
Таблица 1
Оценки параметров Абсолютная погрешность Относительная погрешность
1x_ly 0.00023059 0.00000306 0.01345293
Оо 0.927598 0.004383 0.00470288
Рис. 6. Искаженное изображение.
Для получения таблицы соответствия параметра о и параметров инверсных фильтров (4-й этап) тестовое изображение "текст" (рис. 1) подвергалось искажениям в соответствии с моделью (1),(2) при различных значениях параметра о. Последовательность параметров о задавалась в диапазоне 0.5-3.5. с шагом 0.1. Полученное таким образом для каждого фиксированного значения о искаженное изображение использовалось в качестве входного для идентификации параметров восстанавливающего фильтра.
Для иллюстрации эффективности описанной схемы идентификации с использованием построенного таким образом набора фильтров осуществлялась обработка изображения другого текста, (рис.7 а) подвергшегося искажениям, модель которых была такой же как и на тестовом (рис 7 б).
Восстановленное изображение приведено на рис. 7 в). На рисунке 7 г) для сравнительной визуальной оценки приведено изображение этого же текста, полученное обработкой "средним" фильтром, т.е. инверсным фильтром построенным по той же методике, но без учета пространственной неинвариантности искажений. Как и следовало ожидать, восстановление произошло лишь на отдельных участках изображения.
8. Заключение
Эффективность предложенной схемы идентификации параметров модели и корректирующих фильтров пространственно-зависимых искажений проиллюстрирована на частном примере, модель неизопланатичности которых обладает радиальной симметрией. Однако эта схема может использоваться для более широкого класса искажений. Успешность решения задачи в каждом конкретном случае будет зависеть от того, насколько удачно подобрано параметрическое семейство двумерных функций, описывающих изменение параметра расфокусировки по пространственным координатам.
Б задачах обработки из сообщениями яелягофся слу-нйпые параметры и оп редэпепие которых и печпой целыо иптерпре жения. Это могут быть,| ма, раз мер ы, о риента ц ге мое расположение дета ний, пара метры .опредег ,
а)
1$CiSfidHfcvtHta шЪ
ft
SfiJfeMPIWP» ОбЦиМ^И-*! ЧЙ1
с-йаДи^талм-и Pt«HWWt# -^«n-mit-i
rupa».tet pun I* ■■ гираылп' ра.
Tin f«Д№1 ОПИО которые Ф= -.¡О рвддоо*чи® к0*Ор»< *•
и целью жтврр-рв i-npMf43rt дел ыф иттврг-fi®
Это могут выТЦ' .<ip JrfO могут быт-р »
f-, ы*р ц о i я
hici* рвлпожлдаи» ■!•:'*■ с^яс-псло*«^«» н^-т
м-*, рял мор ы,о ри«чт#^м
в)
«•».м^вииетфичещ«** ^ ■■ - ■ .......
Рис. 7. Результаты восстановления изображения "текст": а) - неискаженное; б) - искаженное; в) - восстановленное; г) - восстановленное "средним " фильтром.
Литература
1. Даджион, Мерсеро. Цифровая обработка многомерных сигналов 1984, 488 с.
2. Василенко Г.И., Тараторин А.М. Восстановление изображений. М., Радио и связь, 1986.
3. Е.З. Демиденко. Линейная и нелинейная регрессии. - Москва, "Финансы и статистика", 1981, 303 с.
4. V.A. Fursov. Identification of optical distorting systems with selecting image informative fragments. Workshop on Digital Image Processing and Computer Graphics. Proceedings SPIE. - 1994. -. 2363, p 62-68.
5. В.А. Фурсов. Идентификация оптических искажающих систем с отбором информативных фрагментов изображений. Компьютерная оптика. Вып. 14-15, 1995, с. 78-79.
6. Ф.Р. Гантмахер. Теория матриц. М.: Наука. 1967, 575 с.
7. В.А. Фурсов. Восстановление изображений КИХ-фильтрами, построенными путем непосредственной идентификации инверсного тракта. Компьютерная оптика. Вып. 16, 1996.
