2. cvičení: Vylepšení/restaurace obrazu – zaostření, dekonvoluce, FFT

Cílem druhého cvičení je praktické vyzkoušení vylepšení obrazu. Jednoduchou, ale velmi často používanou metodou, je zaostření zvané unsharp masking, které se vyskytuje ve většině obrazových editorů. Základem je aplikace jednoduchého 2D filtru, respektive 2D konvoluce, kombinované s jednoduchou aritmetickou operací. Tato technika se používala již v dobách analogového zpracování filmů. Obrazové vysvětlení procesu v přednaškových slidech, viz část věnovaná unsharp masking a Laplaciánu. V současné době se s výhodou používá například jako předzpracování před tiskem na inkoustové tiskárně, která z principu každý obraz trochu rozmaže.

Komplikovanější problém pak představuje rekonstrukce obrazu výrazněji poškozeného, kde jednoduché ostřící metody selhávájí. Uvidíme, jak se lze s výhodou přesunout konvolučním teorémem do frekvenční oblasti a provést některé operace v ní.

Konvoluce

S konvolucí a konvolučními maskami už jste se seznámili na přednášce. Představte si kameru směřující kolmo na list papíru, pro zjednodušení uvažujte, že všechny body na listu jsou od kamery stejně vzdálené. Víme, co je na listu papíru, a chceme konvolucí simulovat obraz

1. v ideálně zaostřené statické kameře;	2. v kameře zaostřené mimo rovinu listu;	3. v kameře, která se posouvala konstantním směrem a rychlostí.

Zkuste se nejdřív sami zamyslet, jak byste toto rozmazání simulovali, pak se podívejte jak to bylo provedeno skriptem simulateblur.m.

Frekvenční oblast, FFT

Nyní si představte, že stojíte před úkolem restaurace obrazu znehodnoceného rozmazáním rovnoměrným pohybem. Toto rozmazání lze modelovat konvolucí s obdélníkovým jádrem o šířce jednoho pixelu a délce odpovídající délce pohybu kamery. Se znalostí tohoto jádra můžeme provést dekonvoluci. Jedinou neznámou je délka pohybu. V prostorové oblasti bychom délku mohli změřit např. z obrazu vzniklého rozmazáním jediného bílého bodu. Ve frekvenční oblasti ale nepotřebujeme takový speciální případ, aby byla délka pohybu jasně měřitelná:

im	log(im_Pfft + 1)	median(log(im_Pfft + 1))

Jak provést FFT v MATLABu ukazuje demo (kód). Pro vysvětlení tmavých pruhů ve spektru potřebujeme znát konvoluční teorém a odvodit obraz obdélníkového jádra ve frekvenční oblasti. Toto odvození zároveň dává návod, jak spočítat velikost rozmazání. Obrázkový úvod do diskrétní FT a konvoluce.

Dekonvoluce

Rozmazaný obraz g vznikl konvolucí originálního obrazu f s "rozmazávacím" jádrem h: g=f*h. Konvoluční teorém o obrazech ve frekvenční oblasti říká G=F.*H a hledaný obraz f tedy lze vyjádřit jako f=ifft2(G./H). Při implementaci je potřeba ošetřit dělení nulou přičtením např. konstanty eps.

Wienerova filtrace

Pokud vám výsledek přijde neuspokojivý, a jste si jisti, že implementace je správně, nezoufejte. Pro možné vysvětlení doporučujeme se podívat do skript na část 11.5 Obnovení obrazu při známé degradaci. Je nutné si uvědomit, že situace není tak ideální, jak si modelujeme. Vysledný (pokažený) obraz g je ve skutečnosti zatížen také šumem, tedy g = f*h + n. Při dekonvoluci tedy dostaneme F = G./H + N./H, kde N je Fourierův obraz šumu, který obvykle neznáme. Aditivní člen N./H se projeví především tam, kde má inverzní filtr malou amplitudu. Velikost modulu H může klesat rychleji než N, a šum tak převáží nad užitečnou informací. Je nutné si uvědomit, že šum na vstupu vzniká už jen zaokrouhlením při kvantizaci (jas ve vstupním obraze nabývá jen 256 hodnot).

Pro vylepšení výsledku je možné celkem jednoduše rozšířit základní implmentaci na Wienerův filtr. F = G.*H_w, kde H_w = (1./H).*(abs(H).^2 ./ (abs(H).^2 + K)). Konstanta K modeluje poměr šum/signál.

Pokud nemáte po ruce skripta, je možné podívat se na zhuštěný výklad dostupný on-line.

Samostatná práce

1. Napište funkci Iout = imsharpen(I,c). Vstupní parameter c určuje míru zaostření. Čím vyšší hodnota, tím větší efekt, pro $c=0$ pak musí platit I=Iout. Rozumné chování funkce se očekává pro hodnoty . [1 bod]

2. Pro experimenty s restaurováním obrazu/dekonvolucí můžete použít předpřipravený rozmazaný obraz. Budeme uvažovat rozmazání pohybem pouze ve směru řádků nebo sloupců, v celočíselných hodnotách.

3. Implementujte funkci Iout = deconvolve(I,h), kde h je konvoluční jádro, které způsobilo rozmazání. Řádkový vektor znamená rozmazání v horizontálním směru, sloupcový vektor pak ve směru vertikálním. Připomínáme požadavek shodného typu a velikosti Iout a I. Splnění tohoto požadavku bude možná tady vyžadovat trochu péče. [1 bod]

4. Implementujte funkci Iout = removemotionblur(I,length,direction), kde vstupní parametery length a direction jsou nepovinné. V případě, že je funkce volána pouze s jedním parametrem, musí funkce odhadnout směr rozmazání a jeho délku sama. Funkce bude testována postupně ve variantách:

   • pouze I - čili funkce musí spočíst jak délku tak i směr rozmazání
   • pouze I,length - neboli funkce musí ještě rozhodnout, zda rozmazání je vertikálním či horizontálním směrem.
   • se všemi třemi parametry. Funkce tedy pouze potřebuje zkontruuovat správné jádro.

   Počet vstupních parametrů lze zjistit pomocí nargin. Očekávává se, že funkce použije výše požadovanout funkci deconvolve. Parametr length je nezáporný, celočíslený, parametr direction je 'char', může nabývat hodnoty 'h' pro horizontální směr nebo 'v' vertikální směr rozmazání. [2 body]