Zpracování obrazu pomocí PDE

Note	Difúzní filtrování (lineární difuze, nelineární izotropní a nelineární anizotropní difuze). Level set metody (pohyb ve směru normály, pohyb řízený křivostí a pohyb ve vnějším vektorovém poli). Fast marching algoritmus. PA166

Tip	Doporučuju kouknout na shrnutí v zápiscích z předmětu PA166 od xrosecky

Gradient $\nabla$: Vektorové pole ve směru největšího nárůstu.

Standardně ho spočítáme jako derivaci obrazu podle x a y. V praxi ale používáme aproximaci derivace podle Taylorova rozvoje.
Aproximace derivace: Taylorův polynom vypadá takto: $f (x + h) = f (x) + h f^{'} (x) + \frac{h ^{2}}{2 !} f^{''} (x) + \dots + \frac{h ^{n}}{n !} f^{(n)} (x) + O (h^{n + 1})$ .

Z něho můžeme odvodit rovnici pro první derivaci (v našem případě ji nazýváme dopředná diference):

$f (x + h) \approx f (x) + h f^{'} (x) f^{'} (x) \approx \frac{f ( x + h ) - f ( x )}{h}$

Tu můžeme dále zpřesnit, pokud si vypíšeme taylorův rozvoj až do druhé derivace včetně (tím získáme centrální diferenci):

$h_{1} = 1, h_{2} = - 1 f (x + h_{1}) - f (x + h_{2}) f (x + 1) - f (x - 1) \approx f (x) + h f^{'} (x) + \frac{h ^{2}}{2 !} f^{''} (x) - f (x) + h f^{'} (x) - \frac{h ^{2}}{2 !} f^{''} (x) f (x + 1) - f (x - 1) \approx 2 h f^{'} (x) f^{'} (x) \approx \frac{f ( x + 1 ) - f ( x - 1 )}{2 h}$

Podobným stylem získáme i druhou derivaci
Divergence: Operace, která nám říká, jak moc vektorové pole míří ven z daného bodu. Pokud je $div j > 0$ , pak se v daném bodě hodnota časem snižuje, pokud je $div j < 0$ , pak se hodnota zvyšuje.

$div j = \nabla^{T} j = (\partial_{x}, \partial_{y}) \cdot (j_{1}, j_{2}) = \partial_{x} j_{1} + \partial_{y} j_{2}$ .
Laplacian $Δ$: Jak se mění teplota v daném místě s časem, proto to odpovídá derivaci teploty podle času.

Laplacian je vlastně divergence gradientu: $Δ u = div (\nabla f) = (\partial_{x}, \partial_{y}) \cdot (\partial_{x} u, \partial_{y} u) = \partial_{xx} u + \partial_{yy} u = u_{xx} + u_{yy}$ .

Difuze

Imituje šíření tepla v 1D tyči / 2D ploše / … Řídí se také rovnicí pro šíření tepla: $u_{t} = Δ u = u_{xx} + u_{yy}$ (kde $Δ$ je Laplaceův operátor).

+ Teplo se v takové ploše může šířit různými způsoby v závislosti na vlastnostech materiálu. Je možné že:

Celý materiál ve všech směrech vede teplo stejně = lineární difuze (LD).
Materiál vede v různých místech teplo různě, ale vždy ve všech směrech stejně = nelineární izotropní difuze (NID).
Materiál vede teplo různě v různých směrech = nelineární anizotropní difuze (NAD).

Lineární difuze

Pokud aplikujeme lineární difuzi na obrázek, dojde k jeho rozmazání (Gaussovský filtr).

Pokud vezmeme prostor všech různých Gaussovských rozmazání v různých časem, dostáváme "Gaussovský prostor měřítek" (Gaussian scale space). Pro ten platí

Zachování průměrné hodnoty šedé
Princip maxima a minima - se zvyšujícím se t se maximum jedině snižuje a minimum jedině zvyšuje
Řešení je nezávislé na
Posuvu hodnoty šedé - “posun po ose Y”
Translaci - “posun po ose X”
Škálování
Vlastnost porovnání - pokud $u \leq v \rightarrow (T_tu) \leq (T_tv)$
Pro dimenze $\geq 2$ mohou vznikat nové extrémy

Chceme vypočítat $u_{ij}^{k + 1}$ pomocí aproximace Taylorova rozvoje:

u_{t} = \frac{u _{ij}^{k + 1} - u _{ij}^{k}}{Δ t} + O (Δ t) u_{xx} = u_{i + 1, j}^{k} - 2 u_{ij}^{k} + u_{i - 1, j}^{k} + O (Δ x^{2}) u_{yy} = u_{i, j + 1}^{k} - 2 u_{ij}^{k} + u_{i, j - 1}^{k} + O (Δ y^{2})

A z této rovnice už vyjádříme $u_{ij}^{k + 1}$ . ( $Δ t$ je časový krok, $Δ x$ a $Δ y$ jsou prostorové kroky)

Nelineární izotropní difuze

Zobecníme funkci pro difuzi kombinací Fick’s law a Mass preservation: $\partial_{t} u = div (g \cdot \nabla u)$ , kde $\div$ je divergence.

Lineární
- $g = 1$
- Rozmazání ve všech bodech stejně ve všech směrech
Nelineární izotropní
- $g =$ skalární funkce
- Rozmazání stejně ve všech směrech, ale v každém bodě jinak
Nelineární anizotropní
- $g =$ maticová funkce
- Rozmazání v každém bodě a směru jinak

Pro NID chceme typicky zabránit rozmazání na výrazných hranách. Chceme tedy, aby fce $g$ byla na hranách (tam, kde je velký gradient) co nejmenší a jinde co největší. Existuje několik různých vzorců:

Perona-Malik difuzivita: $g (\nabla u) = \frac{1}{1 + ∣\nabla u ∣ ^{2} / λ ^{2}}$ .
Charbonnier difuzivita: $g (\nabla u) = \frac{1}{1 + ∣\nabla u ∣ ^{2} / λ ^{2}}$ .
Exponenciální difuzivita: $g (\nabla u) = e^{\frac{- ∣\nabla u ∣ ^{2}}{2 λ ^{2}}}$ .

Nelineární anizotropní difuze

NID nechával "chlupaté hrany", protože kolem hran nerozmazával vůbec. NAD dokáže kolem hran rozmazat jen v tom správném směru, neboť $g$ je v tomto případě matice.

Symetrická matice

Reprezentuje otočení, roztažení a otočení zpátky:

A = (u_{1} u_{2} v_{1} v_{2}) (λ_{1} 0 0 λ_{2}) (u_{1} v_{1} u_{2} v_{2})

Vlastní vektory kolmé na sebe
Vlastní čísla jsou reálná

Díky symetrickým maticím můžeme sestavit sami maticky s předem danými vlastními vektory a čísly.

Edge-enhancing difuzivita (difuzivita zvýrazňující hrany): My chceme matici, co má následující vlastnosti:

$v_{1} ∥ \nabla u_{σ}, λ_{1} = g (∣\nabla u_{σ} ∣^{2}) v_{2} ⊥ \nabla u_{σ}, λ_{2} = 1$

Kde funkce $g$ je funkce nepřímé úměry. Tím zajistíme, že rozmazání podél hran bude maximální a v ostatních směrech bude minimální.
Coherence-enhancing difuzivita (difuzivita zvýrazňující koherenci): Potřebujeme vypočítat difuzní tenzor a podle něj potom aplikujeme rozmazání. Tenzor je dán jako:

$D = \nabla u_{σ} \nabla u_{σ}^{T}$

My opět vytvoříme matici s vlastními vektory a čísly. Tentokrát vektory vezmeme z tensoru: Vlastní vektor tensoru s větším vl. číslem vede přes strukturu → Dáme mu malé číslo blízko nuly; vlastní vektor s menším vl. číslem vede kolmo → Dáme mu číslo podle rozdílu koherence (pokud je koherentní, nechceme rozmazávat tolik, může to být roh).

Level set metody

Level set metody využívají implicitní reprezentace křivek.

Implicitní reprezentace: Představme si funkci, která má uvnitř křivky záporné hodnoty, na křivce nulu a venku kladné hodnoty. Potom můžeme křivku reprezentovat, jako (nulovou) vrstevnici funkce.

Při takovéto reprezentaci nemáme přístup přímo k hranici, ale můžeme jí získat například pomocí marching squares (cubes) algoritmu. Zároveň musíme mít uloženou hodnotu funkce pro celý obraz, což může být nevýhoda. Výhodou je, že můžeme snadno měnit topologii křivky (přidávat díry, spojovat křivky, …).

Toto všechno funguje i ve 3D, kde se ale bavíme o povrchu.

Vývoj křivky můžeme definovat, jako $\frac{\partial C}{\partial t} = β n$ (kde $n$ je normála a $β$ řídí rychlost evoluce). Pokud ho chceme definovat v rámci obalující funkce $u$ , můžeme ho zapsat, jako $\partial_{t} u = β ∣\nabla u ∣$ . $β$ v této rovnici ovlivňuje směr a rychlost pohybu křivky, $∣\nabla u ∣$ je velikost gradientu.

Level set metody se snaží řešit tuto rovnici.

Typy pohybu

Existují 3 základní typy pohybu křivky:

Pohyb ve směru normály - $β = a$ = dilatace / eroze
Pohyb řízený křivostí - $β = - ϵ κ$ = vyhlazování křivky
Pohyb ve vnějším vektorovém poli - $β = V (x, y, t) \cdot n$ , kde $V$ je vnější vektorové pole

Pohyb ve vnějším vektorovém poli

$β = V (x, y, t) \cdot n$ popisuje pohyb ve vektorovém poli definovaném parametrem $V$ . Aproximujeme opět pomocí Taylorova rozvoje.

Pohyb řízený křivostí

Chceme, aby se křivka vyhlazovala, protože křivost je vlastně druhá derivace. $β = - ϵ κ$ popisuje pohyb křivky ve směru opačném k její křivosti. $κ$ je křivost a $ϵ$ je parametr, který řídí rychlost vyhlazování.

Pohyb ve směru normály

Chceme, aby se křivka rozšiřovala, nebo smršťovala. $β = a$ popisuje pohyb křivky ve směru normály. $a$ je parametr, který řídí rychlost pohybu.

Všechny typy pohybu můžeme zapsat do jedné rovnice: $u_{t} = - V \cdot \nabla u - a ∣\nabla u ∣ + ϵ κ ∣\nabla u ∣$ a aproximovat pomocí Taylorova rozvoje.

Fast marching algoritmus

FMA je specifický případ Level set metody s pouze jedním typem pohybu ve směru normály $\partial_{t} u = a ∣\nabla u ∣$ .

Pro různé případy můžeme volit různé hodnoty $a$ .

Euklidovská vzdálenost: Pro výpočet vzdálenosti od křivky ke všem bodům obrazu zvolíme $a = 1$ .
Geodesická vzdálenost: Chceme-li změřit vzdálenost od daného bodu uvnitř objektu, zvolíme $a = 1$ uvnitř objektu a $a \to 0$ venku.
Segmentace: Můžeme jednoduše segmentovat obraz, pokud zvolíme $a = g (∣\nabla f ∣)$ (kde $f$ je obraz a $g$ je nějaká funkce, například Perona-Malik). Potom $g$ jde k 0 při vysokém gradientu a k 1 při nízkém gradientu.

Samotný výpočet algoritmu je trošku složitější, uděláme si 3 množiny

Trial (co chceme testovat, init. kontura)
Far (neobjevené, init. všechno krom kontury)
Known (vyřešené, init. prázdná množina)

Vždy najdeme v Trial nejmenší hodnotu, přehodíme ji do Known a pro její sousedy vypočítáme nový arrival time minimalizací $T_{n e w} = τ (v_{1}, v_{2}) + v_{1} T_{1} + v_{2} T_{2}$ , kde $τ$ je časový krok, $T_{1}, T_{2}$ jsou hodnoty v Known a $v_{1}, v_{2}$ jsou souřadnice bodu X, které chceme minimalizovat ( $v_{1} + v_{2} = 1, v_{1}, v_{2} \geq 0$ ).

Konkrétní výpočet je ve slidech, nebo v notionu, ale myslím, že je celkem zbytečný.

Složitost je v řádu $O (n lo g n)$ , kde $n$ je počet pixelů v obrazu.