消失矩作用的程序
clear;clc;f=50;
T=0.001;
n=1:50;
y=sin(2*pi*f*n*T);
noise=;
y_noise=y+noise;
figure(1);
subplot(2,1,1);
plot(y);
title('signal');
subplot(2,1,2);
plot(y_noise);
title('noise & signal');
=dwt(y,'db2');
=dwt(y,'db10');
figure(2);
subplot(2,1,1);
plot(yl2);
title('db2 low frequency signal');
subplot(2,1,2);
plot(yh2);
title('db2 high frequency signal');
figure(3);
subplot(2,1,1);
plot(yl10);
title('db10 low frequency signal');
subplot(2,1,2);
plot(yh10);
title('db10 high frequency signal');
=dwt(y_noise,'db2');
=dwt(y_noise,'db10');
figure(4);
subplot(2,1,1);
plot(yl2);
title('db2 low frequency noise & signal');
subplot(2,1,2);
plot(yh2);
title('db2 high frequency noise & signal');
figure(5);
subplot(2,1,1);
plot(yl10);
title('db10 low frequency noise & signal');
subplot(2,1,2);
plot(yh10);
title('db10 high frequency noise & signal');
[ 本帖最后由 simon21 于 2007-4-4 07:08 编辑 ]
平移变换平移法(cycle_spinning)消除gibbs效应
clear;clc; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 1.原始信号 f=50; % 信号频率
fs=800; % 采样频率
N=128; % 采样点 % 信号赋值
n=1:N;
y=sin(2*pi*f*n/fs); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 2.噪声
noise=0.4*rand(1,128); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 3.染噪信号
y_noise=y+noise; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 4.硬消噪采用cycle_spinning技术 % 累加量
z5=zeros(1,N); % 平移变换频移法
for i=1:N;
z=circshift(y_noise.',i-1).'; % 源信号右平移
=lwt(z,'db3'); % 小波正变换
z2=zeros(1,N/2); % 高频分量全部为零(主要噪声,硬消噪)
z3=ilwt(z1,z2,'db3'); % 小波反变换
z4=circshift(z3.',-(i-1)).'; % 变换后信号左平移
z5=z5+z4/N; % 平均
end;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 5.显示 error=norm(y-z5)/norm(y); % 相对误差 figure(1);
subplot(2,1,1)
plot(y,'r');
legend('源信号');
subplot(2,1,2);
plot(y_noise);
legend('染噪信号');
figure(2);
subplot(2,1,1)
plot(y,'r');
legend('源信号');
title(error);
subplot(2,1,2);
plot(z5);
legend('消噪后信号');
[ 本帖最后由 simon21 于 2007-4-4 07:08 编辑 ]
请教一个问题,关于三楼那个程序(2代小波示意程序 )
本帖最后由 wdhd 于 2016-9-6 13:39 编辑2代小波示意程序 (三楼)中
……
f1(:,T+1)=f1(:,1); % 补列
f2(T+1,:)=f2(1,:); % 补行
……
补行是把分裂的偶数的第一行补到偶数的最后,但是补列是怎么进行的?比如图像是256×256,分裂后奇数为128×256,偶数也为128×256,进行f1(:,T+1)=f1(:,1); 补列之后,相当于把原来的第129列替换成第1列这样补行和补列后,奇数和偶数的分裂结果为128×256和129×256还请大侠不吝赐教,小生在此谢过了
[此贴子已经被作者于2006-5-24 16:09:37编辑过]
多谢楼主,都是好东西~
真的不错
其实我一直在努力,但是水平有限,还是先向各位学习,谢谢大家! 很不错。。。是需要有人来总结和整理。 虽然不搞小波,但是为楼主的这种无私所敬佩! 楼主真是太好了!!狂赞!!博学多才! 楼主真是大好人!能加楼主的QQ么?我的QQ573618825!谢谢楼主!:handshake 我也是研究小波的,由于我身边没有人和我研究的一样,有好多问题想咨询楼主,在这里先谢谢楼主了。:@D
信号的分析处理
matlab内容太多了,我是想把从实际中测得的信号,用matlab的分析功能进行分析,各位高手请指点,多谢谢! 多谢楼主的无私奉献,在此感谢了/ 顶一下 使用小波包变换分析信号的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));
end
for t=501:1000;
s2(t)=sin(2*pi*200*t*0.001)+sin(2*pi*120*t*0.001)+rand(1,length(t));
end
subplot(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,);
s131=wprcoef(wpt,);
s132=wprcoef(wpt,);
s133=wprcoef(wpt,);
s134=wprcoef(wpt,);
s135=wprcoef(wpt,);
s136=wprcoef(wpt,);
s137=wprcoef(wpt,);
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=
std1=
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,);
s231=wprcoef(wpt,);
s232=wprcoef(wpt,);
s233=wprcoef(wpt,);
s234=wprcoef(wpt,);
s235=wprcoef(wpt,);
s236=wprcoef(wpt,);
s237=wprcoef(wpt,);
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=
std2=
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');
%fft
figure
y1=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('故障信号的功率谱')
figure
subplot(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');
figure
subplot(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)
基于小波消噪的雷达回波检测
基于小波消噪的雷达回波检测,可以检测雷达回波的有无及其准确的位置close all;clear;clc;
signalw=1200;%信号宽度
m=120000; %数据采样点数
pfa=0.001;%虚警概率
dectectw=signalw/4;%滑窗宽度
a=wgn(1,m,13,50,'dbm','real');%噪声
b=gate(a,m,dectectw,pfa)%门限
signal=real(fmlin(signalw,0,0.25));%信号
%加高斯白噪声信号
mix=zeros(1,m);
signalpos=40000;
for l=signalpos:signalpos+signalw-1
mix(1,l)=signal(l-(signalpos-1),1);
end
mix=mix+a;%加高斯白噪声信号
%检测信号有无及位置
g=0;
for l=1:dectectw
g=g+mix(l)*(mix(l))';
end
g0=g;
for k=(dectectw+1):m
g=g+mix(k)*(mix(k))'-mix(k-dectectw)*(mix(k-dectectw))';
if (g>b)
over=1
index=k
break
end
end%检测信号有无及位置
pnoise=mean(abs(a).^2);
psignal=mean(abs(signal).^2);
snr=10*log10(psignal/pnoise) %信噪比
=wden(mix,'rigrsure','s','sln',3,'db2');
=wden(a,'rigrsure','s','sln',3,'db2');
% figure(2)
% subplot(2,1,1)
% plot(a)
% subplot(2,1,2)
% plot(a1)
=wden(signal,'rigrsure','s','sln',3,'db2');
% figure(3)
% subplot(2,1,1)
% plot(signal)
% subplot(2,1,2)
% plot(signal1)
%检测信号有无及位置
dectectw1=signalw*1/4;
b1=gate(a1,m,dectectw1,pfa)
g=0;
for l=1:dectectw1
g=g+mix1(l)*(mix1(l))';
end
g0=g;
for k=(dectectw1+1):m
g=g+mix1(k)*(mix1(k))'-mix1(k-dectectw1)*(mix1(k-dectectw1))';
if (g>b1)
over=1
index=k
break
end
end%检测信号有无及位置
pnoise1=(sum(abs(mix1).^2)-sum(abs(signal1).^2))/m;
psignal1=mean(abs(signal1).^2);
snr1=10*log10(psignal1/pnoise1) %信噪比
subplot(2,1,1)
plot(mix)
axis()
TITLE('消噪前信号')
str={['SNR=',num2str(snr),'dB']};
text(45000,3,str)
subplot(2,1,2)
plot(mix1)
axis()
TITLE('消噪后信号')
str={['SNR=',num2str(snr1),'dB']};
text(45000,2.5,str)