УДК 004.92 DOI: 10.46960/1816-210X2022 226
АДАПТИВНОЕ УПРАВЛЕНИЕ МЕДИАННОЙ ФИЛЬТРАЦИЕЙ ИЗОБРАЖЕНИЙ В MATHCAD
Е.А. Никулин ORCID: 0000-0002-2520-681X e-mail: nea@mail.nnov.ru
Нижегородский государственный технический университет им. Р.Е. Алексеева
Нижний Новгород, Россия
Представлен обзор существующих в пакете MathCAD 15 средств генерации шумов и их фильтрации в растровых изображениях. Разработаны новые методы медианной фильтрации: с произвольными размерами прямоугольного окна и с оптимальными размерами квадратного окна, адаптивно подбираемыми в каждом пикселе изображения. Предложенные алгоритмы и реализующие их функции скомпонованы в законченные MathCAD-программы, очищающие изображения от различных помех с разными качественными показателями. Показана особая эффективность алгоритма адаптивной медианной фильтрации с настраиваемыми размерами окна.
Ключевые слова: изображение, шум, фильтрация, медиана, адаптивный алгоритм, MathCAD.
ДЛЯ ЦИТИРОВАНИЯ: Никулин, Е.А. Адаптивное управление медианной фильтрацией изображений в MathCAD // Труды НГТУ им. Р.Е. Алексеева. 2022. № 2. С. 26-35. DOI: 10.46960/1816-210X_2022_2_26
ADAPTIVE CONTROL OF MEDIAN FILTERING OF IMAGES IN MATHCAD
Е.А. Nikulin
ORCID: 0000-0002-2520-681X e-mail: nea@mail.nnov.ru
Nizhny Novgorod state technical university n.a. R.E. Alekseev Nizhny Novgorod, Russia
Abstract. Overview of existing in the MathCAD 15 package means of noise generation and filtering thereof in bitmap images, is presented. New median filtering methods are developed: with arbitrary rectangular window sizes and with optimal square window sizes adaptively selected in each pixel of the image. Proposed algorithms and functions implementing thereof are compiled into complete MathCAD programs that clean images from various interferences with different quality indicators. Special efficiency of the adaptive median filtering algorithm with adjustable window sizes, is shown.
Key words: image, noise, filtering, median, adaptive algorithm, MathCAD.
FOR CITATION: Nikulin E.A. Adaptive control of median filtering of images in MathCAD. Transactions of NNSTU n.a. R.E. Alekseev. 2022. № 2. Рр. 26-35. DOI: 10.46960/1816-210X_2022_2_26
Введение
Одной из важнейших задач компьютерной обработки изображений, представленных в виде матриц монохромных либо цветных пикселов, является фильтрация шумов, искажающих цветовые коды, в основном, случайным образом, по различным вероятностным законам. Помехи можно условно разделить на импульсные (искажены отдельные пикселы), линейные (царапины в виде цепочек пикселов) и площадные (сплошные пятна пикселов). Задача компьютерной фильтрации шумов в изображениях, основательно изученная в [1, 2], заключается в обнаружении и устранении этих помех с помощью алгоритмов и реализующих их программ. Наиболее популярным, но не самым эффективным, размывающим границы объектов, является линейный метод усреднения в маске, скользящей по матрице изображения.
© Никулин Е.А., 2022
Настоящая статья, посвященная решению задачи фильтрации в программной среде МаШСАО [3], кратко знакомит с необходимыми встроенными функциями обработки изображений. Раскрываются некоторые эффекты их использования и более обстоятельно представляются усовершенствованные автором нелинейные алгоритмы медианной фильтрации, дающие подчас удивительные результаты очистки сильно зашумленных изображений. Моделью исследований служит монохромный фрагмент портрета знаменитой Джоконды.
Генерация шумов в МаШСАБ
Рассмотрим способы управляемого добавления пиксельных шумов в растровые изображения, после чего те могут служить основой для испытания методов шумоподавления в отсутствие реальных зашумленных рисунков.
Самой первой в категории «Обработка изображений» панели «Вставка функции» стоит функция аёёпо1ве (М,р,£), которая с равными вероятностями р/2 увеличивает либо уменьшает код каждого элемента матрицы М точно на £ единиц. При выходе измененного кода из байтового интервала Л=[0, 255] он самой функцией автоматически сдвигается на его ближайшую границу, тем самым исключая сдвиг кодов шаблоном рисунка в Л на ±256. Добавление импульсного шума с разными вероятностями р и фиксированными смещениями £ реализовано в верхнем ряду изображений на рис. 1. Они демонстрируют уплотнение импульсных помех с ростом вероятности их появления в половине пикселов при р=0.5 и в каждом при 2р=1, а также повышение контраста зашумленных точек с окружающим фоном при увеличении параметра разброса от ¿=50 до удвоенного 2^=100.
В отсутствие в МаШСАО других встроенных функций зашумления матричных изображений создадим пользовательские функции:
иш^с, б) := гишГ(1, с - б , с + б)0, погш(е, б) := гпогш(1, с, б)0 (1)
переменных смещений числа с, использующие встроенные функции генерирования случайных векторов: гишГ для равномерно распределенных элементов векторов в интервале (с-£,с+£) и гпогш для нормально распределенных элементов с матожиданием (МО), равным с, и среднеквадратичным отклонением (СКО), равным По свойству нормального распределения 99,73 % смещенных чисел лежат вокруг значения с в интервале (с-3^,с+3^).
Анализ изображений во втором ряду рис. 1 подтверждает вполне ожидаемые эффекты расположения случайных пиксельных шумов и формирования их интенсивностей:
• плотности равномерно и нормально распределенных пиксельных помех при одинаковых значениях параметра разброса £ заметно отличаются друг от друга из-за практически втрое большей ширины интервала расположения гауссовских случайных чисел в сравнении с интервалом локации равномерных;
• цвет большинства точечных помех в матрицах, полученных векторизацией функций (1) от матрицы М, оказался противоположным цвету окружающих фоновых пикселов по принципу «соль и перец» — светлые точки на темном фоне и темные - на светлом.
При этом с удвоением параметра разброса £ такой контрастный шум почти полностью «забил» исходное изображение. В реальности выводимые матрицы таких помех не имеют, а «перец» на светлом челе Моны Лизы и белая «соль» на ее темных волосах - всего лишь иллюзия. Такие эффекты инверсии цвета пикселов объясняются тем, что в них зашумленные матрицы имеют внебайтовые коды, которые шаблоном рисунка лишь визуально конвертируются в интервал Л путем параллельного сдвига на ±256.
Поручим приведение случайных чисел, сгенерированных в (1), в байтовый диапазон Л следующим функциям:
иш(с, б) := шт(шах(ишГ(с, б) , 0) , 255), пог(с, б) := шт(шах(погш(с, б) , 0) , 255). (2)
В нижнем ряду рис. 1 в изображениях, зашумленных функциями (2), эффекты конвертации, присущие шаблону вывода рисунка, исчезли.
Рис. 1. Импульсный шум в изображениях Fig. 1. Impulse noise in images Функции линейной фильтрации изображений в MathCAD
Переходим к решению задачи фильтрации шума с помощью встроенных функций MathCAD, собранных в категории «Обработка изображений» панели «Вставка функции» и предназначенных для сглаживания импульсных помех в монохромных изображениях. Методы подавления помех разделяются на пространственные и частотные [1, 2], но из-за ограниченности возможностей протестировать здесь все функции будут рассмотрены лишь некоторые из них, выполняющие пространственную фильтрацию.
Как в радиотехнике для устранения пульсаций сигнал пропускается через фильтр низких частот, так и в обработке изображений, зашумленных пиксельными помехами, используется обладающий интегрирующим свойством метод усреднения яркостей пикселов в окружающем обрабатываемую точку прямоугольном окне (ядре, маске, апертуре) WeR х нечетных размеров H и D с полями шириной сверху и снизу по h=(H-1)/2, а слева и справа по
окне
(2й+1)х(2й?+1)
^=(П-1)/2. Общее свойство всех масок, сохраняющее яркости пикселов однотонного фона, -сумма всех элементов Жух равна единице. Яркость (/,у')-го пиксела матрицы МеК(т+1)х(и+1) в ), центрированном в элементе Ыу, рассчитывается по формуле:
тш{й,ш-г} ■ ,, .л ' тт{а ,п-ц}
^ £ЖуМуЦ+х [0,ш],7£[0,п]. (3)
у=тах{-й,-г} х~
=
Результирующая матрица свертки 8 имеет такие же размеры, как и матрица М. Указанные в (3) пределы сумм предотвращают при неполном покрытии окна выход индексов элементов на краях матрицы М за пределы интервалов [0,ш] и [0,п], обеспечивая выполнение неравенств 0< /+у < ш и 0<Ц+х< п. При полном покрытии матрицей окнау=-И...И и х=-й...й.
На рис. 2 продемонстрирована работа имеющихся в МаШСАО функций сглаживания со встроенными в них масками
Ж =1
1 1 1 1 1 1 111
Ж5 =—
2
25
1 1 1 1 , Жо3 =— 16 1 2 1
1 2 1 2 4 2
10 1 1 1 1 2 1
Ж* =—
2
32
11111 112 11 1 2 4 2 1 112 11 11111
• иш8шооШ(М) и иш8тооШ5(М) с однородными масками Ж3 и Ж5;
• сеп18тооШ(М) с центрально-взвешенным 3х3-окном Жс, вносящим вдвое больший вес цвета сглаживаемого пиксела Мц в общую сумму (3) с целью уменьшения размытия границ полезных фигур изображения;
• огШоБтоо^М) и ог1Ьо8тоо1Ь5(М) с ортогонально-взвешенными масками Жо3 и Жо5, также уменьшающими размытие и расфокусировку кромок фигур в сравнении с однородными окнами тех же размеров.
В нижнем ряду рис. 2 показано применение функций сопуо12ё(М,Ж), еопуо1уе3(М, Жс3) и сопуо1уе5(М,Жс5) фильтрации с внешне задаваемыми окнами: произвольных разме-
ров {Н, П} и квадратными масками фиксированных размеров ЖсзеК3х3 и Жс5еЯ5х5.
Работа перечисленных функций оконной фильтрации продемонстрирована на примере изображения М=«Царап.Ьтр», созданного в графическом редакторе путем нанесения на рисунок «Мопа.Ьтр» линейных (царапин разных толщин и направлений) и площадных (пятен разных размеров и форм) помех. Последующие добавления в матрицу М пиксельных шумов функциями аёёпо1ве(М,0.5,5), ип1{М,5) и пог(М,5) с параметром разброса 5=50 сформировали матрицы Ма, Ми и Мп. Анализ этих изображений и результатов их фильтрации привел к следующим выводам:
• из трех типов пиксельного шума наименее плотным является фиксированный, равновероятно добавленный в матрице Ма к элементам Мц как 0 либо ±50; из оставшихся менее зернистым и контрастным выглядит равномерно распределенный в интервале (-50,50) шум в матрице Ми, а наиболее контрастным - нормально распределенный с 99,73 % интервалом (-150,150) гауссовский шум в матрице Мп;
• все методы оконного сглаживания обладают присущим им свойством размытия как отдельных пикселов, так и границ широких царапин и пятен, расширяющим резкие перепады яркостей в пологие блеклые площадки и полоски; размытие границ проявляется тем сильнее, чем больше размеры окна фильтрации;
• с увеличением размера окна появляется и далее расширяется черная окантовка изображения, что объясняется меньшим количеством усредняемых пикселов, накрытых окном на краях рисунка;
• в нижнем ряду рис. 2 показана возможность применения не встроенных в функции свертки, а внешне задаваемых масок. Однострочное Ж=[1 24842 22 и одностолбцовое
Жт окна избирательно сглаживают функцией сопуо12ё вертикальные и горизонтальные
9
тонкие царапины. Там же видно действие двух квадратных окон: крестовидного Ж3 и квадратного Ж5, усредняющего 24 числа, окружающих Му.
Рис. 2. Встроенные функции линейной фильтрации Fig. 2. Built-in Line Filtering Functions
Медианная фильтрация изображений
Возможности оконной фильтрации изображений не исчерпываются линейным алгоритмом (3). Среди других, нелинейных операций с яркостями пикселов, попавших в скользящее по рисунку окно, можно назвать выборку минимального, максимального, средних квадратичного, геометрического, гармонического и других математических величин, например, медианы - срединного по расположению числа в упорядоченной по возрастанию последовательности чисел.
Работа встроенного медианного фильтра теёй11;(М), предназначенного для очистки монохромных изображений М от импульсных (пиксельных) помех, основана на поочередном совмещении центра квадратного окна фиксированного размера 3х3 с каждым пикселом рисунка Му, сортировке яркостей покрытых окном пикселов и запоминании в элементе возвращаемой матрицы медианы этой последовательности.
Желая реализовать свойства ориентированной (горизонтальной или вертикальной) линейной фильтрации, присущей функции сопуо12ё, создадим универсальный аналог медианного фильтра матрицы М с окном произвольных размеров НхП, где Н и П — нечетные (1, 3, 5, 7,...) габариты окна с полями вокруг центральной точки шириной сверху и снизу й=(Н-1)/2, а слева и справа ^=(П-1)/2 пикселов. В приведенной ниже функции швй/ИЖНхП-окно, сканирующее все пикселы матрицы М, можно задавать не только квадратным, но и как в функции сопуо12ё, прямоугольным, в том числе, однострочным (Н=1, И=0) или одностолбцовым (П=1, ^=0). Это позволяет фильтровать не только импульсные, но и координатно ориентированные сплошные протяженные помехи типа царапин и трещин, а именно: горизонтальные царапины рекомендуется сглаживать вертикальными Нх1-масками, а вертикальные - горизонтальными масками размерами 1 хП.
medfilW(M, H, D) :=
"Медианный фильтр в окне размером H*D"
[n^ cols(M) - 1 h^ (H- 1) * 2 d^ (D- 1) * 2]
for i e 0.. m ^ rows (M) - 1
(y ^ max(0, i - h) Y ^ min(m, i + h)) for j e 0.. n
(x ^ max(0,j - d) X ^ min(n,j + d)) Si,j ^ median( submatrix(M, y, Y, x, X)) S
Сверху и справа на рис. 3 показаны изображение N полученное импульсным нормально распределенным зашумлением поцарапанного рисунка М=«Царап.Ьтр» функцией пог(М,5) с СКО 5=50 и результаты его медианной фильтрации МаШСАВ-функцией теёйк в окне 3х3 и универсальной функцией швй/ИЖ в квадратных окнах 3х3^9х9. Очевидно, что в 3х3-окне обе функции работают одинаково. Во втором и третьем рядах сконструированная функция швй/ИЖ испытана на неквадратных окнах.
Сопоставление размеров, плотности, яркости и контраста оставшихся на изображениях шероховатостей позволяет сделать новые выводы и подтвердить ранее сделанные предсказания о свойствах медианной фильтрации монохромных изображений окнами разных размеров и конфигураций:
• импульсные (точечные) помехи расплылись в пятна, размер которых примерно совпадает с размером окна фильтрации; что касается наблюдаемого контраста пиксельных помех с окружающим фоном, с ростом размера окна от 3х3 до 9х9 он неуклонно уменьшается вплоть до еле различимого;
,3,1) ше*Ш№(>Г ,5,1) гонШПУф ,7,1) те<1Ш\¥(Ы ,9,9)
Рис. 3. Универсальная медианная фильтрация
Fig. 3. Universal Median Filtering
• это же свойство размытия точечных помех приводит к расплыванию и ухудшению четкости границ не только площадных шумов (пятен), но и полезных фигур, составляющих сюжет изображения; при этом наиболее сильно размываются границы, ориентированные вдоль максимального размера окна;
• медианная фильтрация справляется с возложенной на нее задачей качественнее, чем функции линейного сглаживания, работа которых показана на рис. 2; особенно наглядно ее преимущества проявляются при устранении шумов в виде отдельных точек и линий толщиной менее половины максимального размера окна (после сортировки экстремальные яркости пикселов окна группируются на краях числовой последовательности («перец» в начале, «соль» в конце), оставляя в середине (медиане) яркости фоновых пикселов, благодаря чему контраст зашумленного пиксела в центре окна существенно понижается);
• горизонтальные царапины хорошо сглаживаются вертикально вытянутыми окнами фильтра, а вертикальные - горизонтальными окнами; наклонные и толстые штрихи следует устранять медианным фильтром с квадратным окном;
• помехи в форме пятен размерами более, чем тах{Н,П}, фильтруются плохо, сглаживаются только их контрастные границы.
Адаптивная медианная фильтрация изображений
Приведенные примеры медианной фильтрации зашумленного изображения M как встроенной функцией medfilt(M с фиксированным 3х3-окном, так и универсальной функцией medfilW(M,H,D) с окнами произвольных размеров HxD демонстрируют разное качество удаления шумов различного типа и ставят закономерные вопросы о выборе оптимального размера окна и правомочности применения одного и того же размера на всех участках изображения, отличающихся разными яркостями, плотностями и типами помех. Ответом на эти вопросы является усовершенствование предложенного в [1] алгоритма адаптивной медианной фильтрации изображения MeR(m+1)x(n+1), подбирающего в каждом пикселе My с индексами ie[0, m] и je[0, n] наиболее подходящий размер квадратного (2^+1)х(2^+1)-окна, начиная с минимальной ширины всех полей вокруг его центральной точки d=1. Перед алгоритмом ставятся следующие цели:
• удалить бинарные (черные и белые) импульсные помехи и максимально сгладить шумы с другими функциями распределения (равномерными, нормальными и т. п.);
• минимизировать размытие границ протяженных объектов, составляющих основу сюжета изображения.
Реализующая эти задачи функция
medfilA(M) := "Адаптивный медианный фильтр"
D ^ floor(min(m ^ rows(M) - 1, n ^ cols(M) - 1) ^ 2) for i е a ^ 0.. m for j е 0.. n for d е 1.. D
(y ^ max(f ^ 0, i - d) Y ^ min(m, i + d)) O ^ submatrix(M, y, Y, max(0,j - d), min(n,j + d)) (mi ^ min(O) me^ median(O) ma ^ max(O)) (f ^ 1 break) if mi < me < ma [ Ai j ^ if [(mi < Mi j < ma)-f, Mi j , me] a ^ max(a, d) ] (A 2a + 1)
настраивает во вложенном цикле i = 0, m, j = 0, n для каждого пиксела изображения My его адаптивную медианную яркость Aj, руководствуясь следующей логикой:
• во внутреннем цикле d = 1, D, где D = Lmin(m,n}/2j, производится поиск ширины поля d<D, при которой размер окна (2d+1)x(2d+1) не превышает минимального размера изображения 2D+1. Для этого из матрицы M выбирается прямоугольный блок O покрытых окном пикселов и вычисляются их экстремальные mi=min{O}, ma=max{O} и медианное me=median{O} значения яркостей;
• при выполнении неравенств mi<me<ma бинарные шумы, попавшие в окно, сгруппируются на краях упорядоченной последовательности sort(O) - «перец» на левом, а «соль» на правом. Если при этом еще и mi<Mj<ma, то центр окна также не будет засорен бинарным шумом, что позволяет сохранить его яркость My в (i, ^-м элементе возвращаемой матрицы Aij=My и минимизировать искажения, вносимые обработкой изображения. В противном случае при mi=Mij либо Mij=ma сохраняется менее зашумленное медианное значение Aj=me, соответствующее обычному медианному фильтру;
• если же mi=me либо me=ma, то экстремальные яркости заполнили не менее половины окна O. Это бывает при большой плотности шумов и недостаточном значении d. После инкре-ментирования d++ выполнение предыдущего шага повторяется до тех пор, пока d < D;
• по окончании цикла подбора величины ё элементу матрицы Лц присваивается значение медианы те как наиболее компромиссного результата фильтрации.
Функция теёА1Л(Ы), работа которой иллюстрируется на рис. 4, возвращает список, включающий матрицу адаптивно-медианных яркостей Л и максимально достигнутый размер а=2ётах+1 квадратного аха-окна, где ётах - максимальный размер поля, полученный во всех циклах подбора значений ё для всех пикселов изображения Мц.
твсНШ(ЬГ) тесШ1\¥(1Ч,<1 Б
Рис. 4. Сравнение функций медианной фильтрации
Fig. 4. Comparison of median filtering functions
В представленной программе проведено сравнение трех медианных алгоритмов фильтрации изображения М=«Царап.Ьшр», содержащего царапины и пятна разных размеров и яркостей. В верхнем ряду встроенный фильтр medfilt(M) смог удалить своим 3х3-окном только тонкие, толщиной в один пиксел, царапины, тогда как адаптивный фильтр medfilA(M) с максимальной шириной окна a=27 удалил их все, сохранив при этом четкость границ. Для сравнения алгоритм medfilW(M,9,9) с вручную подобранными размерами окна 9х9, удалившими следы царапин, сильно размыл контуры фигуры. В двух центральных рядах к изображению M добавлены импульсные шумы addnoise(M,p,s) с параметром s=255, делающим их бинарными. При вероятностиp=0.5 функция medfilA(H) с максимальным размером окна ¿=21 удалила и их, и царапины полностью и без размытия границ сплошных объектов, чего не смог сделать алгоритм medfilW(H,b,b) с фиксированным bxb-окном. И даже при 90 % плотности импульсного шума, полностью скрывшего исходный рисунок в хаосе изображения F, алгоритм medfilA(F) с максимальным размером окна с=31 смог более-менее «вытянуть» знакомый образ по сравнению с сильно смазанным изображением medfilW(F,с,с), полученным при постоянном размере окна схс. Наконец, в нижнем ряду показаны весьма посредственные результаты фильтрации, в том числе, и адаптивной, гауссовского шума в изображении N.
Выводы
Сравнение встроенных в MathCAD функций фильтрации шумов в растровых изображениях и вновь разработанных алгоритмов и функций медианной фильтрации однозначно показало более высокое качество очистки последними разнообразных помех - импульсных, линейных и сплошных.
Библиографический список
1. Гонсалес, Р. Цифровая обработка изображений [Текст] / Р. Гонсалес, Р. Вудс. - М.: Техносфера, 2012. - 1104 с.
2. Яне, Б. Цифровая обработка изображений [Текст] / Б. Яне. - М.: Техносфера, 2007. - 584 с.
3. Кирьянов, Д.В. Mathcad 15 / Mathcad Prime 1/0 [Текст] / Д.В. Кирьянов. - СПб.: БХВ-Петербург, 2012. - 432 с.
Дата поступления в редакцию: 14.03.2022