Waveletová transformace a její použití při zpracování signálů BÍLOVSKÝ, Petr1 1
Katedra elektrických měření, VŠB-TU Ostrava, 17. listopadu, Ostrava - Poruba, 708 33,
[email protected]
Abstrakt: Wavelet transforms provide a tool for decomposing signals into elementary building blocks, called wavelets, much as the Fourier transform represents signals in terms of elementary periodic waves. The fundamental difference between Fourier and wavelet transforms is that wavelets are nonperiodic: the terms in a wavelet expansion are built out of dilates and shifts of a single "mother wavelet," so the expansion is localized in both frequency and time. This makes wavelet transforms extremely useful for accurate coding and analysis of nonstationary or aperiodic data. Klíčová slova: wavelet, časově-frekvenční analýza, zpracování dat
1
Úvod
V tomto příspěvku bude ukázáno použití waveletovy transformace pro zpracování naměřených elektrických signálů. V první části jsou položeny základy pro diskrétní waveletovou transformaci a rozklad signálů pomocí pyramidálního algoritmu. V druhé části jsou předloženy dva typické příklady použití waveletovy transformace. Nejprve vyhlazování a komprese signálů nulováním nejmenších waveletových koeficientů a dále časově-frekvenční analýza signálu. Pro ukázku je v příloze zobrazeno několik základních tvarů škálových funkcí, waveletů a jejich frekvenční charakteristiky. Protože jde o problematiku velmi rozsáhlou a již značně rozpracovanou, lze pro hlubší studium waveletové transformace doporučit literaturu uvedenou v seznamu.
2
Waveletová transformace WT
Waveletová transformace poskytuje oproti Fourierové transformaci (FT) informaci o časové lokalizaci spektrálních složek. FT není vhodná pro analýzu nestacionárních signálů. FT využívá kosinové a sinové funkce pro rozklad signálů, a je nejlepší pro popis periodických signálů. WT nabízí nový přístup k analýze signálů použitím speciálního filtru nazvaného wavelet (vlnka). Každá waveletová funkce osciluje pouze v okolí bodu lokalizace což poskytuje dobrou prostorovou lokalizaci. Cílem WT je rozložit vstupní signál do řady waveletových koeficientů. Toto je dosaženo filtrováním signálu párem ortogonálních filtrů. Jsou označeny jako otcovský wavelet a mateřský wavelet. Otcovský wavelet určuje celkový trend signálu - rozklad na škálové koeficienty, zatímco mateřský wavelet zachycuje doplňkovou informaci o „jemnostech“ na jednotlivých úrovních - waveletové koeficienty. Základní rozdíl mezi Fourierovou a WT je v tom, že wavelety nejsou periodické funkce: termíny v wavelet roztažení jsou stavět rozšířit a posunutí jednotlivé "mateřského waveletu," tak roztažení je lokalizovat jak ve frekvenci tak v času. Toto dělá WT velmi výhodnou pro analýzu nestacionárních nebo aperiodických signálů. WT lze úspěšné použít v tak různých oborech jako je komprese obrazu, analýza řeči, analýza přechodných dějů nebo odhalování poruch.
Waveletová transformace je okenní operace. Jádro waveletové transformace je získáno posunutím a roztažením vybrané bázové funkce. Wavelety jsou speciální okenní funkce ψ ( t ) , které mají nulovou střední hodnotu
ò
+∞
−∞
ψ ( t ) dt = 0
WT funkce f je definována jako WT (α, β) =
ò
f (t )
R
α β
kde je
æ ( t − β) ö ÷÷ dt ψ çç α è α ø
1
α, β ∈ R , α ≠ 0
dilatační parametr translační parametr
Bližší informace k waveletové transformaci viz. [1][2][3].
3
Diskrétní waveletová transformace DWT
Často se objevuje požadavek, aby studovaný signál byl plně rekonstruovatelný z nějaké vhodné diskrétní sítě (vhodné je N = 2 n , n ∈ Z ). Příslušný diskrétní systém waveletů musí -
z důvodu rekonstruovatelnosti - tvořit ortonormální bázi v L2 ( R) .
Jedná se o výpočet škálových a waveletových koeficientů. Škálové koeficienty určují celkový trend vstupního signálu a waveletové koeficienty obsahují doplňkovou informaci. Při volbě vhodného waveletu by se mělo přihlédnout k tomu, že wavelet by měl sledovat průběh zpracovávaného signálu. V důsledku ortogonality koeficienty pro výpočet doplňkové informace při daných filtračních koeficientech budou q L −1−n= ( − 1) hn n
åq ∀n
n
=0
Konkrétně pro L = 4 Daubechieové jsou tyto filtrační koeficienty (L - počet filtračních koeficientů):
h0 =
1+ 3 3+ 3 3− 3 1− 3 , h1 = , h2 = , h3 = , (ostatní hn = 0). 4 2 4 2 4 2 4 2
q 0 = h3 , q1 = -h2 , q 2 = h1 , q 3 = -h0
Obr 1. Výsledek waveletovy analýzy Pyramidální algoritmus
Algoritmus výpočtu se realizuje jako postupný přechod z vyšší hladiny j -1 na nižší j. V první hladině multirozkladu j = 1 transformujeme signál do dvou částí, a to hrubší aproximace a detailní složky, jejíž koeficienty představují nízko- a vysokofrekvenční informaci o signálu. V dalším kroku multirozkladu provádíme stejnou transformaci aproximace signálu získané v předchozí fázi. Detailní složky na daných hladinách již
nerozkládáme. Výpočet je proveden na konečném počtu hladin j = 1,2,3, sekvencí na úrovni j = 0. Jedná se o rozklad signálu na škálové koeficienty y k ( j ) = å y n ( j −1 ) hn − 2 k , n = 0, 1, K , N − 1,
K ,J se vstupní
k = 0, 1, K ,
N −1 2
k = 0, 1, K ,
N −1 2
n
a waveletové koeficienty d k ( j ) = å y n ( j −1) qn − 2 k , n = 0, 1, K , N − 1, n
Mějme vstupní signál y s šířkou pásma π . Tento signál je nejprve filtrován dolnopropustným filtrem h a hornopropustným filtrem q. Tyto signály jsou poté podvzorkovány dvěma protože zabírají jen polovinu původního pásma. y[n]
q[n]
f = 0~π
h[n]
f = π/2~π
f = 0~π/2
2 Úroveň 1 DWT koeficientů
2
q[n]
h[n]
f = π/4~π/2
f = 0~π/4
2 Úroveň 2 DWT koeficientů
2
q[n]
h[n]
f = π/8~π/4
f = 0~π/8
2
2
Úroveň 3 DWT koeficientů
. . .
Obr 2: Pyramidální algoritmus DWT rozkládá vstupní signál do malého počtu „velkých“ dat uložených v koeficientech (J) a velkého počtu „malých“ dat, kterými jsou koeficienty d k ( j ) , j = 1,2, K J. Tak se yk otevírá prostor pro různé formy jejich zpracování, jako jsou například komprese, vyhlazování a filtrace. Paketový rozklad
Paketový rozklad nám umožňuje dále analyzovat i detailní složky a zpřesňovat tak frekvenční lokalizaci výsledných koeficientů. Ve svém důsledku je možné rozkládat libovolnou složku ( aproximaci, detail ) na libovolné hladině.
4
Příklady
Příklad 1: Vyhlazení naměřeného signálu pomocí waveletovy transformace; komprese dat.
Vyhlazování využívající DWT spočívá v provedení tří kroků: 1. pyramidálním algoritmem provedeme rozklad vstupního signálu y ( 0) na škálové y ( j ) a waveletové koeficienty d ( j ) , 2. vynulujeme zadané procento nejmenších waveletových koeficientů d ( j ) v absolutní hodnotě,
3. rekonstrukce signálu rekonstrukčním pyramidálním algoritmem (zpětné DWT).
Obr. 3a - Naměřený signál y ( 0) Obr. 3b - Waveletové koeficienty po rozkladu signálu do dvou bodů ( J = 9). Data jsou zobrazená postupně za sebou v tomto pořadí : y ( 9 ) , d ( 9 ) , d ( 8) , K , d ( 2 ) , d (1) (2 + 2 + 4 +K+256 + 512 = 1024 hodnot ) Obr. 3c - Waveletové koeficienty po vynulování 90% nejmenších hodnot Obr. 3d - Rekonstruovaný signál Na obrázku 3d. je vidět rekonstruovaný vyhlazený signál. Pro rozklad byl použit Daubechiesové wavelet L = 4 - viz příloha. Komprese dat je klasickou aplikací waveletů. Prakticky se jedná o speciální použití vyhlazovacího algoritmu. Výsledkem komprese jsou data získána z jeho druhého kroku. Protože waveletové koeficienty obsahují informace o „jemnostech“, vyskytují se po provedení práhování ve vektorech d ( j ) dlouhé úseky nul. Tohoto faktu lze využít k úspoře paměťového místa. Dekompresi pak představuje třetí krok vyhlazovacího algoritmu. Příklad 2: Časově frekvenční analýza
Pro tento příklad byl vygenerován signál u kterého se postupně mění frekvence (chirp). Jeho časový záznam je zobrazen na prvním obrázku. Tento signál byl zpracován pomocí Fourierovy transformace a pomocí waveletovy transformace - paketový rozklad (best basis) Pro analýzu byl zvolen wavelet Coiflet24 - viz příloha. P2a: Kubický chirp t ∈ <0,1> y = y + sin(2*pi*t*(152*t^3));
Obr. 4. Kubický chirp P2b: Lineární kombinace - lineární a kubický chirp t ∈ <0,1> y = sin(2*pi*t*(500*t)) / 3; //lineární y = y + sin(2*pi*t*(152*t^3)); //kubický
Obr. 5. Lineární a kubický chirp Na obrázku je znázorněn vygenerovaný signál, amplitudové spektrum získané pomocí rychlé Fourierovy transformace a časově frekvenční reprezentace tohoto signálu pomocí waveletovy transformace viz [5]. Fourierova transformace nám říká které frekvence se v signálu vyskytují, ale nelze zjistit v kterém čase se vyskytují. Zde je velmi vhodná waveletová transformace.
5
Použitá literatura
[1] DAUBECHIES, I. The Wavelet Transform: Time-Frequency Localization and Signal Analysis, IEEE Trans. Inform. Theory, Vol IT-36, No.5 961-1005, 1990 [2] KAISER, G. A Friendly Guide to Wavelets Birkhauser-Boston, 1994 [3] ADHEMAR, B. Wavelets with applications in signal and image processing, 1999 [4] CHUI, C.K, Wavelets: A Matematical Tool for Signal Analysis, Philadelphia, 1997 [5] OJANEN, H. WAVEKIT: a Wavelet Toolbox for Matlab, 1998
[6] ČASTOVÁ, N. Časově frekvenční analýza, sylabus pro doktorandské studium, VŠB-TU Ostrava, 1996
Příloha A: Některé známé wavelety, jejich škálové funkce a odpovídající frekvenční charakteristiky: