Folyamatosan változó mennyiségek feldolgozása II. 7. előadás
1
Tartalom • • • • • •
sztochasztikus folyamatok mintavételezés (lásd Fizikai geod.) LKN kollokáció (lásd Fizikai geod.) geostatisztika, krigelés szűrések interpolációk 2
Szűrések – Detrekői 7.5 • célja – a hibával terheltnek feltételezett jelek „szűrt” értékeinek meghatározása a mért pontokban – esetenként a jelek értékének becslése nem mért pontokban – zajfüggvény leválasztása az alapfolyamatról • feltételezzük, hogy a trend zérus
• típusok: – – – –
konvolúciós szűrők, medián szűrők, legkisebb négyzetek módszerén alapuló szűrők, Kálmán-szűrők. 3
A szűrések spektrális jellemzői • Lineáris szűrést feltételezve: η (t) = O (t) ξ (t) – ξ (t) eredeti (nem szűrt) függvény – η (t) kimeneti (szűrt) függvény – O (t) a szűrés módját előíró operátor.
• Fourier transzformációval áttérve az ω frekvencia tartományba a szűrés hatását a C(ω) átviteli függvény írja le. – Az átviteli függvény írja le az amplitúdó (és ha C komplex, a fázis) viselkedését a szűrés során. 4
Konvolúciós szűrők • Az egy és kétváltozós konvolúciós szűrés alapösszefüggései: – általános esetben, – szimmetrikusan elhelyezkedő jelek esetében.
• A konvolúciós szűrők fajtái: – aluláteresztő szűrők, – felüláteresztő szűrők – sáváteresztő szűrők5
Konvolúciós szűrők – aluláteresztő szűrők – felüláteresztő szűrők – sáváteresztő szűrők – sávkorlátozó szűrők AF aluláteresztő
OK
0
AF
AF
fL
fH
OK
stop
OK
stop
0
sávkorlátozó
sáváteresztő
felüláteresztő
stop
AF
0
fL
OK
stop
fH
0
OK
stop
fL
fH
6
Egyváltozós aluláteresztő szűrő • ismertek az adott t1, t2, . ..‚ tN pontokhoz tartozó s'1, s'2, . ..‚ s'N értékek (hibával terheltek): s'i = si + ni • az si szűrt értéket a szimmetrikusan elhelyezkedő 2M + 1 jel súlyozott számtani közepeként állítjuk elő: si = Σ(wksi+k), k Î (–M, +M) • ez a (diszkrét) konvolúció 7
Súlyok felvétele • futó átlagolás: wk = 1/(2M + 1) • simító átlagolás: wk = w-k = ((M + 1 – |k|)/ (M + 1)2 • pl. ha M = 2, futó átlagolás: 0.2, 0.2, 0.2, 0.2, 0.2 simító átlagolás: 0.111, 0.222, 0.333, 0.222, 0.111
• kiugró értékek kimutatása: | s'i – si | > e 8
Néhány gyakran alkalmazott kétváltozós konvolúciós szűrés • átlagoló szűrők, • egyéb aluláteresztő szűrők, • felüláteresztő szűrők (Laplace-féle szűrő, LoG, élkiemelés).
9
Medián szűrők • zajszűrés, élmegtartás (alacsony zajszint esetén) • ismertek az adott t1, t2, . ..‚ tN pontokhoz tartozó s'1, s'2, . ..‚ s'N értékek (hibával terheltek): s'i = si + ni • az si szűrt értéket a szimmetrikusan elhelyezkedő 2M + 1 jel mediánjaként állítjuk elő: si = medián(si+k), k Î (–M, +M) • Például x = [2 80 6 3] • szűrés: y[1] = medián[2 2 80] = 2 y[2] = medián[2 80 6] = medián[2 6 80] = 6 y[3] = medián[80 6 3] = medián[3 6 80] = 6 y[4] = medián[6 3 3] = medián[3 3 6] = 3 • tehát y = [2 6 6 3]. 10
Medián szűrő, példa
eredeti kép
M=1px medián szűrő
11 M=3px medián szűrő
M=10px medián szűrő
A legkisebb négyzetek módszerén alapuló szűrők • A legkisebb négyzetek módszerén alapuló szűrés a legkisebb négyzetek elvén alapuló kollokáció speciális esete (trend= 0)
12
A Kálmán-szűrő • Rudolf Emil Kalman publikálja az elméletet (1960) • Első nevezetes alkalmazása: Apolló űrhajó fedélzeti számítógépe (Armstrong) • Navigációs és mérnöki alkalmazások • Szenzor adatok fúziója (radar, lézerszkenner, kamera, INS,…) • GNSS vevők, okostelefonok, számítógépek, játékok, stb. elmaradhatatlan része • Tőzsde előrejelzés, idősor elemzés, függvény approximáció Rudolf Emil Kalman (Kálmán Rudolf) (1930-) magyar származású villamosmérnök
13
A Kálmán-szűrés definíciója • „Olyan algoritmus, amely valamely lineáris dinamikus rendszerben egzakt következetést tesz lehetővé, amely a rejtett Markov-modellhez hasonló Bayesféle modell, azonban a rejtett változók állapottere folytonos, és minden rejtett és megfigyelhető változó (gyakran többváltozós) normális eloszlású.” • UGYE MINDEN VILÁGOS!? 14
Meg lehet érteni? • a Kálmán-szűrés elgondolása egészen egyszerű és közvetlen: • időben változó rendszer, amelynek állapotát leíró változókat (x) nem (teljesen) ismerjük • ezt az állapotot szeretnénk minden lépésben a lehető legjobban megismerni • rendelkezünk folyamatos mérésekkel (y), melyek kapcsolatba hozhatók a rendszer állapotával (x) 15
A Kálmán szűrő jellemzői • lépésenkénti számítás, csak az előző és a jelenlegi lépés függvényében = rekurzív • a kapcsolatok lineárisak (állapot változása, mérések-állapot kapcsolata) • sztochasztikus (mérési hibák, állapot és a hibája valószínűségi változók) • legkisebb hibájú optimális becslés (csak lineáris esetben bizonyítható) 16
Rekurzív becslés – példa • egyváltozós rendszer állapotát mérjük egymás utáni időpontokban • mérések: x1, x2, … xn, • a legjobb állapot becslés az átlag: 1 n m n = å xi n i =1
• új mérésünk van: xn+1, a következő becslés (számológépek STAT funkciói): 1 n +1 n æ1 n 1 ö m n +1 = xi = ç å xi + xn +1 ÷ å n + 1 i =1 n + 1 è n i =1 n ø
17
Erősítés, új információ • a frissített becslés kis korrekcióval számítható: n 1 m n +1 = mn + xn +1 = m n + K ( xn +1 - m n ) n +1 n +1
• K az erősítési tényező („nyereség”) – ez dönti el, mekkora lesz a korrekció
• xn+1 – µn az új információ („innováció”) – ha nincs új információ, nem változik az előző becslés 18
Rekurzív átlag becslés
19
Különböző információk optimális felhasználása • egy távolságot két különböző eszközzel mértünk (mérőszalag, lézer távmérő) • mérések: x1, x2, szórások: σ1, σ2 • az optimális becslést keressük x-re –ez a súlyozott átlag:
xˆ = wx1 + (1 - w) x2
• Mi az optimális w súly? 20
Optimális súly meghatározása • mivel x1, x2 normális eloszlású, a lineáris kombinációjuk is normális eloszlású, melynek szórása 2 2 2 2 2
sˆ = w s 1 + (1 - w) s 2
• az optimális w súly ennek a kifejezésnek a minimuma, ahol a w szerinti derivált zérus lesz 2
¶sˆ 2 2 = 2wopts 1 - 2(1 - wopt )s 2 = 0 ¶w
• az optimális súly tehát
wopt
s = 2 2 s1 + s 2 2 2
21
Rekurzív optimális becslés • Két lépésben mérünk, először x1-et, azután x2 -t. Az első lépésben 2 2 ˆ s1 = s1 xˆ1 = x1 • A második lépésben frissítjük a becslést sˆ12 új információ ˆ ( ) xˆ2 = xˆ1 + 2 x x 2 1 2 (innováció) sˆ1 + s 2 2 æ ö 2 ˆ s 2 1 ÷ sˆ sˆ 2 = çç1 - 2 2 ÷ 1 è sˆ1 + s 2 ø
K erősítés (nyereség)
22
Rekurzív optimális becslés • a becslési egyenletek xˆ2 = xˆ1 + K ( x2 - xˆ1 ) 2 2 ˆ ˆ s 2 = (1 - K )s 1
sˆ12 K= 2 sˆ1 + s 22
becslés hibája (szórás)
K: erősítés
• ismerős az egyenlet alakja? m n +1 = m n + K (xn +1 - m n )
23
Többváltozós eset • egy változó esetén: varianciák sˆ12 K= 2 sˆ1 + s 22
K: erősítés (szám)
• több változó esetén: kovariancia mátrixok -1 ˆ ˆ ˆ K = M1 ( M1 + M 2 )
K: erősítés (mátrix) (nyereségmátrix)
24
Erősítés mértéke • ha jó a mérés, kicsi a σ2 hiba (szórás): sˆ 12 sˆ 12 K= 2 » 2 =1 2 sˆ 1 + s 2 sˆ 1 + 0
egységnyi erősítés
• következménye: xˆ2 = xˆ1 + K ( x2 - xˆ1 ) » xˆ1 + ( x2 - xˆ1 ) = x2 a második mérés lesz az új becslés 25
Erősítés mértéke • ha rossz a mérés, nagy a σ2 hiba (szórás): ˆs12 ˆs12 K= 2 » 2 =0 2 sˆ 1 + s 2 sˆ 1 + ¥
zérus erősítés
• következménye: xˆ2 = xˆ1 + K ( x2 - xˆ1 ) » xˆ1 + 0 = xˆ1 az első mérés lesz az új becslés 26
Közvetett mérések • Mi van, ha x-et nem tudjuk közvetlenül mérni (rejtett változó, rendszerállapot)? • Geodéziában ez jól ismert: közvetett mérések kiegyenlítése (II. kiegyenlítési csoport, paraméteres kiegyenlítés), Gauss-Markov modell • Lineáris (közvetítő) egyenlet az y mérés és x rendszerállapot között (z a mérési hiba, „zaj”):
y = Cx + z
• a mérési zajt zérus átlag ( z = 0 ) és az Sz kovariancia mátrix jellemzi (E[ ] a várható érték):
[
S z = E ( z - z )( z - z )T
]
27
Kálmán szűrők elve Állapotegyenlet: xk+1 = Axk + Buk + wk Átviteli egyenlet: yk = Cxk + zk A, B, C mátrixok; k az időváltozó (index) x a rendszer állapota (közvetlenül nem mérhető) u a rendszer ismert bemenete y a mért kimenete w és z a zaj w folyamat zaj z-t a mérési zaj
28
Rekurzív állapotbecslés állapot terjedése
korrekció
xˆ k +1 = ( Axˆ k + Bu k ) + K k ( y k +1 - Cxˆ k ) K mátrix: Kálmán-féle nyereségmátrix 29
Állapot hiba terjedése Átviteli egyenlet: yk = Cxk + zk Hibaterjedés a mért kimenetre Py = CPkCT+ Sz C mátrix; k az időváltozó (index) x a rendszer állapota (közvetlenül nem mérhető) y a mért rendszer kimenet Pk a becslési hiba kovariancia mátrixa Py a rendszer kimenet kovariancia mátrixa Sz a mérési zaj kovariancia mátrixa
30
K nyereségmátrix Ha nagy a mérési zaj, Sz is nagy lesz, tehát K kicsi, és így nem adunk nagy hitelt az y méréseinknek a következő állapot becslésénél. Viszont ha kicsi a mérési zaj, Sz is kicsi lesz, tehát K nagy, és így több hitelt adunk a méréseinknek a következő állapot becslésénél.
K k = APk C (CPk C + S z ) T
T
-1 31
Kálmán szűrők elve Folyamat zaj kovariancia: Sw = Mérési zaj kovariancia: T S z = E( zk zk )
T E ( wk wk )
Kálmán-szűrő algoritmusa: -1
K k = APk C (CPk C + S z ) xˆ k +1 = ( Axˆ k + Bu k ) + K k ( y k +1 - Cxˆ k ) T
T
Pk +1 = APk A + S w - APk C T
T
-1 T S z CPk A
32
P kovariancia terjedés
A kovariancia terjedés összefüggése: -1 z
Pk +1 = APk A + S w - APk C S CPk A T
T
T
33
Kálmán szűrés
34
Kálmán-szűrő alkalmazásai • • • • • • • • • • • • •
objektumkövetés, mesterséges látás dinamikus helymeghatározás gazdaságtan, makroökonómia inerciális navigáció pályameghatározás radar követés GNSS szeizmológia beszédjavítás időjárás előrejelzés navigációs rendszerek 3D modellezés stb. 35
Példa: járműnavigáció
piros : valódi helyzet kék : mért helyzet zöld : a Kálmán-szűrő által becsült helyzet
36
Interpolációk – Detrekői 7.4 • Célja: a jelek meghatározása olyan pontokban, ahol jel nem áll rendelkezésünkre • trend = zérus • a felhasznált jeleket nem terheli zaj • típusai: – – – – – – –
interpoláció súlyozott számtani középként, legkisebb négyzetek módszerén alapuló interpoláció, legközelebbi szomszéd módszere, lineáris és bilineáris interpoláció, spline-interpolációk polinomiális interpolációk klasszikus és általánosított interpoláció • eltolt lineáris interpoláció 37
Interpoláció súlyozott számtani középként • Az ismeretlen ponthoz tartozó jelet az ismert pontokhoz tartozó jelek (általában a távolságoktól függő) súlyozott számtani közepeként számítjuk. 1D, 2D, 3D esetben egyaránt alkalmazható.
38
A legkisebb négyzetek módszerén alapuló interpoláció • A legkisebb négyzetek módszerén alapuló interpoláció a legkisebb négyzetek elvén alapuló kollokáció speciális esete (trend= 0, zaj=0).
39
Interpoláció feladata • alappontok: xk,, k = 0, 1, ..., n ahol a függvény értékei adottak: fk := f(xk) • approximáció: a közelítő függvény nem halad át az alappontokon • interpoláció: a közelítő függvény áthalad az alappontokon 40
A legközelebbi szomszéd módszere (lépcsős interpoláció) • az ismeretlen ponthoz tartozó jel értékét a hozzá legközelebbi ismert jelű pont jel értékével tekintjük azonosnak. • 1D, 2D, 3D esetben egyaránt alkalmazható.
41
A legközelebbi szomszéd módszere (lépcsős interpoláció) • 12-szeres elforgatás, interpoláció
42
Lineáris interpoláció • számítási összefüggések • Példa az alkalmazásra.
43
Bilineáris interpoláció • számítási összefüggések • Példa az alkalmazásra.
44
Bilineáris interpoláció • 12-szeres elforgatás, interpoláció
45
Spline interpoláció • A spline interpoláció elve. • Egyváltozós harmadfokú spline interpoláció. • Példák az alkalmazásokra.
46
Spline interpoláció elve • A módszer lényege, hogy szakaszosan, több alacsony fokszámú polinommal közelítjük az adatokat. Például harmadfokú spline esetén,
si
3
2
a i× x + b i× x + c i× x + d i
• Az együtthatókat a csatlakozási, simasági és végfeltételekből határozzuk meg
47
Egyváltozós köbös spline • Harmadfokú spline esetén a megoldandó egyenletrendszer a Di deriváltakra (m + 1 adatpont van): é2 1 ù é D0 ù é 3 ( y1 - y0 ) ù ê1 4 1 ú ê D ú ê 3( y - y ) ú 2 0 ú ê úê 1ú ê ú ê 1 4 1 úê . ú ê . ú ê úê ú ê = . . 1 4 1 ú ê úê ú ê ú ê úê . ú ê . . . . . ú ê úê ú ê 1 4 1ú ê . ú ê3 ( ym - ym- 2 )ú ê ê 1 2úû êë Dm úû êë 3 ( ym - ym-1 ) úû ë
• m+1 egyenlet m+1 ismeretlenre ai = y i • Ha megkaptuk a Di-ket, ezekből ai, bi = Di bi, ci, di könnyen számítható: ci = 3( y i +1 - y i ) - 2 Di - Di +1 d i = 2( y i - y i +1 ) + Di + Di +1 48
Köbös spline interpoláció • 12-szeres elforgatás, interpoláció
49
Interpoláció polinommal • interpolációs polinom: olyan p(x) , max. n-edfokú polinom amelyre teljesül, hogy p(xk) = fk , k = 0, 1, ..., n • Tétel: Pontosan egyetlen olyan polinom létezik, amelynek fokszáma ≤ n és az n + 1 különböző xk pontban az adott fk értékeket veszi fel (k = 0, 1, ..., n). 50
Lagrange-féle interpolációs polinom • speciális eset: fm = 1, az összes többi interpolált érték 0 lm(xk) = 1, ha k = m és lm(xk) = 0, ha k ≠ m gyöktényezős alakban: lm(x) = cm(x-x0)(x-x1)···(x-xm-1)(x-xm+1) ···(x-xn)
• a c állandó meghatározása: 1 = cm(xm-x0)···(xm-xm-1)(xm-xm+1) ···(xm-xn) • a keresett Lagrange-féle interpolációs polinom:
x - xk lm ( x ) = Õ k = 0 xm - xk n
k ¹m
51
Lagrange-féle interpoláció • interpoláció a Lagrange-féle polinommal n
f int ( x) = å lm ( x) f m m =0
• interpolációs feltétel n
f int ( xk ) = å lm ( xk ) f m = f k
k = 0, 1, ... , n
m =0
• akkor teljesül, ha ì1, k = m l m ( xk ) = í î0, k ¹ m 52
53
Interpolációk osztályozása • Klasszikus interpoláció • Általánosított interpoláció – Eltolt lineáris interpoláció
54
Klasszikus interpoláció • interpoláció a φint(x) magfüggvénnyel f int ( x) = å f n jint ( x - xn ) nÎZ
• interpolációs feltétel f int ( xk ) = å f njint ( xk - xn ) = f k
k = 0, 1, ... , n
nÎZ
• akkor teljesül, ha ì1, k = n jint ( xk - xn ) = í î0, k ¹ n 55
Egyenközű klasszikus interpoláció • egyenlő, egységnyi xi = i -k f k = å f i jint (k - i ) iÎZ
• diszkrét konvolúció f k = ( f * p) k
k ÎZ
pk = jint (k )
• interpolációs feltételből következik, hogy ì1, k = 0 pk = d k = í î0, k ¹ 0 56
Általánosított interpoláció • interpoláció a φ(x) bázisfüggvénnyel f int ( x) = å cn j ( x - xn ) nÎZ
• interpolációs feltétel f int ( xk ) = å cnj ( xk - xn ) = f k
k = 0, 1, ... , n
nÎZ
• előny: most a bázisfüggvénynek nem kell kielégítenie az interpolációs feltételt, ezért kedvezőbben választható meg • hátrány: több számítás 57 Hivatkozás: Thévenaz et al.(2000)
Kapcsolat a klasszikus interpolációval • klasszikus interpoláció: f int ( x) = å f njint ( x - n) k Î Z nÎZ
• interpoláló és nem interpoláló magfüggvény jint ( x) = å hnj ( x - n) nÎZ
• az ekvivalens interpoláló magfüggvényt a bázisfüggvények szerinti szűréssel kapjuk meg 58
Példa: köbös spline köbös B-spline (nem interpoláló) -> Gauss-féle függvény 1/6, 4/6, 1/6 sorozat
kardinális köbös spline (interpoláló) -> sinc függvény
59
Lineáris interpoláció másként • φ bázisfüggvény: háromszög vagy „kalap” Λ(x)
60
Lineáris interpoláció • átskálázott és egész értékkel eltolt „kalap” függvények összege f int ( x) = å f n L ( x - xn ) nÎZ
61
Eltolt lineáris interpoláció • nem interpoláló Λ bázisfüggvénnyel (T mintavételezési köz, 0 < τ ≤ 0.5 eltolás) æx ö f int ( x) = å cn Lç - n - t ÷ nÎZ èT ø
• rekurziós egyenlet f n = cn (1 - t ) + cn-1t
62
Interpoláció konstrukciója eltolt lineáris
lineáris
fn-1
cn-1
τT
fn
(1-τ)T
cn
τT
fn+1
(1-τ)T
cn+1
τT
63
Az interpoláló bázisfüggvény
1
64
c együtthatók számítása • cn helyben számítható (fn felülírható): t 1 cn = cn -1 + fn 1-t 1-t
• kiindulás: c0 = c-1 = ... = f 0
• cn számítása hatékony, rekurziós szűréssel lehetséges
65
Példa: Gauss-féle függvény interpolációja eredeti függvény minták
f ( x) =
x2 e 2
66
Optimális eltolás
Gain = 10 log10
2 e lineáris 2 e eltolt lineáris
τopt = 0.21
67
Optimális lineáris interpoláció T. Blu et al. (2004) eredeti függvény egyenközű minták lineáris interpoláció eltolt lineáris interpoláció
optimális eltolás:
eltolás 68
A lineáris interpoláció szabályos hibát tartalmaz
69
Képfeldolgozás lineáris interpolációval eredeti kép
hagyományos lineáris interpoláció
eltolt lineáris interpoláció
12x elforgatás és interpoláció után 70
15 x 24°-os elforgatás és interpoláció
SNR = jel-zaj viszony (signal to noise ratio)
lineáris interpoláció SNR=22 dB, CPU=5.1 s
eredeti kép
eltolt lineáris interpoláció SNR=28.25 dB, CPU=5.1 s
Keys köbös interpoláció SNR=28.15 dB, CPU=9.4 s 71
Itt a vége...?
72