Авторы:   Samia Ainouz¹, Jihad Zallat² и Fabrice Meriaudeau¹

Laboratoire d’Electronique, Informatique et Image (LE2I), UMR CNRS 5158, Burgundy University
¹ LE2I, 12 rue de la fonderie, 71200, Le Creusot, France
² LSIIT, Bd Sébastien Brant, 67400, Strasbourg, France

16^-я Европейская конференция по обработке сигналов (European Signal Processing Conference EUSIPCO 2008), Швейцария, Лозанна, 25–29 августа 2008 г. Публикация EURASIP.

перевод: отдел внешнего обслуживания РНБ и Пасяда А.В. 21.2.2012

КРАТКИЙ ОБЗОР

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

1. Введение
2. Получение поляризационного изображения
  2.1. Получение стоксовых изображений
  2.2. Шум в стоксовых изображениях
3. Параметрическая оценка шума
  3.1. Фрактальный алгоритм фильтрации с преобразованием в векторную форму (FVFA)
  3.2. Применение к изображениям в градациях серого
4. Фильтрация поляриметрических изображений
  4.1. Фильтрация на основе поляризации
  4.2. Иллюстрации и обсуждение
5. Заключение
Литература
Переменые в статье

1. ВВЕДЕНИЕ

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

Системы обработки изображений, которые могут измерять состояние поляризации проходящего через сцену света главным образом основаны на возможности создать эффективный анализатор состояния поляризации перед камерой (датчиком), чтобы получать вектора Стокса [1]. Это поляриметры Стокса создают четыре изображения, называемые «стоксовыми изображениями», которые соответствуют четырём параметрам Стокса. Соответственно поляризационно кодированные изображения имеют многомерную структуру, т. е. с каждым пикселем изображения связана многокомпонентная информация. Кроме того, информационное содержание поляризационно кодированных изображений сложно смешано в поляризационных каналах, что затрудняет их правильную интерпретацию в случае присутствия шума.

Шум присутствует в любых формирующих изображения системах и следовательно есть и на стоксовых изображениях. Шум носит аддитивный характер, когда сцена освещается некогерентным светом, и мультипликативный, когда освещение когерентно [3], [4]. Его присутствие ухудшает возможность интерпретировать данные и мешает анализировать физический потенциал поляриметрической информации. Немногие работы в литературе обращаются к вопросу фильтрации поляриметрических изображений. Тем не менее мы отмечаем использование оптимизационных методов в [5] для того, чтобы оптимизировать параметры системы обработки изображений по критерию отношения сигнала к шуму, или рассмотренное в публикации [3] повышение точности степени поляризации с целью снижения шума на стоксовых изображениях.

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

В этой статье, учитывая, что изображения получаются при некогерентном освещении, аддитивный шум оценивается и исключается так, что физическое содержание поляриметрических изображений сохраняется, насколько это возможно. Методы оценки выбраны в области многоспектральных методов фильтрации. Чтобы воспользоваться преимуществами обоих методов, алгоритм фильтрации использует комбинацию методов графика рассеяния [6] и маскировки данных [4]. Так как информационное содержание поляризационных изображений сложно комбинируется в нескольких поляризационных каналах, к зашумлённому изображению применяется обход фрактальной кривой Пеано-Гильберта для того, чтобы сохранить связанность однородных зон и свести к минимуму воздействие выбросов. Характеристики и систематическая погрешность метода статистически исследованы бутстреп-методом [7].

2. ПОЛУЧЕНИЕ ПОЛЯРИЗАЦИОННОГО ИЗОБРАЖЕНИЯ

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

2.1. Получение стоксовых изображений

Общее состояние поляризации световой волны можно описать так называемым вектором Стокса S, который полностью характеризует усреднённые во времени поляризационные свойства излучения. Оно определяется следующим сочетанием комплекснозначных компонентов E_x и E_y вектора напряженности электрического поля вдоль двух ортогональных направлений x и y по формуле из [1]:

S =	S0	=	ExEx* + EyEy*	(1)
	S1		ExEx* – EyEy*
	S2		2 Re ExEy*
	S3		2 Im ExEy*

Легко увидеть, что

Это условие известно как физическое условие формализма Стокса. Произвольный вектор, который не удовлетворяет этому условию, не является вектором Стокса и не имеет какого-либо физического значения.

Общая схема получения стоксовых изображений показана на Рис.1 [1]. Волна, отражённая от мишени, представленная вектором Стокса S_вх, анализируется анализатором состояния поляризации, измеряя её проекции на 4^-х линейно независимых состояниях. Анализатор состоит из линейного поляризатора (ЛП) и четвертьволнового компенсатора (К), вращающегося по четырём углам (θ_i)_i=1,2,3,4. Камера измеряет входящие интенсивности. Полный набор из 4 измерений можно записать в векторной форме как:

I представляет измеренную камерой матрицу интенсивности 4×1. Вектор Стокса S_вх может быть легко получен из матрицы необработанных данных I, при условии, что матрица модуляции анализатора A известна по калибровке. Для идеального случая (теории) матрицу A можно представить как:

A(θi) =	1	–cos2(2θi)	–½ sin (4θi)	sin2(2θi)	(4)
	–1	cos2(2θi)	½ sin(4θi)	–sin2(2θi)
	0	0	0	0
	0	0	0	0

Углы θ_i выбраны так, что матрица A является обратимой. Стоксово изображение определяется как измеренные значения параметров вектора Стокса, связанные с каждым пикселем.

Рис. 1 — Устройство для получения стоксовых изображений.

На Рис.2 показан пример стоксова изображения, полученного в условиях при наличии шума для проверки действия нового алгоритма. Изображение сделано из 4 небольших образцов A, B, C и D, наклеенных на картон. Объекты A и D прозрачны.

Рис. 2 — Стоксово изображение четырех небольших образцов, наклеенных на картон.

2.2. Шум в стоксовых изображениях

Шум, достигающий поляриметрических систем формирования изображения, освещаемых некогерентным источником света, является аддитивным и независимым [4]. Этот тип шума можно моделировать нулевым средним произвольным гауссовым распределением, плотность распределения вероятности которого выражена следующим образом [6]:

f_x(x) = exp            (5)

где σ_ш² — дисперсия шума. Влияние аддитивного шума n_a на цифровое изображение g в пикселе с координатами (i, j) выражается как сумма изображения I, свободного от шума, и шума в форме

В условиях идеального получения [изображения] вектор Стокса рассчитывают по каждому пикселю из уравнения (3) следующим образом

В присутствии шума эта формула приобретает следующий вид:

Рассчитываемый вектор Стокса Ŝ представляет собой независимую сумму из теоретического вектора Стокса S и члена δS, отвечающего за воздействие аддитивного шума. Затем, чтобы восстановить зашумлённое изображение δS и исключить его из вычисляемого стоксова изображения Ŝ, потребуется оценка распределения шума. Однако непосредственная фильтрация может внести нефизический смысл фильтрованного стоксова изображения Ŝ – δS. Это означает, что для значительного числа пикселей вектор Ŝ – δS не удовлетворяет условию, сформулированному в уравнении (2) и, следовательно, не может рассматриваться как вектор Стокса. Без этого условия стоксовы изображения не представляют для нас никакого интереса. По этим причинам будут предприняты дополнительные шаги для достижения компромисса между фильтрацией и физическим смыслом для как можно большего числа пикселей.

3. ПАРАМЕТРИЧЕСКАЯ ОЦЕНКА ШУМА

Так как рассчитываемое стоксово изображение с помехами представляет собой независимую сумму шума и изображения без шума, оценки распределения δS достаточно для того, чтобы иметь информацию об аддитивном шуме. В контексте многомерности стоксовых изображений в данной работе нас интересуют два многоспектральных метода фильтрации: метод графика рассеяния (SP) [6] и метод маскировки данных (DM) [4]. Предложенный здесь алгоритм фильтрации получает выгоду от обоих методов.

Чтобы исключить воздействие отклоняющихся точек, изображение преобразовано во фрактальную кривую Пеано-Гильберта. Новый метод проверен на полутоновом изображении, и результаты сравнены с методами графика рассеяния и маскировки данных. Затем этот метод применялся для фильтрации поляриметрических изображений.

3.1. Фрактальный алгоритм фильтрации с преобразованием в векторную форму (FVFA)

Метод графика рассеяния [6] не пригоден для работы со значительно текстурированными изображениями. Кроме того, этот метод эффективен в отношении времени благодаря расчёту локальной статистики внутри изображения. Метод маскировки данных [4] основан на оценке параметров шума на остаточном изображении, полученном с помощью свёртки зашумлённого изображения с ядром лапласиана. Но из-за оставшихся контуров остаточное изображение также содержит полезную информацию, ведущую к завышению параметров шума.

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

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

L =	1	–2	1	(9)
	–2	4	–2
	1	–2	1

Оставшиеся контуры выявляются посредством применения детектора выделения контуров Собела [8] к обработанному лапласианом свернутому изображению. Заданное пороговое значение отсекает пиксели, которые считаются принадлежащими к оставшимся контурам [4]. Остаточное изображение без контуров представляет изображение, обработанное лапласианом. Это изображение преобразуют в вектор, следуя по пути фрактальной кривой, как показано на примере изображения 9×9 на Рис.3. Используется обход растра по фрактальной кривой Пеано-Гильберта, и вектор получаем так же, как нумерация на Рис.3 от 1 до 81.

Рис. 3 — Обход растра путём фрактальной кривой Пеано-Гильберта.

Локальные средние значения μ и среднеквадратические отклонения σ рассчитываются по результирующему вектору с использованием интервала сдвига (этот интервал представляет на растре размер текущего апертурного окна, в котором ведутся вычисления). Локальное среднее значение и несмещённое среднеквадратическое отклонение определяются, соответственно, по следующим формулам:

μν(i) = ν(i + k)	(10)
σν(i) = (ν(i + k – μν(i))2½

где ν — часть вектора обработанного лапласианом изображения, ограниченного интервалом сдвига величиной 2m + 1 единиц растра. Эта локальная статистика — векторная версия рассчитанного алгоритмом по методу графика рассеяния, примененного к остаточному изображению без контуров.

Затем плоскость (μ, σ) наносится на график. На этой плоскости формируются два типа точек [6]. Первый тип — плотное облако, а второй — выделенный набор точек, соответствующих пикселям оставшимся контуров, которые не исключены фильтрацией по Собелу. Пересечение линейной регрессии точек облака и оси y даёт лучшую оценку среднеквадратического отклонения шума. Аналогично, среднее значение шума может быть оценено путём применения тех же инструкций к плоскости, нанесённой на график (μ, σ).

Фрактальная кривая имеет два преимущества: она облегчает задачу расчёта локальной статистики внутри изображения и сохраняет не более чем окрестность пикселей изображения. Кроме того, преобразование изображения в векторную форму так рассеивает выделенные точки плоскости (μ, σ), что это препятствует их учёту во время процесса регрессии. Это снижает завышение параметров шума, которым отличается метод маскировки данных.

3.2. Применение к изображениям в градациях серого

Чтобы доказать эффективность предложенного алгоритма, на двух полутоновых изображениях было проведено три эксперимента (Рис.4). На первом изображении текстурированное изображение шаблона в клетку, а на втором — срез трёхмерного изображения головы, полученного методом обработки магнитно-резонансных изображений. Оба изображения состоят из 256×256 пикселей. К этим двум изображениям добавлены аддитивные шумы нулевого среднего с дисперсией 5, 10, 30 и 60. Дисперсия шума оценивается затем тремя методами: методом графика рассеяния, методом маскировки данных и предложенным алгоритмом (FVFA). Чтобы получить статистически правильные результаты, расчётное среднеквадратическое отклонение σˆ берётся как эмпирическое среднее для 200 оценок, найденных бутстреп-методом передискретизации [7]. Затем рассчитывается дисперсия σ_d, которая показывает, что точное значение среднеквадратического отклонения лежит в интервале [σˆ – 3σ_d, σˆ + 3σ_d] с вероятностью 99,73% [прим. перев. — приведён статистический критерий, известный как "правило трёх сигм"].

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

Результаты расчётов обобщаются в Табл.1 для Рис.4 (а) и в Табл. 2 для Рис.4 (б). Дисперсия σ_d указана в скобках.

Табл. 1 — Сравнение параметров оценки шума с помощью методов графика рассеяния, маскировки данных и фрактального алгоритма фильтрации FVFA для изображения шахматных клеток с Рис.4 (а).

Табл. 2 — Сравнение параметров оценки шума с помощью методов графика рассеяния, маскировки данных и фрактального алгоритма фильтрации FVFA для изображения с Рис.4 (б).

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

4. ФИЛЬТРАЦИЯ ПОЛЯРИМЕТРИЧЕСКИХ ИЗОБРАЖЕНИЙ

4.1. Фильтрация на основе поляризации

Как видно из параграфа 2, измеренный вектор Стокса Ŝ, связанный с пикселем (i, j), получен как независимая сумма идеального вектора Стокса S_p и связанной с ним шумовой составляющей δS по формуле:

Составляющая δS равна A^-1δI (δI = n_a). В соответствии с гауссовым допущением аддитивного шума, параметры шума рассчитывают рассматриваемым методом FVFA. Восстанавливается интенсивность шума δI, а затем рассчитывается поляриметрический шум δS.

Конечно, реальный вектор Стокса S_p можно вывести из уравнения (11) следующим образом:

Кроме того, информационность наших изображений сильно предопределена физическим содержанием, т.е. вектор S_p должен удовлетворять уравнению (2). Однако прямое применение уравнения (12) к фильтрованному поляриметрическому изображению вызывает нефизическое поведение большого числа пикселей изображения.

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

Вектор Стокса S_p должен удовлетворять уравнению (2), которое эквивалентно формуле:

Объединив уравнения (13) и (14), необходимо найти параметр α, который удовлетворяет условию:

Раскрывая уравнение (15), необходимо найти параметр α, который удовлетворяет неравенству:

Если предположить, что a = δS^T GδS, b = δS^T GŜ – Ŝ^T GŜ и c = Ŝ^T GŜ, то уравнение (16) можно написать в упрощенном виде следующим образом:

Два действительных решения квадратного уравнения, как известно:

Предположим, что α₁ всегда соответствует наименьшему решению, а α₂ — наибольшему решению. Для бесконечно малой величины ε, такой, что, если α₁ положительно, α₁– ε всё ещё положительно, и после классического разрешения неравенства (17) возникают три случая, в зависимости от знака перед a:

• sign(a) > 0
  Если 1 ∈ [α₁, α₂], то α = 1,
  Если 1 ∈ [α₁, α₂] и α₁ > 0, то α = α₁ – ε,
  В противном случае между 0 и 1 нет решений α,
  
  
• sign(a) < 0
  Если 1 ∈ [– ∞, α₁] ∪ [α₂, + ∞], то α = 1,
  Если 0 < α₂ < 1, то α = α₂ – ε,
  В противном случае между 0 и 1 нет решений α,
  
  
• sign(a) = 0
  Существует одно решение α₁ = α₂ = c/b
  Если 0 < α₁ < 1, то α = α₁ – ε,
  Если α₁ > 1, то α = 1,
  В противном случае между 0 и 1 нет решений α.

Эти три случая учитывают тот факт, что значение α должно принадлежать промежутку [0,1].

Наконец, алгоритм можно обобщить следующим образом:

1. Создайте шумовую составляющую δS, связанную с каждым пикселем изображения, методом FVFA
2. Рассчитайте уравнение (10) для каждого пикселя
3. Если вектор Стокса, относящийся к пикселю, физически осуществим, задайте значение α равным 1
4. В противном случае, найдите параметр α, следуя ранее изложенным инструкциям
5. Если дискриминант Δ отрицателен, или если в промежутке [0,1] нет действительного значения α, то выберите вектор Стокса, удовлетворяющий максимуму между Ŝ^TGŜ и S_p^TGS_p.

4.2. Иллюстрации и обсуждение

Наш алгоритм действует на два изображения: на синтезированное стоксово изображение и действительное стоксово изображение. Синтезированное изображение создано следующим образом: вектор Стокса в центре изображения привязывается к S = [1 1/√3 1/√3 1/√3] (белая округлая часть), а в противном случае — к нулю. Стоксовы изображения умножаются на матрицу модуляции A, как в уравнении (3), чтобы получить соответствующие каналы интенсивности. Относящийся к значениям интенсивности изображения, гауссиальный шум с нулевым средним и дисперсией 0,2 добавляется к каналам интенсивности. Изображения с шумом инвертируются, как в уравнении (7), чтобы получить стоксово изображение с шумом (Рис.5).

Рис. 5 — Каналы Стокса с шумом. К соответствующим каналам интенсивности на изображениях добавляется гауссиальный шум с нулевым средним и дисперсией 0,2.

Дисперсии, рассчитанные с помощью алгоритма фильтрации FVFA на четырех каналах интенсивности, соответствующих S₀, S₁, S₂ и S₃, соответственно равны 0,18; 0,19; 0,194 и 0,187. Эти значения очень близки к синтезированной дисперсии. Результаты также показывают, что шум, воздействующий на четыре поляриметрических измерения, примерно тот же самый. Параметр регуляризации α рассчитан для каждого пикселя. На Рис.6 показаны бинарные значения этого параметра. Пиксели, для которых значение α найдено между 0 и 1, отмечены белым, в противном случае они отображаётся чёрными пикселями.

Рис. 6 — Бинарные значения параметра регуляризации α.

Использование параметра регуляризации позволяет обработать значительное количество пикселей на изображении. Действительно, без этого параметра (уравнение (12)) изображение хорошо фильтруется, но только 10% пикселей в S_p имеют физический смысл, тогда как с применением параметра α (уравнение 14) — 70% пикселей. Количество оставшихся пикселей, не имеющих физического смысла, не существенно. Это связано с тем фактом, что условия алгоритма от 1 до 4 не удовлетворяются, и ни Ŝ^TGŜ и S_p^TGS_p не являются положительными (этап 5). Результат двух фильтраций показан на Рис.7. Очевидно, что полностью фильтрованное изображение лучше, но физически параметр регуляризации предпочтительнее, и фильтрация в этом случае всё ещё допустима.

а	б

Рис. 7 — Фильтрованные стоксовы изображения с: (а) параметром регуляризации (уравнение 10), (б) без параметра регуляризации (уравнение 12).

Как показано на Рис.7. (a), наш алгоритм обеспечивает улучшение поляриметрической информации, которую несёт стоксово изображение, сохраняя её физическое ограничение для большей части пикселей.

В случае реального измерения на Рис.2 используется стоксово изображение. Как видно из этого изображения, много информации потеряно, особенно в каналах S₁, S₂ и S₃. Действительно, поляриметрическая информация лежит на трёх последних каналах. Характерные неустранимые дисперсии шума, рассчитанные по соответствующим изображениям интенсивности, доходят соответственно до следующих значений: 10,75; 10,52; 9,78 и 9,88. Эти значения существенны, в сравнении с синтезированными изображениями. Это связано с тем фактом, что значения интенсивности находятся в диапазоне между 0 и 255 для действительного изображения и в диапазоне от 0 до 1 для синтезированного изображения. Результаты физической фильтрации на зашумлённых стоксовых изображениях с шумом представлены на Рис.8. Новая фильтрация обеспечивает 64% физических пикселей, тогда как классическая фильтрация обеспечивает лишь 7% физических пикселей. Предложенный алгоритм, таким образом, является компромиссом между полностью фильтрованным изображением и физически ограниченным (для большей части пикселей) изображением.

а	б

Рис. 8 — (а) Стоксово изображение с шумом, (б) Фильтрованное стоксово изображение.

5. ЗАКЛЮЧЕНИЕ

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

ЛИТЕРАТУРА

[1] R.A. Chipman, "Polarimetry" изд. в Handbook of Optics, M. Bass, ed (McGraw-Hill, New York), стр. 22.1-22.36, 1993.

[2] J. Zallat, P. Graebling, Y.Takakura, "Using polarimetric images for material classification," изд. в ICIP’03, Spain,, 2003.

[3] A. Bénière, F. Goudail, M. Alouini, D. Dolfi, "Precision of degree of polarization estimation in the presence of additive Gaussian detector noise," Optics Communication, том 278, стр. 264–269, Oct. 2007.

[4] A. B. Corner, R. M. Narayanan, S. E. Reichenbach, "Noise estimation in remote sensing imagery using data masking," изд. в Int. J. Remote sensing, том 4, 2003, стр. 689–702.

[5] J. Zallat, S. Ainouz, M-P. Stoll, "Optimal configuration for imaging polarimeters : impact of image noise and systematic errors," изд. в J. Opt., том 8, 2006, стр. 807–814.

[6] B. Aiazzi, L. Alparone, A. Barducci, S. Baronti, I. Pippi, "Estimating noise and information of multispectral imagery ," изд. в Opt. Eng, том .41, 2002, стр. 656–668.

[7] R. C. H. Cheng, "Bootstrap methods in computer simulation experiments," изд. в Simulation Conference Proc., том 2, 1995, стр. 171–177.

[8] N. Kazakova, M. Margala, N. G. Durdle, "Sobel edge detection processor for a real-time volume rendering system," изд. в ISCAS’04, том 2, 2004, стр. 913–916.

ПЕРЕМЕННЫЕ В СТАТЬЕ