ЛЕКЦІЯ 9. ВЕЙВЛЕТ-ПЕРЕТВОРЕННЯ ТА ЙОГО ЗАСТОСУВАННЯ
Вейвлет-перетворення (wavelet transformation) – це сучасний і перспективний метод обробки даних. Англійське слово wavelet (від французького "ondelette") дослівно перекладається як "коротка (маленька) хвиля". Апарат вейвлет-аналізу одержав свій розвиток на початку 1980-x років у роботах Морле, Гроссмана й деяких інших авторів . Результаты, отримані у різних областях за допомогою вейвлет-аналізу, підсилили інтерес до цього напрямку та сприяють його безупинному розвитку.
Методи вейвлет-аналізу можливо застосувати до даних різної природи. Це можуть бути, наприклад, одномірні функції або двовимірні зображення. Грубу класифікацію вейвлет-алгоритмів можна зробити, виділивши безперервне (CWT – Contіnuous Wavelet Transform) і дискретне (DWT – вDіscrete Wavelet Transform) вейвлет-перетворення. Одержати набір вейвлет-коефіцієнтів у випадку дискретного перетворення швидше, і воно дає досить точне представлення сигналу при меншому обсязі одержуваних у результаті даних. Безперервне перетворення вимагає більших обчислювальних витрат, але, разом із цим, дозволяє детальніше роздивитися структуру сигналу.
На відміну від звичайних спектральних перетворень, вейвлет-аналіз дозволяє з однаковою точністю апроксимувати як гладкі функції, так і функції з різкими випадами, що дає можливість визначати незначні об’єкти. Застосування вейвлетів розглянуто в багатьох роботах, наприклад. Оскільки вейвлет-перетворення представляє собою згортку сигналу зображення з вейвлет-функцією, що зміщується відносно координат, даний метод близький за змістом до методів на основі фільтрації.
Вибір того або іншого методу залежить від поставленого завдання й типу наявних даних, які необхідно обробити, від можливостей обчислювальної техніки і від того, у якому вигляді необхідно представити результат.
9.1 ПОНЯТТЯ ПРО ВЕЙВЛЕТИ
Термін вейвлет-перетворення об’єднує два види перетворень – пряме і обернене, які, відповідно, переводять досліджувану функцію f(x) в набір вейвлет-коефіцієнтів Wψ(a,b)f і назад. Розділяють безперервне та дискретне перетворення
Пряме вейвлет-перетворення здійснюється відповідно до правила
(1)
де a і b – параметри, що визначають відповідно масштаб і зсув функції ψ, яка називається аналізуючим вейвлетом, Cψ - нормувальний множник. Інтегрування ведуть по всій числовій осі.
Базисний, або материнський вейвлет (parent wavelet) ψ створює за допомогою розтягнень та зсувів сім’ю ψ(х-b)/a.
Маючи відомий набір коефіцієнтів Wψ(a,b), можна відновити первинний вигляд функції f(x):
(2)
Пряме (1) і обернене (2) перетворення залежать від деякої функції ψ(х)єL2(R), яку називають базисним вейвлетом. Практично єдиним обмеженням на його вибір є умова обмеженості нормувального множника
(3)
де Ψ(ω) – Фур'є- образ вейвлета ψ(х): 
.
Безліч функцій задовольняють дану умову, тому можливо підібрати вид вейвлета, який найбільш підходить для рішення конкретного завдання.
Умова (3) означає, що Фур'є-образ вейвлета дорівнює нулю при нульовій частоті, тобто. Якщо це не так, то знаменник в інтегралі (3) збігається до нуля, у той час як чисельник має відмінне від нуля значення, і коефіцієнт Cψ перестає бути обмеженим.
У свою чергу, цю вимогу можна представити в іншому вигляді. Оскільки Фур'є-образ Ψ(ω) при нульовій частоті має вигляд 
, ми можемо очікувати рівність нулю інтеграла від вейвлета по всій осі:
(4)
9.2 ГОЛОВНІ ОЗНАКИ ВЕЙВЛЕТА
В якості базисних функцій, що утворюють ортогональний базис, можна використовувати широкий набір вейвлетів. Для практичного застосування важливо знати ознаки, якими неодмінно повинна володіти вихідна функція, щоб стати вейвлетом. Приведемо основні з них.
Обмеженість. Квадрат норми функції повинен бути обмеженим:
(5)
Локалізація. ВП на відміну від перетворення Фур'є використовує локалізовану вихідну функцію і у часі, і по частоті. Для цього досить, щоб виконувалися умови:
і 
, при ɛ>0 (6)
Нульове середнє. Графік вихідної функції повинен бути знакозмінним навколо нуля на осі часу і мати нульову площу
(7)
Рівність нулю площі функції ψ(t), тобто нульового моменту, призводить до того, що Фур'є – перетворення Sψ(ω) цієї функції дорівнює нулю при ω=0 і має вигляд смугового фільтра.
Автомодальність. Характерною ознакою ВП є його самоподібність. Всі вейвлети конкретного сімейства ψa,b (t) мають те ж число осциляцій, що й материнський вейвлет ψ(t), оскільки отримані з нього за допомогою масштабних перетворень (a) і зсува (b).
Основні функції, які створюють вейвлети, так звані материнські вейвлети, наведені в табл. 1.
Найбільш поширені базиси створюються на основі похідних функції Гауса G0(t)=exp(-t2/2)). Це обумовлено тим, що функція Гауса має найкращі показники локалізації як у часовій, так і в частотній областях.
На рис. 9.1 показані вейвлети перших чотирьох порядків і модулі їх спектральної щільності. При n =1 одержуємо вейвлет першого порядку, який називається WAVE - вейвлетом з нульовим моментом рівним нулю. При n = 2 одержуємо MHAT – вейвлет, що називається "мексиканський капелюх" (mexіcan hat – схожий на сомбреро). У нього нульовий і перший моменти дорівнюють нулю.
Таблиця 1. Материнські вейвлети
Найбільш простий приклад дискретного вейвлета – це HAAR-вейвлет. Його недоліком є несиметричність форми та не гладкість – різкі границі в t-області. Серед комплексних вейвлетів найбільше часто використовується базис, заснований на добре локалізованому у часовій та в частотної областях вейвлеті Морле. Характерний параметр ω(0) дозволяє змінювати вибірковість базису. Дійсна та уявна частини ψ(t) – це амплітудно-модульовані коливання.
Рисунок 9.1 – Гаусові вейвлети перших чотирьох порядків (зліва) та їх спектральні щільності у частотній області (справа)
Вище був представлений невеликий перелік типів вейвлетів, що аналітично описуються в явному вигляді. Однак більшість типів вейвлетів не мають аналітичного опису у вигляді однієї формули, а задаються ітераційними виразами, які легко обчислюються комп’ютерами. Прикладом таких вейвлетів є функції Добеши (Daubechіes), одна з яких (db4) використовується в якості вбудованої в Mathcad.
У наш час вибір вейвлетів досить великий. Тільки в пакеті Wavelet Toolbox (MATLAB ) представлено півтора десятка материнських вейвлетів; при цьому для ряду з них дано ще безліч варіантів. Для одержання довідки про будь-який тип вейвлета при роботі в командному режимі MATLAB досить виконати команду waveіnfo ('type'), указавши тип вейвлета.
Для перегляду вейвлетів достатньо виконати команду wavemenu і у вікні, що з'явилося, зі списком вейвлет перетворень нажати кнопку Wavelet Dіsplay. Натискання цієї кнопки виводить вікно перегляду вейвлетів, у якому можна переглядати: обраний вейвлет (з ім'ям 'Name') та інформацію про нього.
Вибір конкретного материнського вейвлета цілком залежить від характеру поставленого завдання та від конкретного аналізованого сигналу. А з точки зору обробки зображень вибір конкретного типу материнського вейвлета залежить від не лише від суті завдання обробки, а й від конкретного вигляду зображення.
9.4 ДВОВИМІРНІ ВЕЙВЛЕТИ
При обробці зображень доводиться мати справу з двовимірними масивами S(x, y) . Вони задаються в просторі як функції двох змінних x і y . У цьому випадку двовимірна вейвлет-функція має вигляд:
(8)
де a1 і a2, b1 і b2 – значення a і b по кожному виміру.
Для вейвлет-перетворень (ВП) дискретних зображень батьківський та материнський вейвлети будують наступним чином:
(9)
(10)
де індекси H і L означають реалізацію фільтрів високих частот та низьких частот.
Тоді двовимірні вейвлети запишуться в такому вигляді:
(11)
Таким чином, на двовимірній площині відбувається аналіз по горизонталі, вертикалі й діагоналі з однаковим розширенням відповідно до трьох наведених вище вейвлетів. Формули дискретного ВП (ДВП) двовимірних сигналів і зображень, створені з врахуванням наведених вище співвідношень (11), використані в пакеті Wavelet Toolbox.
З точки зору обробки зображень вибір конкретного типу двовимірного вейвлета залежить від не лише від суті завдання обробки, а й від конкретного вигляду зображення. Тому потрібно обчислювати якийсь кількісний показник для різних типів вейвлетів, так щоб він досяг екстремуму, тобто мінімуму чи максимуму залежно від суті задачі.
Для задачі фільтрації таким показником може бути норма зашумленого зображення, яка повинна бути мінімальна, що означає, що енергія залишкового шуму мінімальна. Якщо ж існує початкову незашумлене зображення, тоді доцільно оцінювати норму різниці відфільтрованого та початкового зображення, або її квадратний корінь (див. вираз (12)).
9.5 АЛГОРИТМ ПРЯМОГО ВП
Пряме ВП зображення відбувається таким способом. Припустимо, що маємо зображення розміром NxN (рис. 9.2,а). Спочатку кожний з N рядків зображення ділиться (фільтрується) на низькочастотну (НЧ) і високочастотну (ВЧ) половини. У результаті виходить два зображення розміром NxN/2 (рис. 9.2, б). Далі кожний стовпець ділиться так само, у підсумку виходить чотири зображення розміром N/2xN/2 (рис. 9.2, в): НЧ по горизонталі та вертикалі (НЧНЧ1), ВЧ по горизонталі та вертикалі (ВЧВЧ1), НЧ по горизонталі та ВЧ по вертикалі (НЧВЧ1) і ВЧ по горизонталі та НЧ по вертикалі (ВЧНЧ1). Перше із зазначених вище зображень ділиться аналогічним чином на наступному кроці (рівні) перетворення (рис. 9.2, г) і т.д.
Для прикладу розглянемо як можна здійснити видалення шумів та компресію зображень за допомогою ВП. Для вирішення даної задачі здійснюється граничне обмеження рівня коефіцієнтів, що деталізують. Задавши певний поріг і "відтинаючи" коефіцієнти нижче цього порога, можна значно знизити рівень шуму і стиснути зображення.
Рисунок 9.2 – Загальна схема фільтрації зображень за допомогою ДВП.
Нижче наведений фрагмент програми фільтрації зображення від шумів, що завантажене з файлу noіse. При цьому використані функції wpbmpen (установка глобального порога) і wpdencmp (видалення шумів і стиск зображень):
load noissi2d; nbc = size(map,1); wname = 'db8';
lev = 2; tree = wpdec2(X,lev,wname);
det1 = [wpcoef(tree,2) wpcoef(tree,3) wpcoef(tree,4)];
sigma = median(abs(det1(:)))/0.6745; alpha = 1.1;
thr = wpbmpen(tree,sigma,alpha); keepapp = 1;
xd = wpdencmp(tree,'s','nobest',thr,keepapp);
colormap(pink(nbc));
subplot(221), image(wcodemat(X,nbc));
title('Вихідне зображення')
subplot(222), image(wcodemat(xd,nbc));
title('Відфільтроване зображення')
end
Рисунок 9.3 – Схема 2-етапного прямого ВП зображення
Рисунок 9.4 – Вихідне і відфільтроване зображення
На рисунку 9.4 наведені вихідне і відфільтроване зображення, отримані при виконанні програми. Ряд прикладів по аналізу та реконструкції зображень, їх компресії та фільтрації від шуму наведено у розділі GUІ Wavelet Toolbox, що активізується кнопкою Wavelet-2D.
9.6 ПРИКЛАД ЗАСТОСУВАННЯ ДВП ДЛЯ СТИСНЕННЯ
На рис. 9.5 наведено приклад застосування Wavelet-2D: у верхньому лівому куті наведено реальне зображення; у нижньому правому куті дано його вейвлет-розклад (dwt) на прямокутні сегменти; у лівому нижньому куті реконструкція зображення (іdwt); у правому верхньому куті ілюструється ділянка декомпозиції зображення.
До стиску зображень проявляється значний інтерес в усьому світі. Це обумовлено стрімким розвитком цифрової техніки обробки зображень, кольорових принтерів, графічних моніторів, цифрових фото- і відеокамер тощо. Зображення, представлене в цифровому вигляді, має досить великий обсяг у бітах. Наприклад, кольорове зображення розміром 512х512 вимагає для свого зберігання 768 кбайт, а якщо передавати відео послідовність таких зображень зі швидкістю 25 кадрів у секунду, тоді необхідна швидкість складе 188 Мбіт/с.
Рисунок 9.5 – Розпізнавання об’єктів методами 2-вимірного ДВП
Принципова відмінність процедури компресії за допомогою ДВП від широко розповсюдженого стиску за стандартом JPEG полягає в тому, що він працює з усім зображенням, у той час як в JPEG зображення розбивається на блоки, які стискуються незалежно. У випадку ВП можна підібрати таку базисну вейвлет-функцію, що адаптована до найбільш інформативних особливостей зображення. Під адаптивністю мають на увазі, що елементи або ділянки зображення з досить плавною зміною яскравості представляються невеликим числом вейвлет- коефіцієнтів.
Рисунок 9.6 – Вигляд початкового (зліва) та стисненого (справа) зображень.
9.7 ОПТИМАЛЬНА ФІЛЬТРАЦІЯ НА ОСНОВІ ДВП
У п. 5 наведено приклад застосування ДВП для фільтрування шумів (див. рис. 9.3 та програмний код). Проте, це дуже простий приклад, бо має такі недоліки:
- зображення є змодельованим та являє собою 2D синусоїду, тобто Фур’є гармоніку,
- вся процедура фільтрації схована у функції wpdencmp,
- всі параметри процедури thr (поріг) та keepapp статичні і ніяк не підбираються, тобто неочевидно, що вони найкращі
Іншими словами, для кожного конкретного зображення потрібно оптимізувати всі кількісні параметри. Поріг обчислюється за допомогою функції wpbmpen, але теж не підбирається і в свою чергу, залежить від параметрів sigma та alpha. Перший з них – це по суті оцінка СКВ шуму. Також, всі ці параметри неявно залежать від типу вейвлету, за який взято db8.
Раніше у Лекції 7 було дано визначення оптимального фільтра як фільтра, при якому досягається мінімум середньоквадратичної похибки (СКП). Для задачі фільтрації СКП, це по суті СКВ між відфільтрованим та початковим зображенням, яке не містить інформативної компоненти, тобто чисто шумове зображення. Але у реальних умовах у нас немає такого зобра-ження, тому СКП оцінюють між відфільтрованим зображенням та початковим, на яке накладено шум, тобто зашумленим, згідно (12).
(12)
де RMSE – root mean square estimation (тобто СКП), Sj,n, S’j,n - значення пік-селів, відповідно зашумленого та початкового зображень, M та N – кількість рядків та стовпчиків пікселів, тобто розмір зображення.
Для синтезу оптимального фільтру на основі ДВП потрібно оптимізувати найбільш критичні параметри - тип вейвлету, величину порогу, а також функцію, яка описує правило, згідно якого вибраковуються спектральні компоненти ДВП (аналогічна АЧХ для звичайного фільтра), яку будемо називати «порогова функція». У наведеному прикладі нижче буде вибрано 4 типи вейвлетів (Хаара, дуально-комплексний та вейвлет в орієнтованому базисі, ОБ) – та дві порогові функції (м’яку та жорстку). Вигляд порогових функцій дається виразом (13), де С – поріг, (13а), зліва – жорстка функція, (13б), справа – жорстка функція.
(13а,б)
1-вимірне перетворення Хаара, тобто елементи НЧ та ВЧ компонент ДВП (див Рис. 3б) даються виразами
(14)
Зворотне перетворення, або реконструкція зображення S виконується згідно виразів (15)
(15)
Дуальне ДВП реалізується паралельним обчисленням двох 2D ДВП, потім для кожної пари субсмуг береться сума та різниця. 2D комплексне складається із 2-х паралельно застосованих дійсних ДВП, одне з яких аналогічно 1-мірному випадку та інтерпретується як уявна частина комплекс-ного вейвлета, а інше – як дійсна. У дуально-комплексному ДВП придушуються не самі коефіцієнти, а їх модулі.
Особливістю ДВП в ОБ є використання 3-х різних типів фільтрів, тому розкладання зображення виконується за 3-ма видами коефіцієнтів. Розрахунок є ітераційним, за якої аналіз функції проводять на різних рівнях розкладу від мінімального інтервалу, який включає три відліки, до максимального, рівного інтервалу визначення.
Результати порівняння різних типів фільтрів на основі 2-х порогових функцій та 4-х типів вейвлетів, тестового зображення, поданого на Рис. 9.7, зашумленого гаусовим шумом, наведено на Рис. 9.8 та в Таблицях 3 та 4.
Таблиця 2. ОБ вейвлет перетворення
|
Пряме ОБ вейвлет перетворення |
Зворотне ОБ вейвлет перетворення |
|
|
|
Рисунок 9.7 – Вигляд зашумленого тестового зображення.
Рисунок 9.8 – Залежності СКП фільтрації від величини порогу для м’якої (зліва) та жорсткої (справа) порогової функції для 4-х типів вейвлетів.
Таблиця 3. Відфільтроване зображення
|
Задання порогу |
ДВП Хаара |
ОБ |
Дуальне ДВП |
Дуально-комплексне ДВП |
|
м’яке |
|
|
|
|
|
жорстке |
|
|
|
|
Таблиця 4. Порівняння методів фільтрації
З Табл. 4 видно, що досліджені методи вейвлет-фільтрації мають майже однакову ефективність, тобто значення мінімальної СКП, але для дуального та дуально-комплексного методів вони досягаються за менших значеннях порогу. Отже, вейвлет-фільтрація слабко залежить від виду вейвлету, бо мінімальна СКП коливається в межах 7-9%. При цьому найгіршим є застосування вейвлетів Хаара (СКП=8-9%), а найкращим випадком є дуально-комплексне ДВП (СКП=7,5%).
На противагу цьому видно, що величина порогу дуже чутлива як від типу вейвлету, так і від виду порогової функції. А саме, – для жорсткого способу задання порогу, оптимальне значення порогу приблизно в 2 рази більше, ніж для м’якого способу. Також застосування більш складних вейвлетів також приблизно в 2 рази зменшує величину порогу. Тут знову найпростіший вейвлет Хаара дає найбільший поріг (С=35/70 для м’якого /жорсткого способу), а найменшу величину забезпечує дуально-комплексне ДВП (С=20/35 для м’якого/жорсткого способу). Отже, для синтезу оптимального вейвлет-фільтру потрібно досліджувати його залежність від порогу і знайти оптимальне значення, при якому досягається мінімум СКП.
Деталі цього прикладу описано в статті: Терещенко В.О., Ямненко Ю.С., Мельниченко О.Л., Панченко М.В. Вейвлет-перетворення для фільтрації зображень із відеокамер спостереження. Вчені записки ТНУ імені В.І. Вернадського. Серія: Технічні науки. 2018, Том 29(68), Ч.2, №3, с. 14-18.
