小波包分解python代碼（小波包分解算法）

• 1、python如何實現類似matlab的小波濾波？
• 2、急求用matlab編小波包3級分解及輸出其係數的程序，那位高手幫助我一下，將不勝感激，先謝謝啦！
• 3、最佳小波包分解層數的獲取
• 4、怎麼用matlab小波包分解
• 5、小波包分解
• 6、matlab中獲得小波樹的末層結點

T=wpdec(y,5,’db40′);

%信號y進行波包解層數5T波樹plot看

a10=wprcoef(T,[1,0]);

%a10節點[1,0]進行重構信號貌似沒層重構說吧能某層某節點進行重構節點編號波樹

%以下為濾波程序（主要調節參數c的大小）

c=10;

wn=0.1;

fs=50000; %採樣頻率；

b=fir1(c,wn/(fs/2),hamming(c+1));

y1=filtfilt(b,1,y);%對y濾波。

load ir; %將信號裝入MATLAB工作環境

s=X106_BA_time(1:1000); %取採樣信號的前1~1000個採樣點

is=length(s); %計算採樣序列長度

subplot(321);plot(s);title(‘原始信號’);%畫出原始信號波形

xlabel(‘樣本序列號’);

ylabel(‘幅值A’);

wpt=wpdec(s,3,’db1′,’shannon’);% 用db1小波包對信號x3層分解，用shannon熵作為熵標準

plot(wpt);%繪製小波包樹

N=allnodes(wpt);% 計算小波包分解樹的結點

%提取各節點的小波包係數

for i=1:length(N)

X=wpcoef(wpt,i-1);

subplot(floor((length(N)+1)/2),2,i);plot(X);

title([‘節點’,num2str(i) ‘的小波包係數’]);

end

figure;

%各節點小波包重構係數

for i=1:length(N)

rcfs=wprcoef(wpt,i-1);

subplot(floor((length(N)+1)/2),2,i);plot(rcfs);

title([‘重構節點’,num2str(i) ‘小波包係數’]);

end

因為小波包分解高頻係數的奇異值在每一分解層都不同，即分解層數的變化會引起奇異值的變化，因此，可以通過分析奇異值變化規律來確定分解層數。定義一個判定函數C（r，k），用其函數值的變化來表徵小波包分解每層信息的變化，設C（r，k）函數為

高光譜遙感影像信息提取技術

式中：r為特徵矩陣維數；k為正整數。

通過大量的實驗分析，表明小波包的分解層數與各層數對應C（r，k）函數值的變化存在一定的內在聯繫。隨着各層C（r，k）函數值的變化，有用信號和噪聲的比值也發生着顯著的變化。當在某分解層處C（r，k）的函數值到達某一程度時，這就使目標地物更好地分離出來，而使噪聲減弱，這樣就可以確定最佳分解層數。

在AVIRIS實驗數據1～9層的小波包分解中，根據高頻係數矩陣分別計算出各層處C（r，k）函數值。隨着分解層數的增加，其C（r，k）函數值在不斷減小，其變化率也發生着顯著變化。圖4.6顯示了各分解層數下的C（r，k）函數值的實際變化趨勢和這種趨勢變化所反演的回歸函數曲線，其回歸曲線方程採用接近C（r，k）函數值變化趨勢的最低次多項式函數來表示。從圖4.6上可以看出，兩條曲線非常接近地重合在一起，說明了回歸函數曲線能很好地表徵C（r，k）的變化特徵。

根據各分解層下離散的C（r，k）函數值，找出最接近其變化趨勢的最低次多項式回歸函數方程y（x）：

高光譜遙感影像信息提取技術

然後求解方程y（x）的二階導數y″（x），即找出曲線拐點的位置，此點位置就代表着小波包分解的最佳分解層數。在表4.1中，可以看到在分解層數達到5層時，y″（x）的值開始從負值變為正值，說明y″（x）在第4層和第5層之間有拐點，從而確認第四層為最佳分解層數。

圖4.6 基於影像AVIRIS的判定函數C（r，k）值變化曲線和回歸函數曲線

表4.1 高光譜影像AVIRIS不同小波包分解層的SVD特徵值分析

小波圖像去噪的方法大概分為3類

1:基於小波變換摸極大值原理

2：基於小波變換係數的相關性

3：基於小波閾值的去噪。

基於小波閾值的去噪方法3個步驟：

1：

計算含噪聲圖像的小波變換。選擇合適的小波基和小波分解層數j，運用matlab

分解算法將含有噪聲圖像進行j層小波分解，得到相應的小波分解係數。

2：對分解後的高頻係數進行閾值量化，對於從1

到j的每一層，選擇一個適當的閾值和合適的閾值函數，將分解得到的高頻係數進行閾值量化，得到估計小波係數。

3:進行小波逆變化，根據圖像小波分解後的第j層，低頻

係數（尺度係數）和經過閾值量化處理的各層高頻係數（小波係數），運用matlab重構算法進行小波重構，得到去噪後的圖像。

1.裝載信號

在MATLAB命令行中輸入

load noisbloc

s=noisbloc(1:1024);

ls=length(s);

plot(s);

2.完成信號的單尺度一維離散小波分解

採用db4基本小波分解信號

在命令窗口中輸入

[cA1,cD1]=dwt(s,’db4′);

這就產生了低頻係數cA1和高頻係數cD1。

可以通過whos命令查看cA1和cD1的長度：

whos

Name Size Bytes Class

cA1 1×515 4120 double array

cD1 1×515 4120 double array

ls 1×1 8 double array

noisbloc 1×1024 8192 double array

s 1×1024 8192 double array

Grand total is 3079 elements using 24632 bytes

3.從係數中重構低頻部分和高頻部分

從第二步產生的係數cA1和cD1構造第一層的低頻和高頻（A1和D1）係數；

A1=upcoef(‘a’,cA1,’db4′,1,ls);

D1=upcoef(‘a’,cD1,’db4′,1,ls);

或用下面兩個函數：

A1=idwt(cA1,[],’db4′,ls);

D1=idwt(cD1,[],’db4′,ls);

4.顯示高頻和低頻部分

為了顯示第一層分解結果，輸入

subplot(211); plot(A1); title(‘低頻A1’);

subplot(212); plot(D1); title(‘低頻D1’);

5.由小波逆變換恢複信號

使用idwt函數很容易實現，在命令窗口輸入

subplot(211); plot(s); title(‘原始信號’);

subplot(212); plot(A0); title(‘重構信號’);

6.多層一維分解

為了完成一個5層的分解，輸入：

[C,L]=wavedec(s,5,’db4′);

7.提取係數的低頻和高頻部分

為了從上面的C中提取第3層的低頻係數，輸入：

cA3=appcoef(C,L,’db4′,3);

8.重構第3層的低頻係數

為了從上面的C中重構第3層的低頻係數，輸入：

A3=wrcoef(‘a’,C,L,’db4′,3);

9.重構第1、2、3、4、5層的高頻信號

其方法是：

cA5=appcoef(C,L,’db4′,5);

A5=wrcoef(‘a’,C,L,’db4′,5);

D1=wrcoef(‘d’,C,L,’db4′,1);

D2=wrcoef(‘d’,C,L,’db4′,2);

D3=wrcoef(‘d’,C,L,’db4′,3);

D4=wrcoef(‘d’,C,L,’db4′,4);

D5=wrcoef(‘d’,C,L,’db4′,5);

下面顯示多尺度一維分解的結果：

subplot(322); plot(D1); title(‘低頻D1’);

subplot(323); plot(D2); title(‘低頻D2’);

subplot(323); plot(D3); title(‘低頻D3’);

subplot(325); plot(A5); title(‘低頻D4’);

subplot(326); plot(A5); title(‘低頻D5’);

10.重構原始信號並顯示

A0=waverec(C,L,’db4′);

subplot(311); plot(s); title(‘原始信號’);

subplot(312); plot(A0); title(‘重構信號’);

subplot(313); plot(s-A0); title(‘誤差信號’);

第一行：將wave 用 meyr小波進行3層小波包分解，獲得一個小波包樹 t

第二行：將小波包樹的第二行的四個節點收起來，也就是讓第二行的節點變為樹的最底層節點。因為第一行中小 波包樹的節點個數是 第一層2個，第二層4個，第三層8個。現在將t2就是將第三層的節點再聚合回第二 層。

第三行：讀取第二層四個節點係數的size

第四~七行：將所有四個節點的小波包係數變為0

第八行：將四個節點的係數重組到t3小波樹中。

第九行：對t3小波樹進行重構，獲得信號wave2