小波变换的程序

numberone00100

% a trous algorithmus (undecimated wavelet transform)

clear;clc;
%%  1.生成信号
f=50;   %  频率
fs=800; %  采样率
T=128;  %  信号长度
n=1:T;
y=sin(2*pi*f*n/fs)+2*exp(-f*n/(4*fs));  %  信号
% y=circshift(y.',3).';
%%  2.Transformation
l1=wfilters('db4','l')*sqrt(2)/2;     %  参考低通滤波器
l1_zeros=[l1,zeros(1,T-length(l1))];  %  低通滤波器1
h1=wfilters('db4','h')*sqrt(2)/2;     %  参考高通滤波器
h1_zeros=[h1,zeros(1,T-length(h1))];  %  高通滤波器1
low1=ifft(fft(y).*fft(l1_zeros));     %  低频分量1
high1=ifft(fft(y).*fft(h1_zeros));    %  高频分量1
l2=dyadup(l1);  %  原滤波器插值
l2_zeros=[l2,zeros(1,T-length(l2))];  %  低通滤波器2
h2=dyadup(h1);  %  原滤波器插值
h2_zeros=[h2,zeros(1,T-length(h2))];  %  高通滤波器2
low2=ifft(fft(low1).*fft(l2_zeros));  %  低频分量2
high2=ifft(fft(low1).*fft(h2_zeros)); %  高频分量2
%%  3.Ruektransformation
lr2=circshift(l2_zeros(end:-1:1).',1).';  %  重构低通滤波器2
hr2=circshift(h2_zeros(end:-1:1).',1).';  %  重构高通滤波器2
lr1=circshift(l1_zeros(end:-1:1).',1).';  %  重构低通滤波器1
hr1=circshift(h1_zeros(end:-1:1).',1).';  %  重构高通滤波器1
lowr=(ifft(fft(low2).*fft(lr2))+ifft(fft(high2).*fft(hr2)));  %  重构低频分量1(lowr=low1)
r_s=(ifft(fft(lowr).*fft(lr1))+ifft(fft(high1).*fft(hr1)));   %  重构源信号(r_s=y)

%%  4.绘图
figure(1);
plot(y);
title('源信号');
figure(2);
plot(low1,'r');
hold on;
plot(low2,'b');
legend('第一层低频','第二层低频');
figure(3);
plot(high1,'r');
hold on;
plot(high2,'b');
legend('第一层高频','第二层高频');
figure(4);
plot(low1,'r');
hold on;
plot(lowr,'b.');
legend('第一层低频','重构第一层低频');
figure(5);
plot(y,'r');
hold on;
plot(r_s,'b.');
legend('源信号','重构信号');
disp(norm(low1-lowr))
disp(norm(y-r_s))

[ 本帖最后由 eisenstange 于 2006-12-8 10:23 编辑 ]

熊猫羊

不错
这里的高手似乎很多
http://bbs.matwav.com/post/page?bid=5&sty=1&age=0