離散小波變換

1、1.離散小波變換長(zhǎng)期以來(lái)，離散小波變換(Discrete Wavelet Transform)在數(shù)字信號(hào)處理、石油勘探、地 震預(yù)報(bào)、醫(yī)學(xué)斷層診斷、編碼理論、量子物理及概率論等領(lǐng)域中都得到了廣泛的應(yīng)用。各種 快速傅氏變換(FFT)和離散小波變換(DWT)算法不斷出現(xiàn)，成為數(shù)值代數(shù)方面最活躍的一個(gè) 研究領(lǐng)域，而其意義遠(yuǎn)遠(yuǎn)超過了算法研究的范圍，進(jìn)而為諸多科技領(lǐng)域的研究打開了一個(gè)嶄新的局面。本章分別對(duì)FFT和DWT的基本算法作了簡(jiǎn)單介紹，若需在此方面做進(jìn)一步研究, 可參考文獻(xiàn)2。1.1離散小波變換DWT離散小波變換DWT及其串行算法先對(duì)一維小波變換作一簡(jiǎn)單介紹。設(shè)f(x)為一維輸入信號(hào)，記jk(x)

2、2 j/2 (2 jx k)，jk(x) 2 j/2 (2 jx k)，這里(x)與(x)分別稱為定標(biāo)函數(shù)與子波函數(shù)， jk(x)與 jk(x)為二個(gè)正交基函數(shù)的集合。記Pof=f,在第j級(jí)上的一維離散小波變換DWT(DiscreteWavelet Transform)通過正交投影 Pjf與Qf將Pj-1f分解為:Pj 1 fPjf Qjfckkjkdk jkkjp 1j 1其中：冷/，dP 1j 1n 嚴(yán))c2kn(j1,2，L,k0,1,N 2j 1),這里,h(n)與g( n)分別為低通與高通權(quán)系數(shù)，它們由基函數(shù) jk(x)與 jk(x)來(lái)確定，p為權(quán)系數(shù)的長(zhǎng)度。C：為信號(hào)的輸入數(shù)據(jù)，N

3、為輸入信號(hào)的長(zhǎng)度，L為所需的級(jí)數(shù)。由上式可見，每級(jí)一維DWT與一維卷積計(jì)算很相似。 所不同的是：在DWT中，輸出數(shù)據(jù)下標(biāo)增加1時(shí), 權(quán)系數(shù)在輸入數(shù)據(jù)的對(duì)應(yīng)點(diǎn)下標(biāo)增加2,這稱為“間隔取樣”。算法22.3 一維離散小波變換串行算法輸入：c0=d 0(C0°, C1°,CN-10)h=(h 0, h1,，hL-1)g=(go, g1,，gL-1)輸出：cij , dij (i=0, 1,，N/2j-1, j > 0)Begi n(1) j=0, n=N(2) While (n> 1) do(2.1) for i=0 to n-1 do(2.1.1) c+1 =0, d

4、 ij+1 =0(2.1.2) for k=0 to L-1 docj1cj1dij1dij1gkd(jk2i)mod nend for end for(2.2) j=j+1, n=n/2end whileEnd顯然，算法22.3的時(shí)間復(fù)雜度為O(N*L)。在實(shí)際應(yīng)用中，很多情況下采用緊支集小波(Compactly Supported Wavelets),這時(shí)相應(yīng)的尺度系數(shù)和小波系數(shù)都是有限長(zhǎng)度的，不失一般性設(shè)尺度系數(shù)只有有限個(gè)非零值：h1,hN, N為偶數(shù)，同樣取小波使其只有有限個(gè)非零值：gw ,gN。為簡(jiǎn)單起見，設(shè)尺度系數(shù)與小波函數(shù)都是實(shí)數(shù)。對(duì)有限長(zhǎng)度的輸入數(shù)據(jù)序列：c0 Xn,n 1,2

5、, ,M (其余點(diǎn)的值都看成0)，它的離散小波變換為dj k1Cn gnn Z2kj 0,1,J 1其中J為實(shí)際中要求分解的步數(shù)，最多不超過log2M，其逆變換為cj 1Ck hn 2k k ZCk hn 2k k Zj J,1注意到尺度系數(shù)和輸入系列都是有限長(zhǎng)度的序列,上述和實(shí)際上都只有有限項(xiàng)。若完全按照上述公式計(jì)算，在經(jīng)過J步分解后，所得到的J+1 個(gè)序列 dkj, j 0,1, ,J 1 和 ck 的非零項(xiàng)的個(gè)數(shù)之和一般要大于 過程。j=0時(shí)計(jì)算過程為M,究竟這個(gè)項(xiàng)目增加到了多少？下面來(lái)分析一下上述計(jì)算1 Mckxnhn 2kn 1Md kxngn 2kn 1的非零值范圍為：kN 1,

6、, 1,0, , 1,即有2 2M2現(xiàn)的規(guī)律相似。1個(gè)非零值。d：的非零值范圍相同。M N 12分解多步后非零項(xiàng)的個(gè)數(shù)可能比輸入序列的長(zhǎng)度增加較多。繼續(xù)往下分解時(shí)，非零項(xiàng)出例如，若輸入序列長(zhǎng)度為100, N=4，則d：有51項(xiàng)非零，d2有27項(xiàng)非零，d3有15項(xiàng)非零，d：有9項(xiàng)非零，dk有6項(xiàng)非零，d6有4項(xiàng)非零，C?有4項(xiàng)非零。這樣分解到 6步后得到的序列的非 零項(xiàng)個(gè)數(shù)的總和為116，超過了輸入序列的長(zhǎng)度。在數(shù)據(jù)壓縮等應(yīng)用中，希望總的長(zhǎng)度基本 不增加，這樣可以提高壓縮比、減少存儲(chǔ)量并減少實(shí)現(xiàn)的難度。可以采用稍微改變計(jì)算公式的方法，使輸出序列的非零項(xiàng)總和基本上和輸入序列的非零項(xiàng)數(shù)相等，并且可以

7、完全重構(gòu)。這種方法也相當(dāng)于把輸入序列進(jìn)行延長(zhǎng)（增加非零項(xiàng)），因而稱為延拓法。只需考慮一步分解的情形，下面考慮第一步分解（j=1）。將輸入序列作延拓，若M為偶數(shù)，直接將其按 M為周期延拓；若 M為奇數(shù)，首先令xM 1 0。然后按M + 1為周期延拓。 作了這種延拓后再按前述公式計(jì)算，相應(yīng)的變換矩陣已不再是H和G,事實(shí)上這時(shí)的變換C1HCD1GC矩陣類似于循環(huán)矩陣。例如，當(dāng)M=8，N=4時(shí)矩陣H變?yōu)椋篽3h40000昨h2hh2h3h4000000h1h2h3h4000000hh2h3h4h3h40000昨h2當(dāng)M=7, N=4時(shí)矩陣H變?yōu)椋篽3h40000hhh2h3h400000h1h2h3h

8、400000h1h2h3h3h40000h從上述的矩陣表示可以看出，兩種情況下的矩陣內(nèi)都有完全相同的行，這說(shuō)明作了重復(fù)計(jì)算，因而從矩陣中去掉重復(fù)的那一行不會(huì)減少任何信息量，也就是說(shuō)，這時(shí)我們可以對(duì)矩陣進(jìn)行截短（即去掉一行），使得所得計(jì)算結(jié)果仍然可以完全恢復(fù)原輸入信號(hào)。當(dāng)M=8，N=4時(shí)截短后的矩陣為：h400000h2Hhh2h3h4000000h1h2館h4000000hh2h3h4當(dāng)M=7, N=4時(shí)截短后的矩陣為：h3h40000hHhh2h3h4000000h2h3h400000h1h2h3這時(shí)的矩陣都只有節(jié)行。分解過程成為:向量C1和D1都只有 個(gè)元素。重構(gòu)過程為:2C 0 H *

9、C1 G* D1可以完全重構(gòu)。矩陣H , G有等式H*H + G*G = I一般情況下，按上述方式保留矩陣的M 行，可以完全恢復(fù)原信號(hào)。2這種方法的優(yōu)點(diǎn)是最后的序列的非0元素的個(gè)數(shù)基本上和輸入序列的非0元素個(gè)數(shù)相同，特別是若輸入序列長(zhǎng)度為2的幕，則完全相同，而且可以完全重構(gòu)輸入信號(hào)。其代價(jià)是得到的變換系數(shù)Dj中的一些元素已不再是輸入序列的離散小波變換系數(shù)，對(duì)某些應(yīng)用可能是不適合的，但在數(shù)據(jù)壓縮等應(yīng)用領(lǐng)域，這種方法是可行的。Begi n對(duì)所有處理器 my_rank(my_rank=0,p-1)同時(shí)執(zhí)行如下的算法：(1) j=0；(2) while (j<r) do(2.1) 將數(shù)據(jù)c,(n 0,1,N/2j 1)按塊分配給p臺(tái)處理器(2.2) 將處理器i+1中前L-1個(gè)數(shù)據(jù)發(fā)送給處理器i(2.3) 處理器 i 負(fù)責(zé) C* 1 和 d/t n i Nj-y,(i1) N-71)的計(jì)算p2j 1p2 j 1(2.4) j=j+1 end whileEnd這里每一步分解后數(shù)據(jù) 材1和dnj 1已經(jīng)是按塊存儲(chǔ)在P臺(tái)處理器上，因此算法第一步中的數(shù)據(jù)分配除了j=0時(shí)需要