搜索
小波滤波器构造和消噪程序（2个）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)，低通分解滤波器fori_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_wavelet.m%此函数用于研究Mallet算法及滤波器设计%此函数用于消噪处理%角度赋值%此处赋值使滤波器系数恰为db9%分解的高频系数采用db9较好,即它的消失矩较大%分解的有用信号小波高频系数基本趋于零%对于噪声信号高频分解系数很大，便于阈值消噪处理[l,h]=wfilters('db10','d');low_construct=l;L_fre=20;%滤波器长度low_decompose=low_construct(end:-1:1);%确定h0(-n)，低通分解滤波器fori_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('小波分解的高频系数');%消噪处理fori_decrease=31:44;ifabs(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;%消噪后信号重构%平滑处理forj=1:2fori=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);holdon;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));fort=1:500;s2(t)=sin(2*pi*50*t*0.001)+sin(2*pi*120*t*0.001)+rand(1,length(t));endfort=501:1000;s2(t)=sin(2*pi*200*t*0.001)+sin(2*pi*120*t*0.001)+rand(1,length(t));endsu