小波包分解

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%小波包分解程序%%

m=load('300_30.txt'); 导入文件名为300_30的txt文件N=length(m);

for i=1:N-1 ;

q(i,1)=m(i,1);

end;

d=q';

s1=d;

change=1000;

[c,l] = wavedec(d,3,'db4'); %提取小波分解后的低频系数

ca3=appcoef(c,l,'db4',3); %提取各层小波分解后的高频系数cd3=detcoef(c,l,3);

cd2=detcoef(c,l,2);

cd1=detcoef(c,l,1); %对信号强制消噪

cdd3=zeros(1,length(cd3));%第三层高频系数cd3全置0 cdd2=zeros(1,length(cd2));%第二层高频系数cd2全置0 cdd1=zeros(1,length(cd1));%第一层高频系数cd1全置0

c1=[ca3,cdd3,cdd2,cdd1];%建立新的系数矩阵

s2=waverec(c1,l,'db4')%为新的分解结构

%[thr,sorh,keepapp]=ddencmp('den','wv',d);

%s2=wdencmp('gbl',c,l,'db4',4,thr,sorh,keepapp);

%subplot(413)

%plot(1:change,s2(1:change));

%title('默认软阈值消噪后信号')

figure(1)

subplot(9,2,1)

plot(1:change,s1(1:change))

title('原始信号')

ylabel('S1')

subplot(9,2,2)

plot(1:change,s2(1:change))

title('强制消噪后信号')

ylabel('S2')

wpt=wpdec(s1,3,'db1','shannon');

%plot(wpt);

%重构第三层个节点小波系数

s130=wprcoef(wpt,[3,0]);

s131=wprcoef(wpt,[3,1]);

s132=wprcoef(wpt,[3,2]);

s133=wprcoef(wpt,[3,3]);

s134=wprcoef(wpt,[3,4]);

s135=wprcoef(wpt,[3,5]);

s136=wprcoef(wpt,[3,6]);

s137=wprcoef(wpt,[3,7]);

%计算第三层个节点小波能量

s10=norm(s130);

s11=norm(s131);

s12=norm(s132);

s13=norm(s133);

s14=norm(s134);

s15=norm(s135);

s16=norm(s136);

s17=norm(s137);

%计算方差

st10=std(s130);

st11=std(s131);

st12=std(s132);

st13=std(s133);

st14=std(s134);

st15=std(s135);

st16=std(s136);

st17=std(s137);

disp('正常信号的特征向量');

snorm1=[s10,s11,s12,s13,s14,s15,s16,s17]; std1=[st10,st11,st12,st13,st14,st15,st16,st17]; %显示三层个节点小波系数

subplot(9,2,3);plot(1:change,s130(1:change)); ylabel('S130');

subplot(9,2,5);plot(1:change,s131(1:change)); ylabel('S131');

subplot(9,2,7);plot(1:change,s132(1:change)); ylabel('S132');

subplot(9,2,9);plot(1:change,s133(1:change)); ylabel('S133');

subplot(9,2,11);plot(1:change,s134(1:change)); ylabel('S134');

subplot(9,2,13);plot(1:change,s135(1:change)); ylabel('S135');

subplot(9,2,15);plot(1:change,s136(1:change)); ylabel('S136');

ylabel('S137');

wpt=wpdec(s2,3,'db1','shannon');

%plot(wpt);

s230=wprcoef(wpt,[3,0]);

s231=wprcoef(wpt,[3,1]);

s232=wprcoef(wpt,[3,2]);

s233=wprcoef(wpt,[3,3]);

s234=wprcoef(wpt,[3,4]);

s235=wprcoef(wpt,[3,5]);

s236=wprcoef(wpt,[3,6]);

s237=wprcoef(wpt,[3,7]);

s20=norm(s230);

s21=norm(s231);

s22=norm(s232);

s23=norm(s233);

s24=norm(s234);

s25=norm(s235);

s26=norm(s236);

s27=norm(s237);

st20=std(s230);

st21=std(s231);

st22=std(s232);

st23=std(s233);

st24=std(s234);

st25=std(s235);

st26=std(s236);

st27=std(s237);

disp('故障信号的特征向量');

snorm2=[s20,s21,s22,s23,s24,s25,s26,s27]; std2=[st20,st21,st22,st23,st24,st25,st26,st27]; subplot(9,2,4);plot(1:change,s230(1:change)); ylabel('S230');

subplot(9,2,6);plot(1:change,s231(1:change)); ylabel('S231');

subplot(9,2,8);plot(1:change,s232(1:change)); ylabel('S232');

ylabel('S233');

subplot(9,2,12);plot(1:change,s234(1:change)); ylabel('S234');

subplot(9,2,14);plot(1:change,s235(1:change)); ylabel('S235');

subplot(9,2,16);plot(1:change,s236(1:change)); ylabel('S236');

subplot(9,2,18);plot(1:change,s237(1:change)); ylabel('S237');

%fft

N=1024;

figure(2)

y1=fft(s1,N);

py1=abs(y1)

%py1=y1.*conj(y1)/N;

y2=fft(s2,N);

py2=abs(y2)

%py2=y2.*conj(y2)/N;

y130=fft(s130,N);

py130=y130.*conj(y130)/N;

y131=fft(s131,N);

py131=y131.*conj(y131)/N;

y132=fft(s132,N);

py132=y132.*conj(y132)/N;

y133=fft(s133,N);

py133=y133.*conj(y133)/N;

y134=fft(s134,N);

py134=y134.*conj(y134)/N;

y135=fft(s135,N);

py135=y135.*conj(y135)/N;

y136=fft(s136,N);

py136=y136.*conj(y136)/N;

y137=fft(s137,N);

py137=y137.*conj(y137)/N;

y230=fft(s230,N);

py230=y230.*conj(y230)/N;

y231=fft(s231,N);

py231=y231.*conj(y231)/N;

y232=fft(s232,N);

py232=y232.*conj(y232)/N;

y233=fft(s233,N);

py233=y233.*conj(y233)/N; y234=fft(s234,N);

py234=y234.*conj(y234)/N; y235=fft(s235,N);

py235=y235.*conj(y235)/N; y236=fft(s236,N);

py236=y236.*conj(y236)/N; y237=fft(s237,N);

py237=y237.*conj(y237)/N; para=512;

f=1000*(0:para-1))/N; subplot(1,2,1);

plot(f,py1(1:p ara));

ylabel('P1');

title('原始信号的功率谱') subplot(1,2,2);

plot(f,py2(1:p ara));

ylabel('P2');

title('故障信号的功率谱') figure

subplot(4,2,1);

plot(f,py130(1:p ara)); ylabel('P130');

title('S130的功率谱') subplot(4,2,2);

plot(f,py131(1:p ara)); ylabel('P131');

title('S131的功率谱') subplot(4,2,3);

plot(f,py132(1:p ara)); ylabel('P132');

subplot(4,2,4);

plot(f,py133(1:p ara)); ylabel('P133');

subplot(4,2,5);

plot(f,py134(1:p ara)); ylabel('P134');

subplot(4,2,6);

plot(f,py135(1:p ara)); ylabel('P135');

subplot(4,2,7);

plot(f,py136(1:p ara));

ylabel('P136');

subplot(4,2,8);

plot(f,py137(1:p ara));

ylabel('P137');

figure

subplot(4,2,1);

plot(f,py230(1:p ara));

ylabel('P230');

title('S230的功率谱')

subplot(4,2,2);

plot(f,py231(1:p ara));

ylabel('P231');

title('S231的功率谱')

subplot(4,2,3);

plot(f,py232(1:p ara));

ylabel('P232');

subplot(4,2,4);

plot(f,py233(1:p ara));

ylabel('P233');

subplot(4,2,5);

plot(f,py234(1:p ara));

ylabel('P234');

subplot(4,2,6);

plot(f,py235(1:p ara));

ylabel('P235');

subplot(4,2,7);

plot(f,py236(1:p ara));

ylabel('P236');

subplot(4,2,8);

plot(f,py237(1:p ara));

ylabel('P237');

figure

E=s10+s11+s12+s13+s14+s15+s16+s17

E0=s10/E;E1=s11/E;E2=s12/E;E3=s13/E;E4=s14/E;E5=s15/E;E6=s16/E;E7 =sE;

EE=[E0,E1,E2,E3,E4,E5,E6,E7]

subplot(211)

bar(EE)

M=s20+s21+s22+s23+s24+s25+s26+s27

M0=s20/M;M1=s21/M;M2=s22/M;M3=s23/M;M4=s24/M;M5=s25/M;M6=s26

/M;M7=s27/M;

MM=[M0,M1,M2,M3,M4,M5,M6,M7] subplot(212)

bar(MM)

小波包能量谱程序

wpt4=wpdec(y4,n,'db30'); %对数据进行小波包分解 for i=1:2^n %wpcoef(wpt4,[n,i-1])是求第n层第i个节点的系数 disp('每个节点的能量E1(i)'); E4(i)=norm(wpcoef(wpt4,[n,i-1]),2)*norm(wpcoef(wpt4,[n,i-1]),2)%求第i个节点的范数平方，其实也就是平方和 end 请教各位，小波包能量如何求？ 我的理解 假设信号x，对齐进行n层分解： wpt=wpdec(x,n,wname); 然后各小波包系数重构分量信号： dp(i,: )=(wprcoef(wpt,i)); 小波包能量为： Edp(i)=sum(dp(i,: ).^2); 这样对吗，谢谢大虾指点！ 1.小波分析中，原始信号被分解为逼近部分和细节部分。逼近部分再分解为另一层的逼近和细节，这样的过程重复进行，直到设定的分阶层。然而，在小波包分解中，细节部分也进行相同的分解。小波包分解具有任意多尺度特点，避免了小波变换固定时频分解的缺陷（如高频段频率分辨率低），为时频分析提供了极大的选择余地，更能反映信号的本质和特征。你理解也算是对的。 2. s%为已知信号源 for i=1:4 wpt=wpdec(s,i,'db3'); e=wenergy(wpt); E=zeros(1,length(e)); for j=1:2^i E(j)=sum(abs(wprcoef(wpt,[i,j-1])).^2); end figure(5) subplot(4,1,i); bar(e); axis([0 length(e) 0 130]); title(['第',num2str(i), ' 层']); for j=1:length(e) text(j-0.2,e(j)+20,num2str(e(j),'%2.2f')); end end 这段程序也是从网上下载的，一起学习一下吧。

05 基于小波包能量谱分析的电机故障诊断要点

应用天地 2008年 2月第 27卷第 2期 基于小波包能量谱分析的电机故障诊断 唐友怀张海涛罗珊姜喆 (工程兵工程学院南京 210007 摘要 :小波包是继小波分析后提出的一种新型的多尺度分析方法 [1], 解决了小波分析在高频部分分辨率差的缺点 , 体现了比小波分析更好的处理效果。文章描述了小波包分析方法的基本原理及其实现算法 , 并从能量分布的角度出发 , 阐述了在电机故障诊断中 , 利用小波包分析方法将模糊故障信号进行量化、分解 , 从而便于用单片机进行处理、判断的一种新的应用途径 , 在实验室中模拟各种电机故障进行了实验验证 , 实验进一步表明基于小波包能量谱分析的电机故障诊断方法是一种方便灵活并且准确度很高的故障诊断方法。关键词 :小波包 ; 故障诊断 ; 能量谱 ; 电机中图分类号 :TP182文献标识码 :A B ased on w avelet p acked energy motor fau lt diagnosis Haitao L uo Shan Jiang Zhe (College of Engineering Corps , Nanjing 210007 Abstract :The wavelet packed is presented as a new kind of multiscale analysis technique followed Wavelet analysis. it re 2solved t he wavelet analysis disadvantage on t he part of high frequency resolution lower , showed better treat ment effect t han wavelet analysis. The f undamental and it s realization arit hmetic of t he wavelet packed analysis met hod are described in t his paper. A new application approach of t he wavelet packed met hod on t he motor fault diagnosis from energy distrib 2uting angle is expatiated. And given t he experimental met hod and t he conclusion. and a new application approach which is convenient for t he microchip to process and judge by using t he wavelet packed analysis met hod to make the f uzzy motor fault diagnosis signals quantized and analyzed

小波分析及小波包分析

小波分析及小波包分析 在利用matlab做小波分析时，小波分解函数和系数提取函数的结果都是分解系数。我们知道，复杂的周期信号可以分解为一组正弦函数之和，及傅里叶级数，而傅里叶变换对应于傅里叶级数的系数；同样，信号也可以表示为一组小波基函数之和，小波变换系数对应于这组小波基函数的系数。 多尺度分解是按照多分辨分析理论，分解尺度越大，分解系数的长度越小（是上一个尺度的二分之一）。我们会发现分解得到的小波低频系数的变化规律和原始信号相似，但要注意低频系数的数值和长度与原始信号以及后面重构得到的各层信号是不一样的。 小波分解：具体实现过程可以分别设计高通滤波器和低通滤波器，得到高频系数和低频系数，并且每分解一次数据的长度减半。小波重构，为分分解的逆过程，先进行增采样，及在每两个数之间插入一个0，与共轭滤波器卷积，最后对卷积结果求和。在应用程中，我们经常利用各层系数对信号进行重构（注意虽然系数数少于原信号点数，但是重构后的长度是一样的），从而可以有选择的观看每一频段的时域波形。从而确定冲击成分所在频率范围。便于更直观的理解，小波分解，利用各层系数进行信号重构过程我们可以认为是将信号通过一系列的不同类型的滤波器，从而得到不同频率范围内的信号，及将信号分解。 小波消噪:运用小波分析进行一维信号消噪处理和压缩处理，是小波分析的两个重要的应用。使用小波分析可以将原始信号分解为一系列的近似分量和细节分量，信号的噪声主要集中表现在信号的细节分量上。使用一定的阈值处理细节分量后，再经过小波重构就可以得到平滑的信号。 小波常用函数 [C,L]=wavedec(s,3,'db1');%用小波函数db1对信号s进行3尺度分解 其中C为分解后低频和高频系数，L存储低频和高频系数的长度。 X=wrcoef(‘type’,C,L,’wname’,N)%对一维小波系数进行单支重构，其中N表示对第几层的小波进行重构 X=wrcoef(‘a’,C,L,’wname’,3)%对第三层的低频信号进行重构，如果a变为d的话，表示对低频分量进行重构。注意重构后数据的长度于原来数据的长度一致。 ca1=appcoef(C,L,'db1',1);%从前面小波3尺度分解结构[C,L]中提取尺度1的低频系数 高频系数提取类似。 选择合适的阈值，小波分解后，重构可以达到去除噪声的目的。 小波包分解，可以将信号分在不同的频带，且不同的频带宽度是一样的。小波分析，只将低

太阳黑子数时间序列的奇异谱分析和小波分析

第32卷第6期 2007年11月 测绘科学 Science of Surveying and M app ing Vol 132No 16 Nov 1 作者简介:徐克红(19822),女,山东泰安人,辽宁工程技术大学与中国测绘科学研究院联合培养硕士研究生,主要研究方向为卫星轨道确定。E 2mail:xukehong0719@1631com 收稿日期:2007206228 太阳黑子数时间序列的奇异谱分析和小波分析 徐克红 ①② ,程鹏飞①,文汉江 ① (①中国测绘科学研究院,北京　100039;②辽宁工程技术大学,辽宁阜新　123000) 【摘　要】本文对小波变换和奇异谱分析方法进行了简要介绍,对离散小波的分解和重构、奇异谱分析的重构进 行了详细阐述。结合太阳黑子数1749年至2007年3月期间的月平均值时间序列进行了小波变换的分解和重构及SS A 方法的重构,提取了其主要的周期特性,并对两种分析方法进行了比较。【关键词】小波分析;离散小波的分解与重构;奇异谱分析;太阳黑子数【中图分类号】P228 【文献标识码】A 【文章编号】100922307(2007)0620035204 1　引言 太阳黑子是太阳光球上经常出现的阴暗斑点,是太阳活动的羁绊标志,是反映太阳辐射变化的重要指标,一般用太阳黑子数表示。太阳黑子数反映了太阳活动强弱的变化,对地球的影响很大,诸如地磁变化、大气运动、气候异常、海洋变化等,都和太阳黑子数变化有着不同程度的关系。因此研究太阳黑子数的变化有利于深入了解它对卫星轨道、定位等方面的影响。 对太阳黑子数变化的研究已有很多,韩延本,韩刚用小波分析的方法对太阳黑子数变化进行研究,验证了小波分析方法的可行性,并得到太阳黑子数变化包含多种周期分量的结论。郝立生,李新,李月英利用Morlet 小波变换对太阳活动变化进行了研究,得到太阳活动存在141和106a 的变化周期。 小波变换的概念是1984年法国地球物理学家J 1Morlte 在分析处理地球物理勘探资料时提出来的。其数学基础是19世纪的傅里叶变换,其后理论物理学家A 1Gr oss man 采用平移和伸缩不变性建立了小波变换的理论体系。1989年S 1Mallat 提出了多分辨率分析概念,统一了在此之前的各种构造小波的方法,特别是提出了二进小波变换的快速算法,使得小波变换完全走向实用性[8]。 奇异谱分析(SS A )是对一维的时间序列进行分析的主成分分析方法。该方法适用于从短噪声时间序列中提取信息。SS A 在时空域中,通过将序列分解成元素行为模式的方法,将含在延迟坐标相空间的信息拆开,通过使用数据适应滤波器来帮助将时间序列分开为统计的独立成分,这些成分可以当作趋势、振动或噪声来进行分类。 本文选用太阳黑子数月平均值,采用小波变换和奇异谱分析的方法对该时间序列进行分析,同时对两种分析方法进行比较。 2　奇异谱分析 主成分分析(PCA,Princi pal Component Analysis ),也称为经验正交函数(E OF,E mp irical O rthogonal Functi on ), 可以由多维的时间序列中获取时间序列的主要成分,是常用的多元统计分析方法之一,主要将多个彼此相关的指标变换为少数几个彼此独立的综合指标即主成分,并要求主成分能反映原始数据的几乎全部信息,其中,常用于对一维的时间序列进行分析的方法称为奇异谱分析(SS A,Sin 2gular s pectru m analysis )。 奇异谱方法(SS A )是一种特别适合于研究周期振荡行为的分析方法,它是从时间序列的动力重构出发,并与经验正交函数相联系的一种统计技术,是E OF 分解的一特殊应用。分解的空间结构与时间尺度密切相关,可以较好地从含噪声的有限尺度时间序列中提取信息,目前已应用于多种时间序列的分析中。 SS A 的具体操作过程是,将一个样本量为n 的时间序列按给定嵌套空间维数(即窗口长度)构造一资料矩阵。当这一个资料矩阵计算出明显成对的特征值,且相应的E OF 几乎是周期性或正交时,通常就对应着信号中的振荡行为,可见SS A 在数学上相应于E OF 在延滞坐标上的表达。 对给定的X 1,X 2,…,X n 的时间序列,给定嵌套维数M ,M

小波包及能量频谱的MATLab算法

一根断条： >> %采样频率 fs=10000; nfft=10240; %定子电流信号 fid=fopen('duantiao.m','r');%故障 N=2048; xdata=fread(fid,N,'int16'); fclose(fid); xdata=(xdata-mean(xdata))/std(xdata,1); %功率谱 figure(1); Y=abs(fft(xdata,nfft)); plot((0:nfft/2-1)/nfft*fs,Y(1:nfft/2)); xlabel('频率f/Hz'); ylabel('功率谱P/W'); %3层小波包分解 T=wpdec(xdata,3,'db4'); %重构低频信号 y1=wprcoef(T,[3,1]); %y1的波形

figure(2); subplot(2,2,1); plot(1:N,y1); xlabel('时间t/n'); ylabel('电流I/A'); %y1的功率谱 Y1=abs(fft(y1,nfft)); subplot(2,2,2); plot((0:nfft/2-1)/nfft*fs,Y1(1:nfft/2)); xlabel('频率f/Hz'); ylabel('功率谱P/W'); 图像的频率是表征图像中灰度变化剧烈程度的指标，是灰度在平面空间上的梯度。如：大面积的沙漠在图像中是一片灰度变化缓慢的区域，对应的频率值很低；而对于地表属性变换剧烈的边缘区域在图像中是一片灰度变化剧烈的区域，对应的频率值较高。傅立叶变换在实际中有非常明显的物理意义，设f是一个能量有限的模拟信号，则其傅立叶变换就表示f的谱。从纯粹的数学意义上看，傅立叶变换是将一个函数转换为一系列周期函数来处理的。从物理效果看，傅立叶变换是将图像从空间域转换到频率域，其逆变换是将图像从频率域转换到空间域。换句话说，傅立叶变换的物理意义是将图像的灰度分布函数变换为图像的频率分布函数，傅立叶逆变换是将图像的频率分布函数变换为灰度分布函数。 这样通过观察傅立叶变换后的频谱图，也叫功率图，我们首先就可以看出，图像的能量分布，如果频谱图中暗的点数更多，那么实际图像是比较柔和的（因

小波包分解

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%小波包分解程序%% m=load('300_30.txt'); 导入文件名为300_30的txt文件N=length(m); for i=1:N-1 ; q(i,1)=m(i,1); end; d=q'; s1=d; change=1000; [c,l] = wavedec(d,3,'db4'); %提取小波分解后的低频系数 ca3=appcoef(c,l,'db4',3); %提取各层小波分解后的高频系数cd3=detcoef(c,l,3); cd2=detcoef(c,l,2); cd1=detcoef(c,l,1); %对信号强制消噪 cdd3=zeros(1,length(cd3));%第三层高频系数cd3全置0 cdd2=zeros(1,length(cd2));%第二层高频系数cd2全置0 cdd1=zeros(1,length(cd1));%第一层高频系数cd1全置0 c1=[ca3,cdd3,cdd2,cdd1];%建立新的系数矩阵 s2=waverec(c1,l,'db4')%为新的分解结构 %[thr,sorh,keepapp]=ddencmp('den','wv',d); %s2=wdencmp('gbl',c,l,'db4',4,thr,sorh,keepapp); %subplot(413) %plot(1:change,s2(1:change)); %title('默认软阈值消噪后信号') figure(1) subplot(9,2,1) plot(1:change,s1(1:change)) title('原始信号') ylabel('S1') subplot(9,2,2) plot(1:change,s2(1:change)) title('强制消噪后信号') ylabel('S2') wpt=wpdec(s1,3,'db1','shannon'); %plot(wpt); %重构第三层个节点小波系数 s130=wprcoef(wpt,[3,0]);

基于小波包的图像压缩及matlab实现

基于小波包的图像压缩及matlab实现 摘要：小波包分析理论作为新的时频分析工具，在信号分析和处理中得到了很好的应用，它在信号处理、模式识别、图像分析、数据压缩、语音识别与合成等等许多方面都取得了很有意义的研究成果。平面图像可以看成是二维信号，因此，小波包分析很自然地应用到了图像处理领域，如在图像的压缩编码、图像消噪、图像增强以及图像融合等方面都很好的应用。本文将对小波包分析在图像处理中的应用作以简单介绍。 关键词：小波包图像处理消噪 1.小波包基本理论 1.1 小波包用于图像消噪 图像在采集、传输等过程中，经常受到一些外部环境的影响，从而产生噪声使得图像发生降质，图像消噪的目的就是从所得到的降质图像中去除噪声还原原始图像。图像降噪是图像预处理中一项应用比较广泛的技术，其作用是为了提高图像的信噪比突出图像的期望特征。图像降噪方法有时域和频域两种方法。频率域方法主要是根据图像像素噪声频率范围，选取适当的频域带通过滤波器进行滤波处理，比如采用Fourier变换（快速算法FFT）分析或小波变换（快速算法Mallat 算法）分析。空间域方法主要采用各种平滑函数对图像进行卷积处理，以达到去除噪声的目的，如邻域平均、中值（Median）滤波等都属于这一类方法。还有建立在统计基础上的lee滤波、Kuan滤波等。但是归根到底都是利用噪声和信号在频域上分布不同进行的：信号主要分布在低频区域。而噪声主要分布在高频区域，但同时图像的细节也分布在高频区域。所以，图像降噪的一个两难问题就是如何在降低图像噪声和保留图像细节上保持平衡，传统的低通滤波方法将图像的高频部分滤除，虽然能够达到降低噪声的效果，但破坏了图像细节。如何构造一种既能够降低图像噪声，又能保持图像细节的降噪方法成为此项研究的主题。在小波变换这种有力工具出现之后，这一目标已经成为可能。 基于小波包变换消噪方法的主要思想就是利用小波分析的多尺度特性，首先对含有噪声的图像进行小波变换，然后对得到的小波系数进行阈值化处理，得到

小波包变换 matlab

小波包分解与重构 function wpt= wavelet_packetdecomposition_reconstruct( x,n,wpname ) %% 对信号进行小波包分解，得到节点的小波包系数。然后对每个节点系数进行重构。% Decompose x at depth n with wpname wavelet https://www.doczj.com/doc/dc10572879.html,ing Shannon entropy. % % x-input signal,列向量。 % n-the number of decomposition layers % wpname-a particular wavelet.type:string. % %Author hubery_zhang %Date 20170714 %% wpt=wpdec(x,n,wpname); % Plot wavelet packet tree (binary tree) plot(wpt) %% wavelet packet coefficients.default：use the front 4. cfs0=wpcoef(wpt,[n 0]); cfs1=wpcoef(wpt,[n 1]); cfs2=wpcoef(wpt,[n 2]); cfs3=wpcoef(wpt,[n 3]); figure; subplot(5,1,1); plot(x); title('原始信号'); subplot(5,1,2); plot(cfs0); title(['结点',num2str(n) ' 1',' 系数']) subplot(5,1,3); plot(cfs1); title(['结点',num2str(n) ' 2',' 系数']) subplot(5,1,4); plot(cfs2); title(['结点',num2str(n) ' 3',' 系数']) subplot(5,1,5); plot(cfs3); title(['结点',num2str(n) ' 4',' 系数']) %% reconstruct wavelet packet coefficients. rex0=wprcoef(wpt,[n 0]); rex1=wprcoef(wpt,[n 1]); rex2=wprcoef(wpt,[n 2]); rex3=wprcoef(wpt,[n 3]); figure; subplot(5,1,1);

小波包分解

一、 首先，小波包的一些基本的基本要弄懂，就是小波包是从原始信号，分级向下分解。如下 图所示。 这就是小波包树，其中节点的命名规则是从（1，0）开始，叫1号，（1，1）是2号，，，，依此类推，（3，0）是7号，（3，7）是14号。每个节点都有对应的小波包 系数，这个系数决定了频率的大小，也就是说频率信息已经有了，但是时域信息在哪里呢？那就是 order。这个order就是这些节点的顺序，也就是频率的顺序。 比如，节点的排序是 1，2，3，，，，14，那么频率就按先1号的频率变化，后2号的， 再3号的，，，然后14号的。 图1 来看一个实例： 采样频率为1024Hz，采样时间是1秒，有一个信号s是由频率100和200Hz的正弦波混合的，我们用小波包来分解。 clear all clc fs=1024; %采样频率 f1=100; %信号的第一个频率 f2=300; %信号第二个频率

t=0:1/fs:1; s=sin(2*pi*f1*t)+sin(2*pi*f2*t); %生成混合信号 [tt]=wpdec(s,3,'dmey'); %小波包分解，3代表分解3层，像图1那样，'dmey'使用meyr小波plot(tt) %这个就是画出图1那个图，可以用鼠标在上面点 wpviewcf(tt,1); %画出时间频率图，如图2 图2 现在开始解释：x轴很简单，就是1024个点，对应1秒，每个点就代表1/1024秒，x轴诡 异一下，最后一个数就是1. y轴上显示的数字对应于图1 中的节点，从下面开始，顺序是 7号节点，8号，10号，9号，，，，11号节点，这个顺序是这么排列的，这是小波包自动排列的，不用管。只要知道怎么查看这个order就可以了。然后，y轴是频率啊，怎么不是100Hz和300Hz呢？原因就是MATLAB这里没有显示频率，显示的是order，频率我们要 自己算，怎么算呢。我们的采样频率是1024Hz，根据采样定理，奈奎斯特采样频率是 512Hz，我们分解了3层，最后一层就是 2^3=8个频率段，每个频率段的频率区间是 512/8=64Hz,对吧，那看图2，颜色重的地方一个是在8那里，一个在13那里，8是第二段，也就是 65-128Hz之间，13是第五段，也就是257-320Hz之间。这样就说通了，正好这个 原始信号只有两个频率段，一个100一个300。如果我们不是分解了3层，而是更多层， 那么每个频率段包含的频率也就越窄，图上有颜色的地方也会更细，也就是说更精细了， 大家可以自己试试。将3改为6试试。由于原始信号的频率在整个1秒钟内都没有改变， 所以有颜色的地方是一个横线。

小波包变换及matlab程序编写

1 小波变换的基本理论 信号分析是为了获得时间和频率之间的相互关系。小波变换（DWT ）是现代谱分析工具，他既能考察局部时域过程的频域特征，又能考察局部频域过程的时域特征，因此即使对于非平稳过程，处理起来也得心应手。傅立叶变换提供了有关频率域的信息，但有关时间的局部化信息却基本丢失。与傅立叶变换不同，小波变换能将图像变换为一系列小波系数，这些系数可以被高效压缩和存储，此外，小波的粗略边缘可以更好地表现图像，因为他消除了DCT 压缩普遍具有的方块效应。通过缩放母小波（Mother wavelet ）的宽度来获得信号的频率特征， 通过平移母小波来获得信号的时间信息。对母小波的缩放和平移操作是为了计算小波系数，这些小波系数反映了小波和局部信号之间的相关程度。 小波变换是当前应用数学中一个迅速发展的领域，是分析和处理非平稳信号的一种有力工具。它是以局部化函数所形成的小波基作为基底展开的，具有许多特殊的性能和优点，小波分析是一种更合理的进频表示和子带多分辨分析。 2小波包变换的基本理论和原理 概论：由于正交小波变换只对信号的低频部分做进一步分解，而对高频部分也即信号的细节部分不再继续分解，所以小波变换能够很好地表征一大类以低频信息为主要成分的信号，但它不能很好地分解和表示包含大量细节信息（细小边缘或纹理）的信号，如非平稳机械振动信号、遥感图象、地震信号和生物医学信号等。与之不同的是，小波包变换可以对高频部分提供更精细的分解，而且这种分解既无冗余，也无疏漏，所以对包含大量中、高频信息的信号能够进行更好的时频局部化分析。 2.1小波包的定义： 正交小波包的一般解释 仅考虑实系数滤波器. {}n n Z h ∈{}n n Z g ∈()11n n n g h -=-( )()( )()22k k Z k k Z t h t k t g t k φφψφ∈∈?=-?? =-??