代码：小波包分解与重构、小波包能量特征提取

1、小波变换的理解

傅里叶变换——短时傅里叶变换——小波变换。

参考文献：以下两篇参考资料讲述得十分清楚，有助于理解小波变换。

但具体的数学角度阐述，请参考其他资料。

（1）知乎专栏：形象易懂讲解算法I——小波变换

https://zhuanlan.zhihu.com/p/22450818

（2）知乎专栏：傅里叶分析之掐死教程。

https://zhuanlan.zhihu.com/p/19763358

2、小波包分解

小波包是为了克服小波分解在高频段的频率分辨率较差，而在低频段的时间分辨率较差的问题的基础上而提出的。

它是一种更精细的信号分析的方法，提高了信号的时域分辨率。

下面是两者的对比图：

3、能量谱

基于小波包分解提取多尺度空间能量特征的原理是把不同分解尺度上的信号能量求解出来，将这些能量值按尺度顺序排列成特征向量供识别使用。

4、Matlab代码

给出两部分代码，写成两个函数。一个是小波包分解与重构，另一个是能量谱函数。

代码名称：wavelet_packetdecomposition_reconstruct

1. function wpt= wavelet_packetdecomposition_reconstruct( x,n,wpname )  
2. %% 对信号进行小波包分解，得到节点的小波包系数。然后对每个节点系数进行重构。   
3. % Decompose x at depth n with wpname wavelet packets.using Shannon entropy.  
4. %     
5. %  x-input signal   
6. %  n-the number of decomposition layers  
7. %  wpname-a particular wavelet.type:string.  
8. %  
9. %Author hubery_zhang  
10. %Date 20170714  
11.   
12. %%  
13. wpt=wpdec(x,n,wpname);  
14. % Plot wavelet packet tree (binary tree)  
15. plot(wpt)  
16. %% wavelet packet coefficients.default：use the front 4.  
17. cfs0=wpcoef(wpt,[n 0]);  
18. cfs1=wpcoef(wpt,[n 1]);  
19. cfs2=wpcoef(wpt,[n 2]);  
20. cfs3=wpcoef(wpt,[n 3]);  
21. figure;  
22. subplot(5,1,1);  
23. plot(x);  
24. title('原始信号');  
25. subplot(5,1,2);  
26. plot(cfs0);  
27. title(['结点 ',num2str(n) '  1',' 系数'])  
28. subplot(5,1,3);  
29. plot(cfs1);  
30. title(['结点 ',num2str(n) '  2',' 系数'])  
31. subplot(5,1,4);  
32. plot(cfs2);  
33. title(['结点 ',num2str(n) '  3',' 系数'])  
34. subplot(5,1,5);  
35. plot(cfs3);  
36. title(['结点 ',num2str(n) '  4',' 系数'])  
37. %% reconstruct wavelet packet coefficients.  
38. rex0=wprcoef(wpt,[n 0]);  
39. rex1=wprcoef(wpt,[n 1]);  
40. rex2=wprcoef(wpt,[n 2]);  
41. rex3=wprcoef(wpt,[n 3]);  
42. figure;  
43. subplot(5,1,1);  
44. plot(x);  
45. title('原始信号');  
46. subplot(5,1,2);  
47. plot(rex0);  
48. title(['重构结点 ',num2str(n) '  1',' 系数'])  
49. subplot(5,1,3);  
50. plot(rex1);  
51. title(['重构结点 ',num2str(n) '  2',' 系数'])  
52. subplot(5,1,4);  
53. plot(rex2);  
54. title(['重构结点 ',num2str(n) '  3',' 系数'])  
55. subplot(5,1,5);  
56. plot(rex3);  
57. title(['重构结点 ',num2str(n) '  4',' 系数'])  
58. end  

代码名称：wavelet_energy_spectrum

1. function E = wavelet_energy_spectrum( wpt,n )  
2. %% 计算每一层每一个节点的能量  
3. %  wpt-wavelet packet tree  
4. %  n-第n层能量  
5. %   
6. % Author hubery_zhang  
7. % Date  20170714  
8.   
9. %%  
10. % 求第n层第i个节点的系数  
11. E(1:2^n )=0;  
12. for i=1:2^n   
13. E(i) = norm(wpcoef(wpt,[n,i-1]),2); %求第i个节点的范数平方，其实也就是平方和  
14. end  
15. %求每个节点的概率  
16. E_total=sum(E);   
17. for i=1:2^n  
18. p_node(i)= 100*E(i)/E_total;  
19. end  
20. % E = wenergy(wpt); only get the last layer  
21. figure;  
22. x=1:2^n;  
23. bar(x,p_node);  
24. title(['第',num2str(n),'层']);  
25. axis([0 2^n 0 100]);  
26. xlabel('结点');  
27. ylabel('能量百分比/%');  
28. for j=1:2^n  
29. text(x(j),p_node(i),num2str(p_node(j),'%0.2f'),...  
30.     'HorizontalAlignment','center',...  
31.     'VerticalAlignment','bottom')  
32. end
33. end