Jelfeldolgozási megközelítés szkenner
A 3D képfeldolgozás és vizualizáció matematikai alapjai
újramintavételezés újrakvantálás
Csébfalvi Balázs E-mail:
[email protected] Irányítástechnika és Informatika Tanszék Budapesti Műszaki és Gazdaságtudományi Egyetem
Tomográfiás vs folytonos rekonstrukció
monitor
diszkrét 2D kép D/A konverzió
2 / 96
Mintavételezés és approximáció elmélete
Tomográfiás rekonstrukció
• A 3D szkennerek 2D vetületeket mérnek • Ezekből kell egy diszkrét 3D adathalmazt előállítani • A rekonstrukciót egy 3D rács rácspontjaiban számoljuk ki
diszkrét 3D adat A/D konverzió mintavételezés kvantálás
Mintavételezés • Feltételezzük, hogy az eredeti jel sávkorlátozott • Melyek a tökéletes rekonstrukció feltételei? • Hogyan lehet az ideális esetet közelíteni?
Folytonos rekonstrukció
Approximáció elmélete • Feltételezzük, hogy az eredeti jel több rendben folytonos • Az eredeti jelet minél magasabb rendben szeretnénk közelíteni • A négyzetes hibát minimalizáljuk
• A diszkrét reprezentáció nem elég a 3D megjelenítéshez • A vetítő sugarak mentén tetszőleges pontban szeretnénk mintát venni • A diszkrét jelből egy folytonos jelet kell rekonstruálni 3 / 96
Matematikai apparátus
Konvolúciós szűrés
Diszkrét Fourier transzformáció
Taylor sor szerinti kifejtés
Approximációs rend
4 / 96
Folytonos rekonstrukció
Az eredeti f(t)folytonos jel csak diszkrét ti=iT pontokban ismert: fi = f(ti)
A folytonos rekonstrukciót egy konvolúciós szűréssel végezhetjük el: ~
f (t ) ≈ f (t ) =
∞
∑
k = −∞
Tér és frekvenciatartománybeli minőségi kritériumok
5 / 96
f kϕ (
t − k) T
A rekonstrukciós szűrő impulzusválasza: φ(t) 6 / 96
1
Szakaszonként lineáris közelítés
Szakaszonként lineáris közelítés
Lineáris sátor szűrő:
Lineáris sátor szűrő:
1− | t | if | t |< 1 ϕ (t ) = otherwise 0
1− | t | if | t |< 1 otherwise 0
ϕ (t ) =
~
~
f (t )
f (t )
t
7 / 96
Szakaszonként lineáris közelítés
8 / 96
t
10 / 96
t
12 / 96
Szakaszonként lineáris közelítés
Lineáris sátor szűrő:
Lineáris sátor szűrő:
1− | t | if | t |< 1 ϕ (t ) = otherwise 0
1− | t | if | t |< 1 otherwise 0
ϕ (t ) =
~
~
f (t )
f (t )
t
9 / 96
Szakaszonként lineáris közelítés
t
Szakaszonként lineáris közelítés
Lineáris sátor szűrő:
Lineáris sátor szűrő:
1− | t | if | t |< 1 ϕ (t ) = otherwise 0
1− | t | if | t |< 1 otherwise 0
ϕ (t ) =
~
~
f (t )
f (t )
t
11 / 96
2
Lineáris rekonstrukció
Fourier analízis
Csak C0 folytonos Könnyen és hatékonyan implementálható Csak a két legközelebbi mintát kell súlyozva összegezni:
Milyen a szűrő átviteli karakterisztikája
Aliasing hatás: A rekonstruált jelben megjelenhet olyan frekvenciakomponens, amely az eredeti jelben nincs benne
Általában az aliasing hatást csak szélesebb hatósugarú (költségesebb) szűrőkkel lehet csökkenteni
Simító hatás (oversmoothing): A szűrő elnyom olyan frekvenciakomponenst, amely az eredeti jelben benne van
~
f (t ) = (t − ti ) f i +1 + (ti +1 − t ) f i =
f i + ( f i +1 − f i )(t − ti ) f i +1
~
f (t )
fi
ti
t
ti +1
13 / 96
Fourier-transzformáció
Fourier-transzformáció (körfrekvenciában)
Fourier-transzformáció (ν: közönséges frekvencia)
F (ν ) =
∞
∫
14 / 96
Fourier-transzformáció (ω = 2πν: körfrekvencia) o
f (t )e −i 2πνt dt
F (ω ) =
Inverz Fourier-transzformáció (szimmetrikus)
f (t ) =
∞
∫ F (ν )e ν
∫ f (t )e
−iωt
dt
t = −∞
t = −∞
∞
i 2πνt
dν
Inverz Fourier-transzformáció (nem szimmetrikus)
f (t ) =
= −∞
∞
1 ω iωt dν F dω = e ∫ dω 2π ω = −∞ 2π
∞ o
∫ F (ω )e ω
iωt
15 / 96
Fourier-transzformált párok Tértartomány
Frekvenciatartomány Frekvenciatartomány (közönséges frekvencia) (körfrekvencia)
f (at )
( f ⊗ g )(t ) 1 δ (t ) box(t )
e −i 2πaν F (ν ) 1 ν F a a F (ν ) ⋅ G (ν )
δ (ν )
1 sin(πν )
πν
16 / 96
Inverzió
a ⋅ f (t ) + b ⋅ g (t ) a ⋅ F (ν ) + b ⋅ G (ν ) a ⋅ F (ω ) + b ⋅ G (ω ) f (t − a )
dω
= −∞
e −iaω F (ω ) 1 ω F a a F (ω ) ⋅ G (ω ) 2πδ (ω ) 1 sin(ω / 2) ω/2
Az eredeti függvény visszaállítható a Fouriertranszformáltból:
f (t ) = =
∞
∫
x = −∞
=
∞
∞ −i 2πνx dx e i 2πνt dν = ∫ f ( x)e ∫ ν = −∞ x = −∞ ∞ f ( x) ∫ e −i 2πν ( x −t ) dν dx = ν =−∞
∞
∫ f ( x)δ ( x − t )dx = f (t )
x = −∞ 17 / 96
18 / 96
3
Eltolás és skálázás Eltolás:
Konvolúciós tétel ∞
∫
f (t − ∆t ) ⇔
Konvolúció a tértartományban:
f (t − ∆t )e −i 2πνt dt =
t = −∞ ∞
∫ f ( x )e
= Skálázás:
− i 2πν ( x + ∆t )
dx = F (ν ) ⋅ e −i 2πν∆t
fáziseltolás
∫ f (at )e
f (at ) =
−i 2πνt
=
ν
∫
f ( x)e
− i 2π x a
x = −∞
dt 1 ν dx = F dx a a
=
19 / 96
∞ if t = 0 0 otherwise
Dirac delta
δ (t ) =
∫ δ ( x)dx = 1
h (t ) =
∫
∞
∞
∞
∑ δ (t − k ) ⇔ ∫ ∑ δ (t − k ) ⋅ e
k = −∞
x = −∞
Impulzusválasz:
∞
Fésűfüggvény:
comb(t ) =
∞
∞ −i 2πνy dy ∫ f ( x) ⋅ h( y − x)dx e ∫ y = −∞ x = −∞
Mintavételezés
Impulzusválasz és átviteli függvény
∞
∞ −i 2πν ( x + z ) dy dz ∫ f ( x) ⋅ h( z )dx e dz z = −∞ x = −∞ ∞ ∞ = ∫ f ( x)e −i 2πνx dx ⋅ ∫ h( z )e −i 2πνz dz = F (ν ) ⋅ H (ν ) x = −∞ z =−∞ 20 / 96
dt =
t = −∞ ∞
Konvolúció a frekvenciatartományban:
R (ν ) =
∞
∫ f ( x) ⋅ h( y − x)dx
x = −∞
x = −∞
∞
r ( y ) = ( f ⊗ h)( y ) =
Mintavételezés:
∞
∫ δ ( x)h(t − x)dx
t = −∞ k = −∞
−i 2πνt
=
∞
∑e
−i 2πνk
=comb(ν )
k = −∞
1 t comb ⋅ f (t ) ⇔ comb(νT ) ⊗ F (ν ) T T diszkrét
periodikus
x = −∞
Frekvenciaválasz, vagy átviteli függvény:
H (ν ) =
∞
∫ δ (t )e
−i 2πνt
dt ⋅ H (ν )
t = −∞
21 / 96
Mintavételezés
DTFT: diszkrét idejű jel Fourier-transzformáltja
A mintavételezett jel Fourier transzformáltja:
∞
∞
∫ ∑ f δ (t − k )e
t = −∞ k = −∞
=
∞
∑fe
k = −∞
k
−i 2πνt
k
−i 2πνk
dt =
∞
∞
∑
k = −∞
fk
∫ δ (t − k )e
22 / 96
−i 2πνt
A mintavételezett jel spektruma az eredeti jel spektrumának periodikus ismétlődése
dt =
t = −∞
= DTFT[ f k ](ν )
A DTFT periodikus:
DTFT[ f k ](ν ) = DTFT[ f k ](ν + 1) 23 / 96
24 / 96
4
Ideális rekonstrukció
Sinc szűrő
Shannon-Nyquist kritérium: Az eredeti jel sávkorlátozott és a mintavételi frekvencia a sávkorlátnak legalább a kétszerese
Egy ideális aluláteresztő (sinc) szűrővel ki tudjuk vágni a középső spektrumot:
impulzusválasz
frekvenciaválasz
1 1 if | ν |< box(ν ) = 2 0 otherwise ∞
box(ν ) ⇔
∫ box(ν )e
i 2πtν
dν =
ν = −∞ 1 2
∫e
=
ν =−
i 2πtν
dν =
1 2
1 2
∫e
i 2πtν
dν =
1 ν =− 2
eiπt − e −iπt sin(πt ) = = sinc(t ) i 2πt πt
25 / 96
Ideális rekonstrukció
∑f
k = −∞
∞
f (t ) = f (t ) =
Inverz DTFT
Konvolúció a sinc szűrővel: ~
26 / 96
k
t ⋅ sinc( − k ) T
Inverz DTFT:
fk =
j = −∞
= box(ν ) ⋅
A sinc szűrő impulzusválasza végtelen kiterjedésű, ezért a gyakorlatban nem alkalmazható A sinc szűrőt véges impulzusválaszú szűrőkkel lehet közelíteni
∞
∑f
∞
∫ ∑f
t = −∞ j = −∞ ∞
∑fe
−i 2πνj
j
j = −∞
∞
sinc(k − j ) ⇔
j
j
sinc(t − j )e −i 2πνt dt =
1/ 2
∫ DTFT[ f
⇔
k
](ν )ei 2πνt dν
ν = −1/ 2
A DTFT-nek csak egy periódusát kell integrálni az eredeti diszkrét jel rekonstruálásához
27 / 96
Diszkrét Fourier transzformáció
Tegyük fel, hogy az eredeti jel periodikus:
Inverz DFT fk = fk +N
DFT:
N −1
DFT[ f k ]u = ∑ f k e −i 2πku / N
A DTFT nem csak periodikus, hanem diszkrét is:
DTFT[ f k ](ν ) = N −1
= ∑ f k e −i 2πνk k =0
28 / 96
∞
N −1
∑∑f
j = −∞ k = 0
∞
∑e
−i 2πνjN
j = −∞
k =0
−i 2πν ( k + jN ) = k + jN e N −1
= ∑ f k e −i 2πνk k =0
fk =
∞
∑ δ (νN − j )
j = −∞
A DFT N darab minta súlyozott összege:
DFT[ f k ]u = ∑ f k e −i 2πku / N k =0
N −1
∞
∑ DFT[ f
j = −∞
k
] j δ (νN − j )
1 N
= ∑ fj 29 / 96
1 N
N −1
∑ DFT[ f u =0
] ei 2πku / N
k u
Reprodukció:
fk =
N −1
DTFT[ f k ](ν ) =
Inverz DFT:
j =0
N −1 N −1
∑∑ f e u = 0 j =0
1 N
N −1
∑e u =0
j
−i 2πju / N i 2πku / N
e
i 2π ( k − j ) u / N
=
N −1
= ∑ f jδ k − j = f k j =0
30 / 96
5
Összefoglalás
Nevezetes szűrőcsaládok
tértartomány
frekvenciatartomány
folytonos függvény tisztán valós függvény tisztán képzetes függvény
folytonos függvény páros függvény páratlan függvény
DTFT
diszkrét függvény (végtelen sor)
folytonos és periodikus függvény
Fourier-sorba fejtés
folytonos és periodikus függvény
diszkrét függvény (végtelen sor)
DFT
diszkrét periodikus függvény
diszkrét periodikus függvény
Fourier-transzformáció
B-spline szűrők BC-spline szűrők Kardinális spline szűrők Catmull-Rom spline Nem polinomiális szűrők Ablakozott sinc szűrők
31 / 96
Smoothing vs. Postaliasing
Szűrők elemzése a frekvenciatartományban
Smoothing: (simítás): A szűrő frekvenciaválasza nem egy az áteresztő sávban (pass band) – finom részletek elkenése Postaliasing: A szűrő frekvenciaválasza nem nulla a levágó sávban (stop band) – a rekonstruált jelben fals frekvenciakomponensek is megjelennek
Prealiasing:
32 / 96
• Nem a szűrő tökéletlenségének a következménye • A nem megfelelő mintavételezés miatt a spektrum periodikus ismétlődései átlapolódnak 33 / 96
Ringing - Hullámzás
34 / 96
Hullámzás 2D-ben
Ha az eredeti jel nem sávkorlátozott, akkor az ideális sinc szűrő hullámzást hoz be (Gibbsjelenség) A drasztikus vágást a frekvenciatartományban célszerű elkerülni
35 / 96
36 / 96
6
B-spline szűrők családja
Elsőrendű B-spline szűrő
Generátorfüggvény (0-ad rendű B-spline):
1 2 0 otherwise
β 0 (t ) = 1
if | t |<
A magasabb rendű B-spline szűrőket rekurzívan definiáljuk:
β i +1 (t ) = β i (t ) ⊗ β 0 (t )
Az elsőrendű B-spline a lineáris sátor szűrő:
1− | t | if | t |< 1 otherwise 0
β 1 (t ) = β 0 (t ) ⊗ β 0 (t ) =
37 / 96
Harmadrendű B-spline szűrő
Elsőrendű B-spline szűrő levezetése ∞
β 1 (t ) = ∫ β 0 ( x) ⋅ β 0 (t − x)dx = −∞
38 / 96
integrálási tartomány: -1/2<x<1/2 -1/2
Köbös B-spline szűrő:
2 1 3 2 if | t |< 1 2 | t | − | t | + 3 1 4 3 1 1 3 2 β (t ) = β (t ) ⊗ β (t ) = − | t | + | t | −2 | t | + if 1 ≤| t |< 2 3 6 otherwise 0 A levezetés hasonló mint a lineáris B-spline esetén
1/ 2
=
∫β
0
(t − x )dx =
−1 / 2
t +1 / 2 if min(1 / 2 ,t +1 / 2 ) ∫ dx =t + 1 −1/ 2 = ∫ dx = 1/ 2 max( −1 / 2 ,t −1 / 2 ) dx =1 − t if t −1∫/ 2 )
t <0 t>0
1− | t | if | t |< 1 = otherwise 0 39 / 96
B-spline szűrők frekvenciaválasza
40 / 96
Köbös B-spline szűrő frekvenciaválasza
A box szűrő (0-ad rendű B-spline) Fourier-transzformáltja sinc(ν)
Az n-ed rendű B-spline szűrő a 0-ad rendű B-spline konvolúciója önmagával n-szer
Minden konvolúciónak megfelel egy szorzás a frekvenciatartományban
Az n-ed rendű B-spline szűrő Fourier transzformáltja sinc(ν)n+1 41 / 96
42 / 96
7
BC-Spline szűrők
BC-Spline szűrők levezetése
A BC-spline szűrők a szakaszonként harmadfokú szűrőknek egy kétparaméteres (B és C) családja:
(12 − 9 B − 6C ) | t |3 + if | t |< 1 2 (−18 + 12 B + 6C ) | t | + (6 − 2 B) 3 2 1 (− B − 6C ) | t | +(6 B + 30C ) | t | − ϕ BC (t ) = if 1 ≤| t |< 2 6 (−12 B − 48C ) | t | +(8B + 24C ) otherwise 0
Szakaszonként harmadfokú szűrők általános alakja nyolc paraméterrel:
P | t |3 +Q | t |2 + R | t | + S 1 ϕ BC (t ) = T | t |3 +U | t |2 +V | t | +W 6 0
Kardinális spline szűrők: interpolációs feltétel teljesül A köbös B-spline a BC spline speciális eseteként adódik B=1, C=0 paraméterezés mellett Catmull-Rom spline: B=0, C=0.5
if if
otherwise
Feltételek: • legyen a szűrő legalább C1 folytonos, azaz a határokon legyen folytonos átmenet a nulladrendű és elsőrendű deriváltak között • a frekvenciaválasz legyen nulla az aliasing spektrumok közepén – ne vigyük be a spektrum ismétlődéseinek DC komponensét a rekonstruált jelbe (ne legyen „sample frequency ripple”)
43 / 96
Kvadratikus konvergencia
44 / 96
Nem polinomiális szűrők
Ha 2C+B=1, akkor teljesül a kvadratikus konvergencia:
Csonkolt Gauss-szűrő
e −t / 2σ 0 2
~
2
ϕ Gauss (t ) =
| f (t ) − f (t ) |= (2C + B − 1)Tf ' (t )r (t ) + O(T ) 2
ahol r(t) egy polinomiális tag, T pedig a mintavételi távolság
| t |< 1 1 ≤| t |< 2
A köbös B-spline (B=1, C=0) és a Catmull-Rom spline (B=0, C=0.5) kielégíti ezt a feltételt
if | t |< b otherwise
Koszinuszos harangszűrő 1 + cos(πt / b) if | t |< b ϕCos (t ) = otherwise 0
Ablakozott (windowing) sinc szűrő
(1 + cos(πt / b)) sinc( 4t / b) if | t |< b 0 otherwise
ϕWindowedSi nc (t ) = 45 / 96
Smoothing és Postaliasing mérése
Smoothing és Postaliasing mérése
Smoothing:
Sϕ = 1 −
46 / 96
1 | ϕˆ (ν ) |2 dν ∫ R N | RN |
Postaliasing:
Pϕ =
1 | ϕˆ (ν ) |2 dν | RN | ∫RN
47 / 96
48 / 96
8
Smoothing és Postaliasing mérése
Kiterjesztés magasabb dimenzióra
Keresztszorzatos szeparálható kiterjesztés • • • •
2D: h(x,y) = h(x) * h(y) 3D: h(x,y,z) = h(x) * h(y) * h(z) Hatékonyan implementálható A kiterjesztett szűrő megőrzi az 1D szűrő előnyös tulajdonságait (pl. interpoláló jelleg, approximációs rend)
Kör- vagy gömbszimmetrikus kiterjesztés • • • •
2D: h(x,y) = h(sqrt(x2+y2)) 3D: h(x,y,z) = h(sqrt(x2+y2+z2)) Kevésbé irányfüggő, de költségesebb kiértékelni A kiterjesztett szűrő nem őrzi meg az 1D szűrő előnyös tulajdonságait (pl. interpoláló jelleg, approximációs rend)
49 / 96
Anisotropy - irányfüggőség
Konvolúció magasabb dimenzióban ~
f ( x, y , z ) ≈ f ( x, y , z ) = =
∞
∞
∞
∑∑ ∑
i = −∞ j = −∞ k = −∞
f i , j , kϕ (
Legközelebbi szomszéd
50 / 96
Szeparálható kiterjesztés • az áteresztő sávban lesznek preferált irányok • a lezáró sávban nem lesz egyenletes a postaliasing hatás elnyomása a különböző irányok mentén • irányfüggő artifaktumok jelennek meg a rekonstrukcióban
Sugarasan szimmetrikus kiterjesztés • irányfüggetlen az áteresztő sávban • irányfüggő a lezáró sávban, mivel a spektrum a duális rács rácspontjai körül ismétlődik
x y z − i, − j , − k ) Tx Ty Tz
Lineáris
Köbös
Figyelembe vett 2D pixelek száma
1
4
16
Figyelembe vett 3D voxelek száma
1
8
64 51 / 96
Vizuális összehasonlítás
52 / 96
Trilineáris interpoláció
Marschner-Lobb benchmark
53 / 96
A lineáris B-spline szeparálható kiterjesztése Erős postaliasing hatás de nem nagyon mossa el a nagyfrekvenciás részleteket
54 / 96
9
Triköbös B-spline szűrés
Trilineáris interpoláció kiértékelése
A nyolc legközelebbi rácspont értékeinek súlyozott összege
A köbös B-spline szűrő szeparálható kiterjesztése Nagyon elmossa (blurring) a nagyfrekvenciás részleteket Jól elnyomja a postaliasing hatást
55 / 96
Catmull-Rom spline
56 / 96
Más BC-spline szűrők
Interpoláló (kardinális BC-spline, B=0, C=0.5) Kvadratikus konvergencia (2C+B=1) Jó kompromisszum a postaliasing és a smoothing között
Előfordulhatnak kockás és hullámzó artifaktumok
B=0.5, C=0.85
B=0.26, C=0.1
57 / 96
Ablakozott sinc szűrő
58 / 96
Osztályozás aszimptotikus hiba szerint
Koszinuszos ablakozó függvénnyel (r=4.8)
59 / 96
Feltételezzük, hogy az eredeti jel kellőképpen sima, nincsenek hirtelen ugrások a magasabb rendű deriváltakban sem
Hasonló a sávkorlátozottság feltételezéséhez
A rekonstruált jel minél gyorsabban közelítse az eredetit, ha a mintavételi távolság tart nullához
Tértartománybeli kritériumok a rekonstrukciós illetve derivált szűrők tervezéséhez 60 / 96
10
Taylor sor szerinti kifejtés
Közelítés konvolúciós szűréssel: ~
f (t ) ≈ f (t ) =
∞
∑
k = −∞
Osztályozás aszimptotikus hiba szerint
f kϕ (
t − k) T
~
f (t ) = a0ϕ (t ) f (t ) + a1ϕ (t ) f ' (t ) + a2ϕ (t ) f '' (t ) + ...
Taylor-sor szerinti kifejtés a kT pont körül
fk =
N
f
∑
(n)
(t ) ( kT − t ) n + O ( T n!
n=0
~
f (t ) =
N
∑
n=0
a nϕ ( t ) =
a nϕ ( t ) f
T n n!
∞
∑
(n)
(k −
k = −∞
( t ) + O (T
N +1
N +1
)
)
t n t ) ϕ( − k) T T
A szűrőre jellemző hibaegyüttható
Függvényrekonstrukció: • a0 = 1 • a1 , a2 , a3 … aN = 0 ⇒ a hibafüggvény rendje O(TN+1) • L-EF rekonstrukció (L = N+1) • L-1 fokú polinomot tökéletesen rekonstruál • L-ed rendű kvázi-interpoláció
61 / 96
Deriváltszűrés aszimptotikus hiba szerint
Taylor-sor szerinti kifejtés:
62 / 96
Osztályozás aszimptotikus hiba szerint
Taylor-sor szerinti kifejtés: ~
f (t ) = a0ϕ (t ) f (t ) + a1ϕ (t ) f ' (t ) + a2ϕ (t ) f '' (t ) + ...
Box szűrő: 1-EF (nem folytonos) Lineáris B-spline: C0 2-EF Köbös B-spline: C2 2-EF Catmull-Rom spline: C1 3-EF
Deriváltszűrés: • a0 = 0 • a1 = 1 • a2 , a3 … aN = 0 ⇒ a hibafüggvény rendje O(TN+1) • L-EF deriváltszűrés (L = N+1) 63 / 96
Szűrők tervezése
64 / 96
Szűrők tervezése
Tegyük fel, hogy k-adik deriváltat akarunk számolni L-EF pontossággal Cm folytonossággal 1. feltétel: an = 0 minden n
65 / 96
A wk(τ) tagokat polinomiális alakban keressük Az 1-3 feltételek egy lineáris egyenletrendszert alkotnak Az egyenletrendszer ismeretlenjei a polinomok együtthatói Ha nem oldható meg az adott feladat, akkor növelni kell a polinomok fokszámát illetve a szűrő kiterjedését (a wk(τ) tagok számát) 66 / 96
11
Nevezetes szűrők reprodukciója
Hibaegyütthatók Fourier-transzformáltja
A szűrőtervező algoritmus speciális esetként reprodukál nevezetes szűrőket:
A hibaegyütthatók tértartományban:
anϕ (t ) =
• C2 2-EF, kiterjedés [-2,2]: Köbös B-spline
Tn n!
∞
t
t
∑ (k − T ) ϕ ( T − k ) n
k = −∞
A hibaegyütthatók a frekvenciatartományban:
• C1 3-EF, kiterjedés [-2,2]: Catmull-Rom spline
Anϕ (ν ) =
Tn comb(ν )(−i ) n Φ ( n ) (ν ) n!
67 / 96
Előszűrt rekonstrukció
Hibaegyütthatók Fourier-transzformáltja
68 / 96
L-EF függvényrekonstrukció (L-ed rendű kvázi-interpoláció): • a0 = 1 • a1 , a2 , a3 … aL-1 = 0
Diszkrét előszűrés – megváltoztatjuk az eredeti mintákat Folytonos utószűrés egy rekonstrukciós kernellel
A frekvenciatartományban: • A0 = δ(ν) ⇒ Φ(k) = δk • A1 , A2 , A3 … AL-1= 0 ⇒ Φ(n)(k) = 0 minden 0
Előszűrt rekonstrukció
70 / 96
Hagyományos interpoláció
Általánosított interpoláció – egy nem interpoláló szűrőt teszünk interpolálóvá Kvázi-interpoláció – az előszűrőt úgy optimalizáljuk, hogy a legmagasabb rendű kváziinterpolációt kapjuk Áteresztő sávban optimális rekonstrukció – az előszűrő egy olyan eredő frekvenciaválaszt eredményez, amely ideális az áteresztő sávban 71 / 96
Hagyományos interpoláció: ~
f (t ) ≈ f (t ) =
∞
∑
k = −∞
f kϕ (
t − k) T
A rekonstrukciós szűrő teljesíti az interpolációs feltételt:
ϕ (k ) = δ k ~
f ( jT ) =
∞
∑
k = −∞
f kϕ ( j − k ) = f
j 72 / 96
12
Általánosított interpoláció
Általánosított interpoláció
Hagyományos interpoláció:
~
∞
f ( jT ) =
t f (t ) ≈ f (t ) = ∑ c kϕ ( − k ) T k = −∞ ~
A rekonstrukciós szűrő nem teljesíti az interpolációs feltételt:
∞
f ( jT ) =
∑
k = −∞
c kϕ ( j − k ) = f
Általánosított interpoláció
Feltételezzük, hogy ck és fk periodikusak N-ben Ekkor az együttható mátrix ciklikus lesz:
N −1 k =0
ckϕ ( j − k ) = f
j
Toeplitz mátrix: ciklikus sávmátrix, iteratív módszerekkel invertálható
74 / 96
Konvolúció: N −1
∑
k =0
ckϕ ( j − k ) = f
= c⊗ϕ
−1
A dekonvolúció egy osztás a frekvenciatartományban
DFT [ c k ] j =
= c⊗ϕ
j
Dekonvolúció:
c = f ⊗ϕ
Az együttható mátrixszal való szorzás valójában egy ciklikus konvolúció:
∑
Lineáris egyenletrendszer a ck ismeretlenekre:
Megoldás a frekvenciatartományban
ϕ (1) c0 f 0 L ϕ (0) ϕ (−1) ϕ (1) ϕ (0) ϕ (−1) c f ⋅ 1 = 1 M ϕ ( −1) M M O ϕ (1) ϕ (0) cN −1 f N −1 ϕ ( −1)
j
j 73 / 96
c kϕ ( j − k ) = f
L ϕ (0) ϕ (−1) c0 f 0 ϕ (1) ϕ (0) ϕ (−1) c f ⋅ 1 = 1 M ϕ ( −1) M M O ϕ (1) ϕ (0) c N −1 f N −1
A ck együtthatókat úgy kell meghatározni, hogy ~
∞
∑
k = −∞
ϕ (k ) ≠ δ k
A ck együtthatókat úgy kell meghatározni, hogy
DFT [ f k ] j DFT [ ϕ k ] j
75 / 96
Példa: Köbös B-spline
Nem interpoláló szűrő:
β 3 (0 ) =
2 , 3
Eredő frekvenciaválasz
β 3 ( ± 1) =
Az előszűrő DTFT-je:
Az előszűrt rekonstrukció eredő frekvenciaválasza:
DTFT[ pk ](ν ) =
W (ω ) = DTFT[ pk ](ν ) β 3 (ν ) =
76 / 96
1 6
Az előszűrés csökkenti a simítást – a finom részleteket nem mossuk el annyira A postaliasing hatás viszont növekszik
1 2 1 + cos(ν ) 3 3
sinc 4 (ν ) 3 sinc 4 (ν ) = 2 1 + cos(ν ) 2 + cos(ν ) 3 3 77 / 96
78 / 96
13
Eredő impulzusválasz
Köbös B-spline előszűrés nélkül
Az előszűrés vizuális hatása
eredeti jel
Köbös B-spline előszűréssel 79 / 96
Implementáció
80 / 96
Kvázi-interpoláció
1.
Az fk eredeti minták DFT-je
2.
Előszűrés a frekvenciatartományban: DFT[ck]j= DFT[fk]j DTFT[φk](j/N)
3.
DFT[ck]j inverz DFT-je
4.
A ck mintákat konvolváljuk a φ szűrővel
L-ed rendű kvázi-interpoláció: legfeljebb L-1 fokú polinomot tökéletesen reprodukál
Feltételek a frekvenciatartományban: • Strang-Fix kritériumok (szükséges feltétel) : • Φ (k) = δk • Φ (n)(k) = 0 minden 0
Előszűrő tervezése rögzített rekonstrukciós kernelhez úgy, hogy a feltételek teljesüljenek
81 / 96
Példa: Köbös B-spline
82 / 96
FIR előszűrő tervezése
Approximációs rend: L = 4 (Strang-Fix)
Előszűrés nélkül csak másodrendű kvázi-interpoláció - a szűrő approximációs rendje nincs teljesen kihasználva
Eredő frekvenciaválasz:
Három pontos diszkrét előszűrő, ismeretlenek: p0, p1 Az előszűrő DTFT-je:
DTFT[ pk ](ν ) = p0 + 2 p1 cos(2πν ) pk p0
∞
W (ν ) = β 3 (ν ) DTFT[ pk ](ν ) = sinc 4 (ν ) ∑ pk e −i 2πνk k = −∞
A pk súlyok az ismeretlenek, ezeket kell úgy beállítani, hogy a kvázi-interpoláció feltételei teljesüljenek
p1
83 / 96
p1
84 / 96
14
FIR előszűrő tervezése
Eredő frekvenciaválasz
1. feltétel: W(k) = δk 2. feltétel: W(n)(k) = 0 minden 0
Az előszűrőre a következő feltételeket kapjuk • [ p0 + 2 p1 cos(2πk )] sinc 4 (k ) = δ k • [[ p + 2 p cos(2πk )] sinc 4 (k )]( n ) = 0 0 1
A megoldás L=4-re: p0=4/3, p1=-1/6 85 / 96
Az előszűrés vizuális hatása köbös B-spline
Catmull-Rom spline
86 / 96
IIR előszűrő tervezése előszűrt köbös B-spline
Három pontos diszkrét dekonvolúciós előszűrő, ismeretlenek: q0, q1 A dekonvolúciós előszűrő DTFT-je:
DTFT[qk ](ν ) = q0 + 2q1 cos(2πν ) qk q0
Felbontás: 40×40×40
q1 Felbontás: 60×60×60
87 / 96
IIR előszűrő tervezése
88 / 96
FIR versus IIR előszűrés
1. feltétel: W(k) = δk 2. feltétel: W(n)(k) = 0 minden 0
FIR • A tértartományban is implementálható • Nem feltételez periodikus kiterjesztést • Kisebb mértékben javítja a frekvenciaválaszt az áteresztő sávban • Kisebb a postaliasing hatása
IIR • Csak a frekvenciatartományban implementálható • Periodikus vagy tükrös kiterjesztést feltételez • Drasztikusan javítja a frekvenciaválaszt az áteresztő sávban • Nagyobb a postaliasing hatása 90 / 96
1 sinc 4 (k ) = δ k q0 + 2q1 cos( 2πk ) 1 sinc 4 (k ) q0 + 2q1 cos(2πk )
q1
( n)
=0
A megoldás L=4-re ugyanaz, mint az általánosított interpoláció esetén: q0=2/3, q1=1/6 89 / 96
15
Implementáció
Áteresztő sávban optimális rekonstrukció
Az előszűrés után az áteresztő sávban egy lesz az eredő frekvenciaválasz Az előszűrő diszkrét, ezért a DTFT-je periodikus kell, hogy legyen A [-0.5, 0.5] intervallumon értelmezett periódus a rekonstrukciós szűrő Fourier transzformáltjának reciproka Az előszűrő DTFT-je:
1 DTFT[ pk ](ν ) = box(ν ) ⊗ comb(ν ) Φ (ν )
1.
Az fk eredeti minták DFT-je
2.
Előszűrés a frekvenciatartományban: DFT[ck]j= DFT[fk]j DTFT[pk](j/N)
3.
DFT[ck]j inverz DFT-je
4.
A ck mintákat konvolváljuk a φ szűrővel
91 / 96
Eredő frekvenciaválasz
Áteresztő sávban optimális előszűrt B-spline rekonstrukció
Relatíve erős postaliasing hatás
92 / 96
Az előszűrés vizuális hatása
ringing
Drasztikus vágás ⇒ ringing Minél magasabb rendű szűrőt használunk annál inkább csökken a postaliasing hatás
ringing tükrös kiterjesztéssel
93 / 96
Tesztelés CT adaton
periodikus kiterjesztéssel 94 / 96
Összefoglalás
95 / 96
A diszkrét előszűréssel felülírjuk az eredeti adatokat Előfeldolgozás – az újramintavételezés költségét nem növeli FIR vagy IIR előszűrés Az IIR előszűrés DFT-vel implementálható Mindegyik előszűrő módszer javítja frekvenciaválaszt az áteresztő sávban de rontja a lezáró sávban Minél magasabb rendű rekonstrukciós kernelt használunk, annál érdemesebb diszkrét előszűrést alkalmazni
96 / 96
16
