关于小波分析的matlab程序
（参考）小波分析及其matlab实现
第8章 小波分析及其MATLAB 实现本文档节选自：Matlab 在数学建模中的应用，卓金武等编著， 北航出版社，2011年4月出版8.1小波分析基本理论8.1.1 Fourier 变换的局限性8.1.1.1变换的含义我们把那些定义域和因变域都不是数值或常量的函数称为变换或算子，它们是定义域和值域本身为函数集的函数，如傅里叶变换（Fourier Transform ）和拉普拉斯变换（Laplace Transform ），其定义域是时间的函数，而因变域是频率的函数。
简单地说，变换的基本思想仍然是映射，变换是函数的函数。
8.1.1.2 Fourier 变换的局限性信号分析的主要目的是寻找一种简单有效的信号变换方法，使信号包含的重要特征能显示出来。
在小波变换兴起之前，Fourier 级数展开和Fourier 分析是刻画函数空间、求解微分方程、进行数学计算的有效方法和有效的数学工具。
从物理直观上看，一个周期性振动的量可以看成是具有简单频率的简谐振动的叠加。
Fourier 级数展开则是这一物理过程的数学描述，Fourier 变化和Fourier 逆变换公式如下：函数)()(1R L t f ∈的连续Fourier 变换定义为t t f ft d e )()(ˆi -⎰+∞∞-=ωω)(ˆωf 的Fourier 逆变换定义为 ωωωd e )(ˆπ21)(i ⎰+∞∞-=t f t f 从公式上看，Fourier 变换是域变换，它把时间域和频率域联系起来，在时间域上难以观察到的现象和规律，在频域往往能十分清楚地显示出来。
频谱分析的本质就是对)(ˆωf的加工、分析和滤波等处理。
Fourier 变换是平稳信号分析的最重要的工具。
然而在实际运用中，所遇到的信号大多数并不平稳，而是时变频率信号，这时人们需要知道信号在突变时刻所对应的频率成分，显然Fourier 变换的积分作用平滑了非平稳过程的突变部分，作为积分核tω-i e的幅值在任何情况下均为1，因此，在频谱)(ˆωf的任一频点值是由时间过程)(t f 在整个时间域上)(∞+-∞，上的贡献决定的；反之，时间过程)(t f 在某一时刻的状态也是由)(ˆωf在整个频率域)(∞+-∞，上贡献决定的。
关于小波分析的matlab程序
关于小波分析的matlab程序（个人收集关于小波分析的matlab程序）小波滤波器构造和消噪程序 (2)小波谱分析mallat算法经典程序 (16)小波包变换分析信号的MATLAB程序 (20)利用小波变换实现对电能质量检测的算法实现 (31)基于小波变换的图象去噪Normalshrink算法 (35)小波滤波器构造和消噪程序1．重构% mallet_wavelet.m% 此函数用于研究Mallet算法及滤波器设计% 此函数仅用于消噪a=pi/8; %角度赋初值b=pi/8;%低通重构FIR滤波器h0(n)冲激响应赋值h0=cos(a)*cos(b);h1=sin(a)*cos(b);h2=-sin(a)*sin(b);h3=cos(a)*sin(b);low_construct=[h0,h1,h2,h3];L_fre=4; %滤波器长度low_decompose=low_construct(end:-1:1); %确定h0(-n)，低通分解滤波器for i_high=1:L_fre; %确定h1(n)=(-1)^n,高通重建滤波器if(mod(i_high,2)==0);coefficient=-1;elsecoefficient=1;endhigh_construct(1,i_high)=low_decompose(1,i_hi gh)*coefficient;endhigh_decompose=high_construct(end:-1:1); %高通分解滤波器h1(-n)L_signal=100; %信号长度n=1:L_signal; %信号赋值f=10;t=0.001;y=10*cos(2*pi*50*n*t).*exp(-20*n*t);figure(1);plot(y);title('原信号');check1=sum(high_decompose); %h0(n)性质校验check2=sum(low_decompose);check3=norm(high_decompose);check4=norm(low_decompose);l_fre=conv(y,low_decompose); %卷积l_fre_down=dyaddown(l_fre); %抽取,得低频细节h_fre=conv(y,high_decompose);h_fre_down=dyaddown(h_fre); %信号高频细节figure(2);subplot(2,1,1)plot(l_fre_down);title('小波分解的低频系数');subplot(2,1,2);plot(h_fre_down);title('小波分解的高频系数');l_fre_pull=dyadup(l_fre_down); %0差值h_fre_pull=dyadup(h_fre_down);l_fre_denoise=conv(low_construct,l_fre_pull); h_fre_denoise=conv(high_construct,h_fre_pull);l_fre_keep=wkeep(l_fre_denoise,L_signal); %取结果的中心部分,消除卷积影响h_fre_keep=wkeep(h_fre_denoise,L_signal);sig_denoise=l_fre_keep+h_fre_keep; %信号重构compare=sig_denoise-y; %与原信号比较figure(3);subplot(3,1,1)plot(y);ylabel('y'); %原信号subplot(3,1,2);plot(sig_denoise);ylabel('sig\_denoise'); %重构信号subplot(3,1,3);plot(compare);ylabel('compare'); %原信号与消噪后信号的比较2．消噪、% 此函数用于研究Mallet算法及滤波器设计% 此函数用于消噪处理%角度赋值%此处赋值使滤波器系数恰为db9%分解的高频系数采用db9较好,即它的消失矩较大%分解的有用信号小波高频系数基本趋于零%对于噪声信号高频分解系数很大，便于阈值消噪处理[l,h]=wfilters('db10','d');low_construct=l;L_fre=20; %滤波器长度low_decompose=low_construct(end:-1:1); %确定h0(-n)，低通分解滤波器for i_high=1:L_fre; %确定h1(n)=(-1)^n,高通重建滤波器if(mod(i_high,2)==0);coefficient=-1;elsecoefficient=1;endhigh_construct(1,i_high)=low_decompose(1,i_hi gh)*coefficient;endhigh_decompose=high_construct(end:-1:1); %高通分解滤波器h1(-n)L_signal=100; %信号长度n=1:L_signal; %原始信号赋值f=10;t=0.001;y=10*cos(2*pi*50*n*t).*exp(-30*n*t); zero1=zeros(1,60); %信号加噪声信号产生zero2=zeros(1,30);noise=[zero1,3*(randn(1,10)-0.5),zero2]; y_noise=y+noise;figure(1);subplot(2,1,1);plot(y);title('原信号');subplot(2,1,2);plot(y_noise);title('受噪声污染的信号');check1=sum(high_decompose); %h0(n),性质校验check2=sum(low_decompose);check3=norm(high_decompose);check4=norm(low_decompose);l_fre=conv(y_noise,low_decompose); %卷积l_fre_down=dyaddown(l_fre); %抽取,得低频细节h_fre=conv(y_noise,high_decompose);h_fre_down=dyaddown(h_fre); %信号高频细节figure(2);subplot(2,1,1)plot(l_fre_down);title('小波分解的低频系数');subplot(2,1,2);plot(h_fre_down);title('小波分解的高频系数');% 消噪处理for i_decrease=31:44;if abs(h_fre_down(1,i_decrease))>=0.000001 h_fre_down(1,i_decrease)=(10^-7);endendl_fre_pull=dyadup(l_fre_down); %0差值h_fre_pull=dyadup(h_fre_down);l_fre_denoise=conv(low_construct,l_fre_pull); h_fre_denoise=conv(high_construct,h_fre_pull);l_fre_keep=wkeep(l_fre_denoise,L_signal); %取结果的中心部分,消除卷积影响h_fre_keep=wkeep(h_fre_denoise,L_signal);sig_denoise=l_fre_keep+h_fre_keep; %消噪后信号重构%平滑处理for j=1:2for i=60:70;sig_denoise(i)=sig_denoise(i-2)+sig_denoise(i+2)/ 2;end;end;compare=sig_denoise-y; %与原信号比较figure(3);subplot(3,1,1)plot(y);ylabel('y'); %原信号subplot(3,1,2);plot(sig_denoise);ylabel('sig\_denoise'); %消噪后信号subplot(3,1,3);plot(compare);ylabel('compare'); %原信号与消噪后信号的比较小波谱分析mallat算法经典程序clc;clear;%% 1.正弦波定义f1=50; % 频率1f2=100; % 频率2fs=2*(f1+f2); % 采样频率Ts=1/fs; % 采样间隔N=120; % 采样点数n=1:N;y=sin(2*pi*f1*n*Ts)+sin(2*pi*f2*n*Ts); % 正弦波混合figure(1)plot(y);title('两个正弦信号')figure(2)stem(abs(fft(y)));title('两信号频谱')%% 2.小波滤波器谱分析h=wfilters('db30','l'); % 低通g=wfilters('db30','h'); % 高通h=[h,zeros(1,N-length(h))]; % 补零（圆周卷积，且增大分辨率变于观察）g=[g,zeros(1,N-length(g))]; % 补零（圆周卷积，且增大分辨率变于观察）figure(3);stem(abs(fft(h)));title('低通滤波器图')figure(4);stem(abs(fft(g)));title('高通滤波器图')%% 3.MALLET分解算法(圆周卷积的快速傅里叶变换实现)sig1=ifft(fft(y).*fft(h)); % 低通(低频分量) sig2=ifft(fft(y).*fft(g)); % 高通(高频分量) figure(5); % 信号图subplot(2,1,1)plot(real(sig1));title('分解信号1')subplot(2,1,2)plot(real(sig2));title('分解信号2')figure(6); % 频谱图subplot(2,1,1)stem(abs(fft(sig1)));title('分解信号1频谱')subplot(2,1,2)stem(abs(fft(sig2)));title('分解信号2频谱')%% 4.MALLET重构算法sig1=dyaddown(sig1); % 2抽取sig2=dyaddown(sig2); % 2抽取sig1=dyadup(sig1); % 2插值sig2=dyadup(sig2); % 2插值sig1=sig1(1,[1:N]); % 去掉最后一个零sig2=sig2(1,[1:N]); % 去掉最后一个零hr=h(end:-1:1); % 重构低通gr=g(end:-1:1); % 重构高通hr=circshift(hr',1)'; % 位置调整圆周右移一位gr=circshift(gr',1)'; % 位置调整圆周右移一位sig1=ifft(fft(hr).*fft(sig1)); % 低频sig2=ifft(fft(gr).*fft(sig2)); % 高频sig=sig1+sig2; % 源信号%% 5.比较figure(7);subplot(2,1,1)plot(real(sig1));title('重构低频信号');subplot(2,1,2)plot(real(sig2));title('重构高频信号');figure(8);subplot(2,1,1)stem(abs(fft(sig1)));title('重构低频信号频谱');subplot(2,1,2)stem(abs(fft(sig2)));title('重构高频信号频谱');figure(9)plot(real(sig),'r','linewidth',2);hold on;plot(y);legend('重构信号','原始信号')title('重构信号与原始信号比较')小波包变换分析信号的MATLAB程序%t=0.001:0.001:1;t=1:1000;s1=sin(2*pi*50*t*0.001)+sin(2*pi*120*t*0.001)+ rand(1,length(t));for t=1:500;s2(t)=sin(2*pi*50*t*0.001)+sin(2*pi*120*t*0.00 1)+rand(1,length(t));endfor t=501:1000;s2(t)=sin(2*pi*200*t*0.001)+sin(2*pi*120*t*0.001)+rand(1,length(t));endsubplot(9,2,1)plot(s1)title('原始信号')ylabel('S1')subplot(9,2,2)plot(s2)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(s130);ylabel('S130');subplot(9,2,5);plot(s131);ylabel('S131');subplot(9,2,7);plot(s132);ylabel('S132');subplot(9,2,9);plot(s133);ylabel('S133');subplot(9,2,11);plot(s134);ylabel('S134');subplot(9,2,13);plot(s135);ylabel('S135');subplot(9,2,15);plot(s136);ylabel('S136');subplot(9,2,17);plot(s137);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(s230);ylabel('S230');subplot(9,2,6);plot(s231);ylabel('S231');subplot(9,2,8);plot(s232);ylabel('S232');subplot(9,2,10);plot(s233);ylabel('S233');subplot(9,2,12);plot(s234);ylabel('S234');subplot(9,2,14);plot(s235); ylabel('S235');subplot(9,2,16);plot(s236); ylabel('S236');subplot(9,2,18);plot(s237); ylabel('S237');%fftfigurey1=fft(s1,1024);py1=y1.*conj(y1)/1024;y2=fft(s2,1024);py2=y2.*conj(y2)/1024;y130=fft(s130,1024);py130=y130.*conj(y130)/1024; y131=fft(s131,1024);py131=y131.*conj(y131)/1024; y132=fft(s132,1024);py132=y132.*conj(y132)/1024; y133=fft(s133,1024);py133=y133.*conj(y133)/1024; y134=fft(s134,1024);py134=y134.*conj(y134)/1024; y135=fft(s135,1024);py135=y135.*conj(y135)/1024; y136=fft(s136,1024);py136=y136.*conj(y136)/1024; y137=fft(s137,1024);py137=y137.*conj(y137)/1024;y230=fft(s230,1024);py230=y230.*conj(y230)/1024; y231=fft(s231,1024);py231=y231.*conj(y231)/1024; y232=fft(s232,1024);py232=y232.*conj(y232)/1024; y233=fft(s233,1024);py233=y233.*conj(y233)/1024; y234=fft(s234,1024);py234=y234.*conj(y234)/1024; y235=fft(s235,1024);py235=y235.*conj(y235)/1024;y236=fft(s236,1024);py236=y236.*conj(y236)/1024; y237=fft(s237,1024);py237=y237.*conj(y237)/1024;f=1000*(0:511)/1024;subplot(1,2,1);plot(f,py1(1:512));ylabel('P1');title('原始信号的功率谱') subplot(1,2,2);plot(f,py2(1:512));ylabel('P2');title('故障信号的功率谱') figuresubplot(4,2,1);plot(f,py130(1:512));ylabel('P130');title('S130的功率谱') subplot(4,2,2);plot(f,py131(1:512)); ylabel('P131');title('S131的功率谱') subplot(4,2,3);plot(f,py132(1:512)); ylabel('P132'); subplot(4,2,4);plot(f,py133(1:512)); ylabel('P133'); subplot(4,2,5);plot(f,py134(1:512)); ylabel('P134'); subplot(4,2,6);plot(f,py135(1:512)); ylabel('P135'); subplot(4,2,7);plot(f,py136(1:512)); ylabel('P136'); subplot(4,2,8);plot(f,py137(1:512)); ylabel('P137'); figuresubplot(4,2,1);plot(f,py230(1:512)); ylabel('P230');title('S230的功率谱') subplot(4,2,2);plot(f,py231(1:512)); ylabel('P231');title('S231的功率谱') subplot(4,2,3);plot(f,py232(1:512)); ylabel('P232'); subplot(4,2,4);plot(f,py233(1:512)); ylabel('P233'); subplot(4,2,5);plot(f,py234(1:512)); ylabel('P234'); subplot(4,2,6);plot(f,py235(1:512)); ylabel('P235'); subplot(4,2,7);plot(f,py236(1:512));ylabel('P236');subplot(4,2,8);plot(f,py237(1:512));ylabel('P237');figure%plottree(wpt)利用小波变换实现对电能质量检测的算法实现N=10000;s=zeros(1,N);for n=1:Nif n<0.4*N||n>0.8*Ns(n)=31.1*sin(2*pi*50/10000*n);elses(n)=22.5*sin(2*pi*50/10000*n);endendl=length(s);[c,l]=wavedec(s,6,'db5'); %用db5小波分解信号到第六层subplot(8,1,1);plot(s);title('用db5小波分解六层：s=a6+d6+d5+d4+d3+d2+d1');Ylabel('s');%对分解结构【c，l】中第六层低频部分进行重构a6=wrcoef('a',c,l,'db5',6);subplot(8,1,2);plot(a6);Ylabel('a6');%对分解结构【c，l】中各层高频部分进行重构for i=1:6decmp=wrcoef('d',c,l,'db5',7-i);subplot(8,1,i+2);plot(decmp);Ylabel(['d',num2str(7-i)]);end%----------------------------------------------------------- rec=zeros(1,300);rect=zeros(1,300);ke=1;u=0;d1=wrcoef('d',c,l,'db5',1);figure(2);plot(d1);si=0;N1=0;N0=0;sce=0;for n=20:N-30rect(ke)=s(n);ke=ke+1;if(ke>=301)if(si==2)rec=rect;u=2;end;si=0;ke=1;end;if(d1(n)>0.01) % the condition of abnormal append.N1=n;if(N0==0)N0=n;si=si+1;end;if(N1>N0+30)Nlen=N1-N0;Tab=Nlen/10000;end;end;if(si==1)fork=N0:N0+99 %testin g of 1/4 period signals tosce=sce+s(k)*s(k)/10000;end;re=sqrt(sce*200) %re indicate the pike value of .sce=0;si=si+1;end;end;NlenN0n=1:300;figure(3)plot(n,rec);基于小波变换的图象去噪Normalshrink算法function[T_img,Sub_T]=threshold_2_N(img,levels)% reference :image denoising using wavelet thresholding[xx,yy]=size(img);HH=img((xx/2+1):xx,(yy/2+1):yy);delt_2=(std(HH(:)))^2;%(median(abs(HH(:)))/0. 6745)^2;%T_img=img;for i=1:levelstemp_x=xx/2^i;temp_y=yy/2^i;% belt=1.0*(log(temp_x/(2*levels)))^0.5;belt=1.0*(log(temp_x/(2*levels)))^0.5; %2.5 0.8%HLHL=img(1:temp_x,(temp_y+1):2*temp_y);delt_y=std(HL(:));T_1=belt*delt_2/delt_y;%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%% %%%%%%%%%%T_HL=sign(HL).*max(0,abs(HL)-T_1);%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%% %%%%%%%%%T_img(1:temp_x,(temp_y+1):2*temp_y)=T_HL;Sub_T(3*(i-1)+1)=T_1;%LHLH=img((temp_x+1):2*temp_x,1:temp_y);delt_y=std(LH(:));T_2=belt*delt_2/delt_y;%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%% %%%%%%%%%T_LH=sign(LH).*max(0,abs(LH)-T_2);%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%T_img((temp_x+1):2*temp_x,1:temp_y)=T_LH;Sub_T(3*(i-1)+2)=T_2;%HHHH=img((temp_x+1):2*temp_x,(temp_y+1):2*te mp_y);delt_y=std(HH(:));T_3=belt*delt_2/delt_y;%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%T_HH=sign(HH).*max(0,abs(HH)-T_3);%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%T_img((temp_x+1):2*temp_x,(temp_y+1):2*tem p_y)=T_HH;Sub_T(3*(i-1)+3)=T_3;end。
利用dwt实现信号的小波分解 matlab代码
利用dwt实现信号的小波分解 matlab代码
一、简介
小波变换是一种信号处理方法，它通过伸缩和平移小波函数来分析信号的局部特性。
离散小波变换（Discrete Wavelet Transform，DWT）是小波变换的一种形式，它可以用于信号的多尺度分析。
在MATLAB中，我们可以使用内置的小波分析工具箱来实现DWT。
二、代码实现
这段代码首先生成一个简单的正弦信号，然后使用Daubechies 1小波对其进行DWT分解。
最后，它绘制出原始信号以及各尺度的小波系数。
请注意，您需要先安装小波工具箱才能运行此代码。
三、说明与注意事项
此代码仅用于演示DWT的基本概念和实现。
在实际应用中，您可能需要根据具体的需求和信号特性选择合适的小波函数和分解层次。
此外，还需要考虑如何处理噪声和小波系数的阈值压缩等问题。
matlab小波分解程序
matlab小波分解程序小波分解是一种信号处理的方法，可以用于信号的分析和压缩。
在MATLAB中，可以使用内置的`wavedec`函数来进行小波分解。
下面是一个简单的MATLAB小波分解程序示例：matlab.% 创建一个示例信号。
x = randn(1,1024);% 选择小波基和分解级别。
wname = 'db4'; % 选择小波基，这里使用db4小波。
level = 3; % 选择分解级别。
% 进行小波分解。
[c, l] = wavedec(x, level, wname);% 从分解系数和长度信息中重构近似和细节系数。
appx = wrcoef('a',c,l,wname,level); % 近似系数。
det1 = wrcoef('d',c,l,wname,1); % 第一层细节系数。
det2 = wrcoef('d',c,l,wname,2); % 第二层细节系数。
det3 = wrcoef('d',c,l,wname,3); % 第三层细节系数。
% 绘制原始信号和重构的近似信号。
t = 1:1024;subplot(2,1,1);plot(t, x);title('Original Signal');subplot(2,1,2);plot(t, appx);title('Approximation Coefficients'); % 显示细节系数。
figure;subplot(3,1,1);plot(t, det1);title('Detail Coefficients Level 1'); subplot(3,1,2);plot(t, det2);title('Detail Coefficients Level 2'); subplot(3,1,3);plot(t, det3);title('Detail Coefficients Level 3');在这个示例中，我们首先生成了一个长度为1024的随机信号。
小波分析中的matlab使用
小波分析中的matlab使用Matlab主窗口File菜单File菜单,弹出如图1所示的菜单选项。
其中,各子菜单选项的功能如下:图 1New选项包含5个选项:M-File,Figure,Varible，Model和gui。
1)M-File选项:打开m文件编辑器;2)Figure选项:将打开一个空白的图形窗口;3)Variable选项：可变因素；4)Model选项:用于创建新模型的窗口;5)Gui选项:创建新的图形用户界面的对话框。
Open选项:打开一个open对话框,可以在对话框中选择相应的文件,然后matlab将用相应的编辑器打开该文件。
Close…选项:跟随某个打开的视窗名。
单击该选项,将关闭该视窗。
Importdata…选项:打开一个import对话框,用户可以选择相应的数据文件,然后将该数据文件中的数据导入到matlab工作空间。
Saveworkspaceas…选项:打开一个savetomat-File对话框,用户需要为保存的工作空间命名。
Setpath…选项:打开设置路径对话框。
通过该对话框可以更改matlab执行命令时搜索的路径。
Preferences：首选参数。
Pagesetup选项:用于设置页面布局,页面的页眉,页面所用的文字。
Print…选项:用于打印预定义好的页面内容,也可以设置一些参数。
Printselection…选项:当选中命令窗口内的一部分内容后,该选项将处于激活状态,此时单击该选项,将打印对话框中选中的内容。
Exitmatlab选项:关闭matlab。
也可以通过快捷键ctrl+O来关闭。
Edit菜单单击edit菜单,会弹出如图2所示的菜单选项。
其中,各子菜单选项的功能如下:Undo选项：取消上一次的操作。
Redo选项：重复上一次的操作。
Cut选项：剪切所选中的部分。
Copy选项：选复制被选中的部分。
Paste选项：把存放在缓冲区中的内容粘贴到光标所在的位置。
Pastespecial选项：打开导入数据向导,该向导引导用户把存放在缓冲区中的内容以特定格式存放到剪贴板变量中。
小波变换的原理及matlab仿真程序
基于小波变换的信号降噪研究2 小波分析基本理论设Ψ(t)∈L 2( R) ( L 2( R) 表示平方可积的实数空间，即能量有限的信号空间) ， 其傅立叶变换为Ψ(t)。
当Ψ(t)满足条件[4,7]：2()Rt dw wCψψ=<∞⎰（1）时，我们称Ψ(t)为一个基本小波或母小波，将母小波函数Ψ(t)经伸缩和平移后，就可以得到一个小波序列：,()()a bt bt aψ-=,,0a b R a ∈≠ （2） 其中a 为伸缩因子，b 为平移因子。
对于任意的函数f(t)∈L 2( R)的连续小波变换为：,(,),()()f a b Rt bW a b f f t dt aψψ-=<>=⎰（3） 其逆变换为：211()(,)()fR R t b f t W a b dadb C a aψψ+-=⎰⎰ （4） 小波变换的时频窗是可以由伸缩因子a 和平移因子b 来调节的，平移因子b,可以改变窗口在相平面时间轴上的位置，而伸缩因子b 的大小不仅能影响窗口在频率轴上的位置,还能改变窗口的形状。
小波变换对不同的频率在时域上的取样步长是可调节的，在低频时，小波变换的时间分辨率较低，频率分辨率较高：在高频时，小波变换的时间分辨率较高，而频率分辨率较低。
使用小波变换处理信号时，首先选取适当的小波函数对信号进行分解，其次对分解出的参数进行阈值处理，选取合适的阈值进行分析，最后利用处理后的参数进行逆小波变换，对信号进行重构。
3 小波降噪的原理和方法3.1 小波降噪原理从信号学的角度看 ,小波去噪是一个信号滤波的问题。
尽管在很大程度上小波去噪可以看成是低通滤波 ,但由于在去噪后 ,还能成功地保留信号特征 ,所以在这一点上又优于传统的低通滤波器。
由此可见 ,小波去噪实际上是特征提取和低通滤波的综合 ,其流程框图如图所示[6]：小波分析的重要应用之一就是用于信号消噪 ,一个含噪的一维信号模型可表示为如下形式：(k)()()S f k e k ε=+* k=0.1…….n-1其中 ,f( k)为有用信号,s(k)为含噪声信号,e(k)为噪声,ε为噪声系数的标准偏差。
小波分析MATLAB实例
小波分析MATLAB实例小波分析是一种信号处理方法，可以用于信号的时频分析和多尺度分析。
在MATLAB中，可以使用Wavelet Toolbox实现小波分析。
这个工具箱提供了丰富的函数和工具，可以方便地进行小波分析的计算和可视化。
小波分析的核心是小波变换，它将信号分解成一组不同尺度和频率的小波基函数。
在MATLAB中，可以使用`cwt`函数进行连续小波变换。
以下是一个小波分析的MATLAB实例，用于分析一个心电图信号的时频特性。
首先，导入心电图信号数据。
假设心电图数据保存在一个名为`ecg_signal.mat`的文件中，包含一个名为`ecg`的变量。
可以使用`load`函数加载这个数据。
```MATLABload('ecg_signal.mat');```接下来，设置小波变换的参数。
选择一个小波基函数和一组尺度。
这里选择Morlet小波作为小波基函数，选择一组从1到64的尺度。
可以使用`wavelet`函数创建一个小波对象，并使用`scal2frq`函数将尺度转换为频率。
```MATLABwavelet_name = 'morl'; % 选择Morlet小波作为小波基函数scales = 1:64; % 选择1到64的尺度wavelet_obj = wavelet(wavelet_name);scales_freq = scal2frq(scales, wavelet_name, 1);```然后，使用`cwt`函数进行小波变换，得到信号在不同尺度和频率下的小波系数。
将小波系数的幅度平方得到信号的能量谱密度。
```MATLAB[wt, f] = cwt(ecg, scales, wavelet_name);energy = abs(wt).^2;``````MATLABimagesc(1:length(ecg), scales_freq, energy);colormap('jet');xlabel('时间（样本）');ylabel('频率（Hz）');```运行整个脚本之后，就可以得到心电图信号的时频图。
第六章基于MATLAB的小波分析
第六章基于MATLAB的小波分析小波分析是一种用来分析和处理信号的数学方法，其基本原理是通过将信号分解成不同频率范围的小波基函数来揭示信号的特征。
MATLAB是一种功能强大的科学计算和数据分析软件，提供了丰富的工具箱和函数，可以方便地进行小波分析。
在MATLAB中，小波分析可以通过使用Wavelet Toolbox来实现。
该工具箱提供了几种常用的小波基函数，如Daubechies、Coiflets、Symlets等，同时还包括了一系列小波分析的函数。
下面将介绍基于MATLAB的小波分析的基本步骤。
首先，需要导入待分析的信号数据。
可以使用MATLAB的数据导入和处理工具来加载信号数据，如load函数、importdata函数等。
加载数据后，可以使用plot函数将信号数据可视化，以便直观地了解信号的特点。
接下来，需要选择合适的小波基函数进行分析。
小波基函数的选择与信号的特征和分析目标相关。
可以使用waveinfo函数来查看Wavelet Toolbox提供的小波基函数的特性和参数，并选择适合的小波基函数。
然后，使用wavedec函数对信号进行小波分解。
wavedec函数可以将信号分解成多个尺度的小波系数。
分解得到的小波系数包括近似系数和细节系数，近似系数反映了信号在低频范围的特征，而细节系数则反映了信号在高频范围的细节特征。
分解后，可以使用可视化函数如plot、imshow等来展示小波系数的分布和变化情况。
通过观察小波系数的变化，可以得到信号的频率特征和局部特征。
除了观察小波系数，还可以根据需要进行小波系数的处理和分析。
例如，可以使用细节系数来提取信号中的细节特征，如边缘、尖峰等，也可以使用近似系数来提取信号的整体趋势。
最后，可以使用waverec函数将处理过的小波系数重构成原始信号。
重构得到的信号可以与原始信号进行对比，以验证分析的结果和提取的特征。
综上所述，MATLAB提供了丰富的工具和函数来实现小波分析，可以方便地进行信号的频率分析和特征提取。
Daubechies小波分析的部分Matlab实现
实验5：Daubechies 小波分析的部分Matlab 实现注：本实验在原本操作过程中在重构部分存在着问题，经思考检验调试之后，排除了隐藏的问题，实验得到成功进行同时达到了预期的结果。
第一部分Daubechies 尺度序列，二进点上的尺度函数图像及小波函数图像第一部分实验的目标有两个：1，给出各阶Daubechies 小波的低通滤波器，即尺度函数两尺度序列。
具体的程序算法由《小波导论》中7.3节中的引理7.16，定理7.17给现，其中定理7.17的证明过程实际上给出了一个显式的计算过程，这是程序实现的重要依据。
2，根据两尺度关系，由第一步骤中得到的两尺度序列逐步计算得到Daubechies 尺度函数及小波函数的二进点值，进而给出近似图像。
第一步骤相关函数程序解释：1、sinj(n)，实现将引理7.16中的2(sin)2n ω转化为形如cos k ω，n 即表示(sin)2ω的次数。
2、sinj(n)中涉及函数cosn(n), cosn(n)实现将cos kω转化为cos()kkw ∑，n 表示三角函数次数。
程序编写的思路就是先将2(sin)2n ω转化为cos k ω再转化为cos()kkw ∑，然后根据定理7.17的证明过程计算出两尺度序列，其中用于提取多项式系数的函数利用了三角函数的正交性，使用了积分函数int()。
daucoef(n)给出n 阶daubechies 尺度函数的两尺度序列。
具体的原程序代码见附件，部分程序运行结果如下：各阶（2--12）Daubechies 尺度函数的两尺度序列如下：）第二步骤：给出n 阶daubechies 尺度函数图像，由函数dauinterp(n,max1)实现，即给出n 阶daubechies 尺度函数的max1次迭代的函数值，所得到的点都是函数的二进点。
Dauinterpwa(n,max1)则给出小波函数的二进点上的函数值。
这一部分内容中，比较难以实现的地方在于下标的变换，这常常引起一些数据的混乱或错误，相关的数学推导总是难以避免的，由于过程过于繁复，这里不给出详细的推导过程及解释。
实验三 MATLAB小波实现
实验三 MATLAB 小波实现实验题目：利用MATLAB 实现一些小波函数，并绘图实验内容：1. Haar 小波函数（1） 定义函数源程序：function Haarfor x=0:0.001:2if x>=0&x<=0.5H=1;elseif x>=0.5&x<=1H=-1;elseH=0;endplot(x,H,'k')axis([0,2,-2,2]) title('Harr 定义函数')hold onendhold offend运行程序得图：⎪⎪⎩⎪⎪⎨⎧<≤-≤≤=其它012112/101x x H ψ（2） 尺度函数源程序：function harrfor x=0:0.001:2if x>=0&x<=1H=1;elseH=0;endplot(x,H,'k')axis([0,2,-2,2]) title('Harr 尺度函数')hold onendhold offend00.20.40.60.81 1.2 1.4 1.6 1.82Harr 定义函数⎩⎨⎧≤≤=其它0101)(x x φ运行程序得图：2. Mexican Hat 小波函数函数定义为：源程序：function MexicanHat[psi,x]=mexihat(-5,5,1000);plot(x,psi,'k') title('Mexican Hat 小波')end运行程序得图：00.20.40.60.81 1.2 1.4 1.6 1.82Harr 尺度函数2/24/12e )1(π32)(x x x ---=ψ3. Meyer 小波函数（1）定义函数v (a )为构造Meyer 小波的辅助函数，且有v (a )＝a 4(35－84a ＋70a 2－20a 3)，a ∈［0，1］（2）尺度函数-5-4-3-2-1012345-0.4-0.20.20.40.60.81Mexican Hat 小波⎪⎪⎪⎩⎪⎪⎪⎨⎧∉≤≤-≤≤-=--]3π8,3π2[,03π83π4)),1π23(2πcos(e π)2(3π43π2)),1π23(2πsin(e π)2()(ˆ2/i 2/12/i 2/1ωωωνωωνωψωω⎪⎪⎪⎩⎪⎪⎪⎨⎧>≤≤-≤=--3π403π43π2))1π23(2πcos(π)2(3π2π)2()(2/12/1ωωωνωωΦ源程序：function Meyer[phi,psi,x]=meyer(-8,8,1024);figure(1);plot(x,psi,'k'),title('Meyer 定义函数') figure(2);plot(x,phi,'k'),title('Meyer 尺度函数') end运行程序得图：-8-6-4-202468-1-0.50.511.5Meyer 定义函数4. Morlet 小波函数函数定义为：源程序：function Morlet[psi,x]=morlet(-4,4,1000);plot(x,psi,'k')title('MorletС²¨')end运行程序得图：-8-6-4-202468-0.4-0.20.20.40.60.811.2Meyer 尺度函数)5cos(e )(2/2x C x x -=ψ5. Gauss 函数族 同一坐标系中ααπα4221)(x e x g -= 161,41,1=α 及其频域形式 源程序：function Gaussx=-5:0.001:5;g1=1/(2*pi^(1/2))*exp(-x.^2/4);g2=1/(pi^(1/2))*exp(-x.^2);g3=2/(pi^(1/2))*exp(-4*x.^2);plot(x,g1,x,g2,x,g3);title('Gauss 函数族')end运行程序得图：-4-3-2-101234-1-0.8-0.6-0.4-0.20.20.40.60.81Morlet 小波-5-4-3-2-101234500.20.40.60.811.21.4Gauss 函数族。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
/forum-viewthread-tid-9141-extra-page%3D1-page-1.html（个人收集关于小波分析的matlab程序）小波滤波器构造和消噪程序 (1)小波谱分析mallat算法经典程序 (7)小波包变换分析信号的MATLAB程序 (9)利用小波变换实现对电能质量检测的算法实现 (15)基于小波变换的图象去噪Normalshrink算法 (17)小波滤波器构造和消噪程序1．重构% mallet_wavelet.m% 此函数用于研究Mallet算法及滤波器设计% 此函数仅用于消噪a=pi/8; %角度赋初值b=pi/8;%低通重构FIR滤波器h0(n)冲激响应赋值h0=cos(a)*cos(b);h1=sin(a)*cos(b);h2=-sin(a)*sin(b);h3=cos(a)*sin(b);low_construct=[h0,h1,h2,h3];L_fre=4; %滤波器长度low_decompose=low_construct(end:-1:1); %确定h0(-n)，低通分解滤波器for i_high=1:L_fre; %确定h1(n)=(-1)^n,高通重建滤波器if(mod(i_high,2)==0);coefficient=-1;elsecoefficient=1;endhigh_construct(1,i_high)=low_decompose(1,i_high)*coefficient;endhigh_decompose=high_construct(end:-1:1); %高通分解滤波器h1(-n) L_signal=100; %信号长度n=1:L_signal; %信号赋值f=10;t=0.001;y=10*cos(2*pi*50*n*t).*exp(-20*n*t);figure(1);plot(y);title('原信号');check1=sum(high_decompose); %h0(n)性质校验check2=sum(low_decompose);check3=norm(high_decompose);check4=norm(low_decompose);l_fre=conv(y,low_decompose); %卷积l_fre_down=dyaddown(l_fre); %抽取,得低频细节h_fre=conv(y,high_decompose);h_fre_down=dyaddown(h_fre); %信号高频细节figure(2);subplot(2,1,1)plot(l_fre_down);title('小波分解的低频系数');subplot(2,1,2);plot(h_fre_down);title('小波分解的高频系数');l_fre_pull=dyadup(l_fre_down); %0差值h_fre_pull=dyadup(h_fre_down);l_fre_denoise=conv(low_construct,l_fre_pull);h_fre_denoise=conv(high_construct,h_fre_pull);l_fre_keep=wkeep(l_fre_denoise,L_signal); %取结果的中心部分,消除卷积影响h_fre_keep=wkeep(h_fre_denoise,L_signal);sig_denoise=l_fre_keep+h_fre_keep; %信号重构compare=sig_denoise-y; %与原信号比较figure(3);subplot(3,1,1)plot(y);ylabel('y'); %原信号subplot(3,1,2);plot(sig_denoise);ylabel('sig\_denoise'); %重构信号subplot(3,1,3);plot(compare);ylabel('compare'); %原信号与消噪后信号的比较2．消噪、% 此函数用于研究Mallet算法及滤波器设计% 此函数用于消噪处理%角度赋值%此处赋值使滤波器系数恰为db9%分解的高频系数采用db9较好,即它的消失矩较大%分解的有用信号小波高频系数基本趋于零%对于噪声信号高频分解系数很大，便于阈值消噪处理[l,h]=wfilters('db10','d');low_construct=l;L_fre=20; %滤波器长度low_decompose=low_construct(end:-1:1); %确定h0(-n)，低通分解滤波器for i_high=1:L_fre; %确定h1(n)=(-1)^n,高通重建滤波器if(mod(i_high,2)==0);coefficient=-1;elsecoefficient=1;endhigh_construct(1,i_high)=low_decompose(1,i_high)*coefficient;endhigh_decompose=high_construct(end:-1:1); %高通分解滤波器h1(-n) L_signal=100; %信号长度n=1:L_signal; %原始信号赋值f=10;t=0.001;y=10*cos(2*pi*50*n*t).*exp(-30*n*t);zero1=zeros(1,60); %信号加噪声信号产生zero2=zeros(1,30);noise=[zero1,3*(randn(1,10)-0.5),zero2];y_noise=y+noise;figure(1);subplot(2,1,1);plot(y);title('原信号');subplot(2,1,2);plot(y_noise);title('受噪声污染的信号');check1=sum(high_decompose); %h0(n),性质校验check2=sum(low_decompose);check3=norm(high_decompose);check4=norm(low_decompose);l_fre=conv(y_noise,low_decompose); %卷积l_fre_down=dyaddown(l_fre); %抽取,得低频细节h_fre=conv(y_noise,high_decompose);h_fre_down=dyaddown(h_fre); %信号高频细节figure(2);subplot(2,1,1)plot(l_fre_down);title('小波分解的低频系数');subplot(2,1,2);plot(h_fre_down);title('小波分解的高频系数');% 消噪处理for i_decrease=31:44;if abs(h_fre_down(1,i_decrease))>=0.000001h_fre_down(1,i_decrease)=(10^-7);endendl_fre_pull=dyadup(l_fre_down); %0差值h_fre_pull=dyadup(h_fre_down);l_fre_denoise=conv(low_construct,l_fre_pull);h_fre_denoise=conv(high_construct,h_fre_pull);l_fre_keep=wkeep(l_fre_denoise,L_signal); %取结果的中心部分,消除卷积影响h_fre_keep=wkeep(h_fre_denoise,L_signal);sig_denoise=l_fre_keep+h_fre_keep; %消噪后信号重构%平滑处理for j=1:2for i=60:70;sig_denoise(i)=sig_denoise(i-2)+sig_denoise(i+2)/2;end;end;compare=sig_denoise-y; %与原信号比较figure(3);subplot(3,1,1)plot(y);ylabel('y'); %原信号subplot(3,1,2);plot(sig_denoise);ylabel('sig\_denoise'); %消噪后信号subplot(3,1,3);plot(compare);ylabel('compare'); %原信号与消噪后信号的比较小波谱分析mallat算法经典程序clc;clear;%% 1.正弦波定义f1=50; % 频率1f2=100; % 频率2fs=2*(f1+f2); % 采样频率Ts=1/fs; % 采样间隔N=120; % 采样点数n=1:N;y=sin(2*pi*f1*n*Ts)+sin(2*pi*f2*n*Ts); % 正弦波混合figure(1)plot(y);title('两个正弦信号')figure(2)stem(abs(fft(y)));title('两信号频谱')%% 2.小波滤波器谱分析h=wfilters('db30','l'); % 低通g=wfilters('db30','h'); % 高通h=[h,zeros(1,N-length(h))]; % 补零（圆周卷积，且增大分辨率变于观察）g=[g,zeros(1,N-length(g))]; % 补零（圆周卷积，且增大分辨率变于观察）figure(3);stem(abs(fft(h)));title('低通滤波器图')figure(4);stem(abs(fft(g)));title('高通滤波器图')%% 3.MALLET分解算法(圆周卷积的快速傅里叶变换实现)sig1=ifft(fft(y).*fft(h)); % 低通(低频分量)sig2=ifft(fft(y).*fft(g)); % 高通(高频分量)figure(5); % 信号图subplot(2,1,1)plot(real(sig1));title('分解信号1')subplot(2,1,2)plot(real(sig2));title('分解信号2')figure(6); % 频谱图subplot(2,1,1)stem(abs(fft(sig1)));title('分解信号1频谱')subplot(2,1,2)stem(abs(fft(sig2)));title('分解信号2频谱')%% 4.MALLET重构算法sig1=dyaddown(sig1); % 2抽取sig2=dyaddown(sig2); % 2抽取sig1=dyadup(sig1); % 2插值sig2=dyadup(sig2); % 2插值sig1=sig1(1,[1:N]); % 去掉最后一个零sig2=sig2(1,[1:N]); % 去掉最后一个零hr=h(end:-1:1); % 重构低通gr=g(end:-1:1); % 重构高通hr=circshift(hr',1)'; % 位置调整圆周右移一位gr=circshift(gr',1)'; % 位置调整圆周右移一位sig1=ifft(fft(hr).*fft(sig1)); % 低频sig2=ifft(fft(gr).*fft(sig2)); % 高频sig=sig1+sig2; % 源信号%% 5.比较figure(7);subplot(2,1,1)plot(real(sig1));title('重构低频信号');subplot(2,1,2)plot(real(sig2));title('重构高频信号');figure(8);subplot(2,1,1)stem(abs(fft(sig1)));title('重构低频信号频谱');subplot(2,1,2)stem(abs(fft(sig2)));title('重构高频信号频谱');figure(9)plot(real(sig),'r','linewidth',2);hold on;plot(y);legend('重构信号','原始信号')title('重构信号与原始信号比较')小波包变换分析信号的MATLAB程序%t=0.001:0.001:1;t=1:1000;s1=sin(2*pi*50*t*0.001)+sin(2*pi*120*t*0.001)+rand(1,length(t));for t=1:500;s2(t)=sin(2*pi*50*t*0.001)+sin(2*pi*120*t*0.001)+rand(1,length(t));endfor t=501:1000;s2(t)=sin(2*pi*200*t*0.001)+sin(2*pi*120*t*0.001)+rand(1,length(t));endsubplot(9,2,1)plot(s1)title('原始信号')ylabel('S1')subplot(9,2,2)plot(s2)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(s130);ylabel('S130');subplot(9,2,5);plot(s131);ylabel('S131');subplot(9,2,7);plot(s132);ylabel('S132');subplot(9,2,9);plot(s133);ylabel('S133');subplot(9,2,11);plot(s134);ylabel('S134');subplot(9,2,13);plot(s135);ylabel('S135');subplot(9,2,15);plot(s136);ylabel('S136');subplot(9,2,17);plot(s137);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(s230);ylabel('S230');subplot(9,2,6);plot(s231);ylabel('S231');subplot(9,2,8);plot(s232);ylabel('S232');subplot(9,2,10);plot(s233);ylabel('S233');subplot(9,2,12);plot(s234);ylabel('S234');subplot(9,2,14);plot(s235);ylabel('S235');subplot(9,2,16);plot(s236);ylabel('S236');subplot(9,2,18);plot(s237);ylabel('S237');%fftfigurey1=fft(s1,1024);py1=y1.*conj(y1)/1024;y2=fft(s2,1024);py2=y2.*conj(y2)/1024;y130=fft(s130,1024);py130=y130.*conj(y130)/1024;y131=fft(s131,1024);py131=y131.*conj(y131)/1024;y132=fft(s132,1024);py132=y132.*conj(y132)/1024;y133=fft(s133,1024);py133=y133.*conj(y133)/1024;y134=fft(s134,1024);py134=y134.*conj(y134)/1024; y135=fft(s135,1024);py135=y135.*conj(y135)/1024; y136=fft(s136,1024);py136=y136.*conj(y136)/1024; y137=fft(s137,1024);py137=y137.*conj(y137)/1024;y230=fft(s230,1024);py230=y230.*conj(y230)/1024; y231=fft(s231,1024);py231=y231.*conj(y231)/1024; y232=fft(s232,1024);py232=y232.*conj(y232)/1024; y233=fft(s233,1024);py233=y233.*conj(y233)/1024; y234=fft(s234,1024);py234=y234.*conj(y234)/1024; y235=fft(s235,1024);py235=y235.*conj(y235)/1024; y236=fft(s236,1024);py236=y236.*conj(y236)/1024; y237=fft(s237,1024);py237=y237.*conj(y237)/1024;f=1000*(0:511)/1024; subplot(1,2,1);plot(f,py1(1:512));ylabel('P1');title('原始信号的功率谱') subplot(1,2,2);plot(f,py2(1:512));ylabel('P2');title('故障信号的功率谱') figuresubplot(4,2,1);plot(f,py130(1:512));ylabel('P130');title('S130的功率谱') subplot(4,2,2);plot(f,py131(1:512));title('S131的功率谱') subplot(4,2,3);plot(f,py132(1:512)); ylabel('P132'); subplot(4,2,4);plot(f,py133(1:512)); ylabel('P133'); subplot(4,2,5);plot(f,py134(1:512)); ylabel('P134'); subplot(4,2,6);plot(f,py135(1:512)); ylabel('P135'); subplot(4,2,7);plot(f,py136(1:512)); ylabel('P136'); subplot(4,2,8);plot(f,py137(1:512)); ylabel('P137'); figuresubplot(4,2,1);plot(f,py230(1:512)); ylabel('P230');title('S230的功率谱') subplot(4,2,2);plot(f,py231(1:512)); ylabel('P231');title('S231的功率谱') subplot(4,2,3);plot(f,py232(1:512)); ylabel('P232'); subplot(4,2,4);plot(f,py233(1:512)); ylabel('P233'); subplot(4,2,5);plot(f,py234(1:512)); ylabel('P234'); subplot(4,2,6);plot(f,py235(1:512)); ylabel('P235'); subplot(4,2,7);plot(f,py236(1:512));subplot(4,2,8);plot(f,py237(1:512));ylabel('P237');figure%plottree(wpt)利用小波变换实现对电能质量检测的算法实现N=10000;s=zeros(1,N);for n=1:Nif n<0.4*N||n>0.8*Ns(n)=31.1*sin(2*pi*50/10000*n);elses(n)=22.5*sin(2*pi*50/10000*n);endendl=length(s);[c,l]=wavedec(s,6,'db5'); %用db5小波分解信号到第六层subplot(8,1,1);plot(s);title('用db5小波分解六层：s=a6+d6+d5+d4+d3+d2+d1');Ylabel('s');%对分解结构【c，l】中第六层低频部分进行重构a6=wrcoef('a',c,l,'db5',6);subplot(8,1,2);plot(a6);Ylabel('a6');%对分解结构【c，l】中各层高频部分进行重构for i=1:6decmp=wrcoef('d',c,l,'db5',7-i);subplot(8,1,i+2);plot(decmp);Ylabel(['d',num2str(7-i)]);end%-----------------------------------------------------------rec=zeros(1,300);rect=zeros(1,300);ke=1;u=0;d1=wrcoef('d',c,l,'db5',1);figure(2);plot(d1);si=0;N1=0;N0=0;sce=0;for n=20:N-30rect(ke)=s(n);ke=ke+1;if(ke>=301)if(si==2)rec=rect;u=2;end;si=0;ke=1;end;if(d1(n)>0.01) % the condition of abnormal append.N1=n;if(N0==0)N0=n;si=si+1;end;if(N1>N0+30)Nlen=N1-N0;Tab=Nlen/10000;end;end;if(si==1)for k=N0:N0+99 %testing of 1/4 period signals to sce=sce+s(k)*s(k)/10000;end;re=sqrt(sce*200) %re indicate the pike value of .sce=0;si=si+1;end;end;NlenN0n=1:300;figure(3)plot(n,rec);基于小波变换的图象去噪Normalshrink算法function [T_img,Sub_T]=threshold_2_N(img,levels)% reference :image denoising using wavelet thresholding[xx,yy]=size(img);HH=img((xx/2+1):xx,(yy/2+1):yy);delt_2=(std(HH(:)))^2;%(median(abs(HH(:)))/0.6745)^2;%T_img=img;for i=1:levelstemp_x=xx/2^i;temp_y=yy/2^i;% belt=1.0*(log(temp_x/(2*levels)))^0.5;belt=1.0*(log(temp_x/(2*levels)))^0.5; %2.5 0.8%HLHL=img(1:temp_x,(temp_y+1):2*temp_y);delt_y=std(HL(:));T_1=belt*delt_2/delt_y;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%T_HL=sign(HL).*max(0,abs(HL)-T_1);%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%T_img(1:temp_x,(temp_y+1):2*temp_y)=T_HL;Sub_T(3*(i-1)+1)=T_1;%LHLH=img((temp_x+1):2*temp_x,1:temp_y);delt_y=std(LH(:));T_2=belt*delt_2/delt_y;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%T_LH=sign(LH).*max(0,abs(LH)-T_2);%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%T_img((temp_x+1):2*temp_x,1:temp_y)=T_LH;Sub_T(3*(i-1)+2)=T_2;%HHHH=img((temp_x+1):2*temp_x,(temp_y+1):2*temp_y);delt_y=std(HH(:));T_3=belt*delt_2/delt_y;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%T_HH=sign(HH).*max(0,abs(HH)-T_3);%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%T_img((temp_x+1):2*temp_x,(temp_y+1):2*temp_y)=T_HH;Sub_T(3*(i-1)+3)=T_3;end。