小孩害羞不自信怎么办:小波分析实验
实验1 连续小波变换
实验目的:
在理解连续小波变换原理的基础上,通过编程实现对一维信号进行连续小波变换,(实验中采用的是墨西哥帽小波),从而对连续小波变换增加了理性和感性的认识,并能提高编程能力,为今后的学习和工作奠定基础。
实验工具:
计算机,matlab6.5
程序附录:
(1) 墨西哥帽小波函数,按照(***)式编程
% mexh.m
function Y=mexh(x)
if abs(x)<=8
Y=exp(-x*x/2)*(1-x^2);
else
Y=0;
End
(2) 实验程序,按照(**)式编程,详细过程请参考“本实验采取的一些小技巧”
%
clc;clear;
load('data.mat');
len=length(dat);
lna=70; % (尺度a)的长度
a=zeros(1,lna);
wfab=zeros(lna,len); %小波系数矩阵
mexhab=zeros(1,len); % 离散化小波系数矩阵
for s=1:lna %s 表示尺度
for k=1:len
mexhab(k)=mexh(k/s);
end
for t=1:len % t 表示位移
wfab(s,t)=(sum(mexhab.*dat))/sqrt(s); %将积分用求和代替
mexhab=[mexh(-1*t/s),mexhab(1:len-1)]; %mexhab修改第一项并右移
end
end
figure(1);
plot(dat);
title('原始数据图');
figure(2); %小波系数谱
image(wfab);
colormap(pink(128));
title('小波系数图');
%surf(wfab);
%title('小波系数谱网格图');
%pwfab=wfab.*wfab; %%瞬态功率谱
%figure(3);
%subplot(1,2,1);
%surf(pwfab);
%title('瞬态功率谱网格图');
%subplot(1,2,2);
%contour(pwfab);
%title('瞬态功率谱等值线');
(3)test函数。
%test 函数
clc;clear;
for i=1:200
dat(i)=sin(2*pi*i*0.05); %正弦波函数
end
len=length(dat);
lna=40;
wfab=zeros(lna,len);
mexhab=zeros(1,len);
for s=1:lna %s 表示尺度
for k=1:len
mexhab(k)=mexh(k/s);
end
for t=1:len % t 表示位移
wfab(s,t)=(sum(mexhab.*dat))/sqrt(s); %将积分用求和代替
mexhab=[mexh(-1*t/s),mexhab(1:len-1)]; %mexhab修改第一项并右移
end
end
figure(1);
plot(dat);
title('orignal dat');
figure(2); %小波系数谱
image(wfab);
colormap(pink(128));
title('正弦波的小波系数图');
(4)用fft实现cwt
%按照圆周卷积定理,原周卷积和线性卷积的关系L>=M+N-1
%按照圆周卷积的定义,相关和线性卷积的关系(原始算法和线性卷积的关系)
%注意画图理解
clc;clear;
t1=cputime;
load('data.mat');
len=length(dat);
lna=70; % a(尺度)的长度
a=zeros(1,lna); % a 表示尺度
b=zeros(1,len); % b 表示位移
wfab=zeros(lna,len); %小波系数矩阵
mexhab=zeros(1,2*len-1);
data=[zeros(1,len-1),dat];
Ydata=fft( data ,4*len);
for s=1:lna
for k=1:2*len-1
mexhab(k)=mexh((k-len)/s);
end
temp=ifft( Ydata.*fft( mexhab,4*len ) ,4*len);
wfab(s,:)=real(temp(2*len-1:3*len-2))/sqrt(s); %为什么要取实部而不是取模,我也不是很清楚,可是有种感觉
end
figure(1);
plot(dat);
title('原始数据图');
figure(2); %小波系数谱
image(wfab);
colormap(pink(128));
title('小波系数谱 ');
cputime-t1
4)fft快速计算cwt
%按照圆周卷积的定义,
%注意画图理解
clc;clear;
t1=cputime;
load('data.mat');
len=length(dat);
lna=70; % a(尺度)的长度
a=5;
data=[dat,zeros(1,len)];
Ydata=fft(dat,2*len);
for s=1:lna
mexhab=zeros(1,2*len);
k=[-a*s:1:a*s];
mexhab(k+len)=mexh2(k./s);
temp=ifft( Ydata.*fft( mexhab,2*len ) ,2*len);
wfab(s,:)=real(temp(len+1:2*len))/sqrt(s); %要取实部而不是取模,呵呵
end
figure(1);
plot(dat);
title('原始数据图');
figure(2); %小波系数谱
image(wfab);
colormap(pink(128));
title('小波系数谱 ');
cputime-t1
5)保存为mexh2.m
function Y=mexh2(x)
Y=exp(-x.*x/2).*(1-x.^2);
Torstan
2005.09.16
实验2 二维离散小波变换(Mallat快速算法)
实验目的:
在理解离散小波变换原理和Mallat快速算法的基础上,通过编程对图像进行二维离散小波变换,从而加深对二维小波分解和重构的理性和感性认识,并能提高编程能力,为今后的学习和工作奠定基础。
实验工具:
计算机,matlab6.5
附录:
(1)二维小波分解函数
%二维小波分解函数
function Y=mallatdec2(X,wname,level)
%输入:X 载入的二维图像像数值;
% level 小波分解次(级)数设定值(如果设定值超过最高可分解次数,按最高分解次数分解)
% wname 小波名字wavelet name
%输出:Y 多极小波分解后的小波系数矩阵
[h,g]=wfilters(wname,'d'); %h,g分别为低通和高通滤波器
X=double(X);
t=1;
hh=size(X,2);
while t<=level
%先进行行小波变换
for row=1:hh
Y(row,1:hh)=mdec1(X(row,1:hh),h,g) ;
end
%再进行列小波变换
for col=1:hh
temp=mdec1( Y(1:hh,col)',h,g);
Y(1:hh,col)=temp';
end
t=t+1;
hh=hh/2;
X=Y;
end
%内部子函数,对一行(row)矢量进行一次小波变换,利用fft实现
function y=mdec1(x,h,g)
%输入:x 行数组
% h为低通滤波器
% g为高通滤波器
%输出: y 进行一级小波分解后的系数
lenx=size(x,2);
lenh=size(h,2);
rh=h(end:-1:1);
rrh=[zeros(1,(lenx-lenh)),rh];
rrh=circshift(rrh',1)';
rg=g(end:-1:1);
rrg=[zeros(1,(lenx-lenh)),rg];
rrg=circshift(rrg',1)';
r1=dyaddown(ifft(fft(x).*fft(rrh,lenx)),1); %use para 1
r2=dyaddown(ifft(fft(x).*fft(rrg,lenx)),1);
y=[r1,r2];
(2)二维小波重构函数
%二维小波重构函数
function Y=mallatrec2(X,wname,level)
%输入:X 载入的小波系数矩阵;
% level 小波分解次(级)数设定值(如果设定值超过最高可分解次数,按最高分解次数分解)
% wname 小波名字wavelet name
%输出:Y 重构图像矩阵
[h,g]=wfilters(wname,'d'); %h,g分别为重构低通滤波器和重构高通滤波器
hz=size(X,2);
h1=hz/(2^(level-1));
while h1<=hz
% 对列变换
for col=1:h1
temp=mrec1(X(1:h1,col)',h,g)';
X(1:h1,col)=temp;
end
%再对行变换
for row=1:h1
temp=mrec1(X(row,1:h1),h,g);
X(row,1:h1)=temp;
end
h1=h1*2;
end
Y=X;
%内部子函数,对一行小波系数进行重构
function y=mrec1(x,h,g)
%输入:x 行数组
% h为低通滤波器
% g为高通滤波器
%输出: y 进行一级小波重构后值
lenx=size(x,2);
r3=dyadup(x(1,1:lenx*0.5),0); %内插零use para 0
r4=dyadup(x(1,(lenx*0.5+1):lenx),0); %use para 0
y=ifft(fft(r3,lenx).*fft(h,lenx))+ ifft(fft(r4,lenx).*fft(g,lenx));