MATLAB2维小波变更典型挨次
% FWT_DB.M;
% 此展现挨次用DWT实现二维小波变更
% 编程光阴2004-4-10，编程人沙威
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;clc;
T=256; % 图像维数
SUB_T=T/2; % 子图维数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 1.调原始图像矩阵
load wbarb; % 下载图像
f=X; % 原始图像
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 2.妨碍二维小波分解
l=wfilters('db10','l'); % db10（消逝矩为10)低通分解滤波器侵略照应（长度为20）
L=T-length(l);
l_zeros=[l,zeros(1,L)]; % 矩阵行数与输入图像不同，为2的整数幂
h=wfilters('db10','h'); % db10（消逝矩为10)高通分解滤波器侵略照应（长度为20）
h_zeros=[h,zeros(1,L)]; % 矩阵行数与输入图像不同，为2的整数幂
for i=1:T; % 列变更
row(1:SUB_T,i)=dyaddown( ifft( fft(l_zeros).*fft(f(:,i)') ) ).'; % 圆周卷积<->FFT
row(SUB_T+1:T,i)=dyaddown( ifft( fft(h_zeros).*fft(f(:,i)') ) ).'; % 圆周卷积<->FFT
end;
for j=1:T; % 行变更
line(j,1:SUB_T)=dyaddown( ifft( fft(l_zeros).*fft(row(j,:)) ) ); % 圆周卷积<->FFT
line(j,SUB_T+1:T)=dyaddown( ifft( fft(h_zeros).*fft(row(j,:)) ) ); % 圆周卷积<->FFT
end;
decompose_pic=line; % 分解矩阵
% 图像分为四块
lt_pic=decompose_pic(1:SUB_T,1:SUB_T); % 在矩阵左上方为低频份量--fi(x)*fi(y)
rt_pic=decompose_pic(1:SUB_T,SUB_T+1:T); % 矩阵右上为--fi(x)*psi(y)
lb_pic=decompose_pic(SUB_T+1:T,1:SUB_T); % 矩阵左下为--psi(x)*fi(y)
rb_pic=decompose_pic(SUB_T+1:T,SUB_T+1:T); % 右下方为高频份量--psi(x)*psi(y)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 3.分解服从展现
figure(1);
colormap(map);
subplot(2,1,1);
image(f); % 原始图像
title('original pic');
subplot(2,1,2);
image(abs(decompose_pic)); % 分解后图像
title('decomposed pic');
figure(2);
colormap(map);
subplot(2,2,1);
image(abs(lt_pic)); % 左上方为低频份量--fi(x)*fi(y)
title('Phi(x)*Phi(y)');
subplot(2,2,2);
image(abs(rt_pic)); % 矩阵右上为--fi(x)*psi(y)
title('Phi(x)*Psi(y)');
subplot(2,2,3);
image(abs(lb_pic)); % 矩阵左下为--psi(x)*fi(y)
title('Psi(x)*Phi(y)');
subplot(2,2,4);
image(abs(rb_pic)); % 右下方为高频份量--psi(x)*psi(y)
title('Psi(x)*Psi(y)');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 5.重构源图像及服从展现
% construct_pic=decompose_matrix'*decompose_pic*decompose_matrix;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
l_re=l_zeros(end:-1:1); % 重构低通滤波
l_r=circshift(l_re',1)'; % 位置调解
h_re=h_zeros(end:-1:1); % 重构高通滤波
h_r=circshift(h_re',1)'; % 位置调解
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
top_pic=[lt_pic,rt_pic]; % 图像上半部份
t=0;
for i=1:T; % 行插值低频
if (mod(i,2)==0)
topll(i,:)=top_pic(t,:); % 偶数行坚持
else
t=t+1;
topll(i,:)=zeros(1,T); % 奇数行动零
end
end;
for i=1:T; % 列变更
topcl_re(:,i)=ifft( fft(l_r).*fft(topll(:,i)') )'; % 圆周卷积<->FFT
end;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
bottom_pic=[lb_pic,rb_pic]; % 图像下半部份
t=0;
for i=1:T; % 行插值高频
if (mod(i,2)==0)
bottomlh(i,:)=bottom_pic(t,:); % 偶数行坚持
else
bottomlh(i,:)=zeros(1,T); % 奇数行动零
t=t+1;
end
end;
for i=1:T; % 列变更
bottomch_re(:,i)=ifft( fft(h_r).*fft(bottomlh(:,i)') )'; % 圆周卷积<->FFT
end;
construct1=bottomch_re+topcl_re; % 列变更重构竣事
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
left_pic=construct1(:,1:SUB_T); % 图像左半部份
t=0;
for i=1:T; % 列插值低频
if (mod(i,2)==0)
leftll(:,i)=left_pic(:,t); % 偶数列坚持
else
t=t+1;
leftll(:,i)=zeros(T,1); % 奇数列为零
end
end;
for i=1:T; % 行变更
leftcl_re(i,:)=ifft( fft(l_r).*fft(leftll(i,:)) ); % 圆周卷积<->FFT
end;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
right_pic=construct1(:,SUB_T+1:T); % 图像右半部份
t=0;
for i=1:T; % 列插值高频
if (mod(i,2)==0)
rightlh(:,i)=right_pic(:,t); % 偶数列坚持
else
rightlh(:,i)=zeros(T,1); % 奇数列为零
t=t+1;
end
end;
for i=1:T; % 行变更
rightch_re(i,:)=ifft( fft(h_r).*fft(rightlh(i,:)) ); % 圆周卷积<->FFT
end;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
construct_pic=rightch_re+leftcl_re; % 重修全副图像
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 服从展现
figure(3);
colormap(map);
subplot(2,1,1);
image(f); % 源图像展现
title('original pic');
subplot(2,1,2);
image(abs(construct_pic)); % 重构源图像展现
title('reconstructed pic');
error=abs(construct_pic-f); % 重构图形与原始图像误值
figure(4);
mesh(error); % 倾向三维图像
title('absolute error display');
clear
clc
%在噪声情景下语音信号的增强
%语音信号为读入的声音文件
%噪声为正态随机噪声
sound=wavread('c12345.wav');
count1=length(sound);
noise=0.05*randn(1,count1);
for i=1:count1
signal(i)=sound(i);
end
for i=1:count1
y(i)=signal(i)+noise(i);
end
%在小波基'db3'下妨碍一维离散小波变更
[coefs1,coefs2]=dwt(y,'db3'); %[低频 高频]
count2=length(coefs1);
count3=length(coefs2);
energy1=sum((abs(coefs1)).^2);
energy2=sum((abs(coefs2)).^2);
energy3=energy1+energy2;
for i=1:count2
recoefs1(i)=coefs1(i)/energy3;
end
for i=1:count3
recoefs2(i)=coefs2(i)/energy3;
end
%低频系数妨碍语音信号清浊音的分说
zhen=160;
count4=fix(count2/zhen);
for i=1:count4
n=160*(i-1)+1:160+160*(i-1);
s=sound(n);
w=ha妹妹ing(160);
sw=s.*w;
a=aryule(sw,10);
sw=filter(a,1,sw);
sw=sw/sum(sw);
r=xcorr(sw,'biased');
corr=max(r);
%为浊音（unvoice）时，输入为1；为浊音（voice）时，输入为0
if corr>=0.8
output1(i)=0;
elseif corr<=0.1
output1(i)=1;
end
end
for i=1:count4
n=160*(i-1)+1:160+160*(i-1);
if output1(i)==1
switch abs(recoefs1(i))
case abs(recoefs1(i))<=0.002
recoefs1(i)=0;
case abs(recoefs1(i))>0.002 & abs(recoefs1(i))<=0.003
recoefs1(i)=sgn(recoefs1(i))*(0.003*abs(recoefs1(i))-0.000003)/0.002;
otherwise recoefs1(i)=recoefs1(i);
end
elseif output1(i)==0
recoefs1(i)=recoefs1(i);
end
end
%对于高频系数妨碍语音信号清浊音的分说
count5=fix(count3/zhen);
for i=1:count5
n=160*(i-1)+1:160+160*(i-1);
s=sound(n);
w=ha妹妹ing(160);
sw=s.*w;
a=aryule(sw,10);
sw=filter(a,1,sw);
sw=sw/sum(sw);
r=xcorr(sw,'biased');
corr=max(r);
%为浊音（unvoice）时，输入为1；为浊音（voice）时，输入为0
if corr>=0.8
output2(i)=0;
elseif corr<=0.1
output2(i)=1;
end
end
for i=1:count5
n=160*(i-1)+1:160+160*(i-1);
if output2(i)==1
switch abs(recoefs2(i))
case abs(recoefs2(i))<=0.002
recoefs2(i)=0;
case abs(recoefs2(i))>0.002 & abs(recoefs2(i))<=0.003
recoefs2(i)=sgn(recoefs2(i))*(0.003*abs(recoefs2(i))-0.000003)/0.002;
otherwise recoefs2(i)=recoefs2(i);
end
elseif output2(i)==0
recoefs2(i)=recoefs2(i);
end
end
%在小波基'db3'下妨碍一维离散小波反变更
output3=idwt(recoefs1, recoefs2,'db3');
%对于输入信号抽样点值妨碍归一化处置
maxdata=max(output3);
output4=output3/maxdata;
%读出带噪语音信号，存为'101.wav'
wavwrite(y,5500,16,'c101');
%读出处置后语音信号，存为'102.wav'
wavwrite(output4,5500,16,'c102');
function [I_W , S] = func_DWT(I, level, Lo_D, Hi_D);
%经由这个函数将I妨碍小波分解，并将分解后的一维向量转换为矩阵方式
% Matlab implementation of SPIHT (without Arithmatic coding stage)
% Wavelet decomposition
% input: I : input image
% level : wavelet decomposition level
% Lo_D : low-pass decomposition filter
% Hi_D : high-pass decomposition filter
% output: I_W : decomposed image vector
% S : corresponding bookkeeping matrix
% please refer wavedec2 function to see more
[C,S] = func_Mywavedec2(I,level,Lo_D,Hi_D);
S(:,3) = S(:,1).*S(:,2); % dim of detail coef nmatrices 求低频以及每一个尺度中高频的元素个数
%st=S(1,3)+S(2,3)*3+S(3,3)*3;%%%%对于前两层加密
%C(1:st)=0;
L = length(S); %a求S的列数
I_W = zeros(S(L,1),S(L,2));%设一个与原图像巨细相同的全零矩阵
% approx part
I_W( 1:S(1,1) , 1:S(1,2) ) = reshape(C(1:S(1,3)),S(1,1:2)); %将LL层从C中复原为S(1,1)*S(1,2)的矩阵
for k = 2 : L-1 %将C向量中复原出HL,HH,LH 矩阵
rows = [sum(S(1:k-1,1))+1:sum(S(1:k,1))];
columns = [sum(S(1:k-1,2))+1:sum(S(1:k,2))];
% horizontal part
c_start = S(1,3) + 3*sum(S(2:k-1,3)) + 1;
c_stop = S(1,3) + 3*sum(S(2:k-1,3)) + S(k,3);
I_W( 1:S(k,1) , columns ) = reshape( C(c_start:c_stop) , S(k,1:2) );
% vertical part
c_start = S(1,3) + 3*sum(S(2:k-1,3)) + S(k,3) + 1;
c_stop = S(1,3) + 3*sum(S(2:k-1,3)) + 2*S(k,3);
I_W( rows , 1:S(k,2) ) = reshape( C(c_start:c_stop) , S(k,1:2) );
% diagonal part
c_start = S(1,3) + 3*sum(S(2:k-1,3)) + 2*S(k,3) + 1;
c_stop = S(1,3) + 3*sum(S(2:k,3));
I_W( rows , columns ) = reshape( C(c_start:c_stop) , S(k,1:2) );
end
%%%%%%%mallat algorithm%%%%% clc; clear;tic; %%%%original signal%%%% f=100;%%frequence ts=1/800;%%抽样距离 N=1:100;%%点数 s=sin(2*ts*pi*f.*N);%%源信号 figure(1) plot(s);%%%源信号s title('原信号'); grid on; %%%%小波滤波器%%%% ld=wfilters('db1','l');%%低通 hd=wfilters('db1','h');%%高通 figure(2) stem(ld,'r');%%%低通 grid on; figure(3) stem(hd,'b')%%%高通 grid on; %%%%% tem=conv(s,ld);%%低通以及原信号卷积 ca1=dyaddown(tem);%%抽样 figure(4) plot(ca1); grid on; tem=conv(s,hd);%%高通以及原信号卷积 cb1=dyaddown(tem);%%抽样 figure(5) plot(cb1); grid on; %%%%%%%% %[ca3,cb3]=dwt(s,'db1');%%小波变更 %%%%%%%% [lr,hr]=wfilters('db1','r');%%重构滤波器 figure(6) stem(lr); figure(7) stem(hr); tem=dyadup(cb1);%%插值 tem=conv(tem,hr);%%卷积 d1=wkeep(tem,100);%%去掉中间的份量 %%%%%%%%% tem=dyadup(ca1);%%插值 tem=conv(tem,lr);%%卷积 a1=wkeep(tem,100);%%去掉中间的份量 a=a1+d1;%%%重构原信号 %%%%%%%%% %a3=idwt(ca3,cb3,'db1',100);%%%小波逆变更 %%%%%%%%% figure(8) plot(a,'.b'); hold on; plot(s,'r'); grid on; title('重构信号以及原信号的比力');toc; %figure(9) %plot(a3,'.b'); %hold on; %plot(s,'r'); %grid on; %title('重构信号以及原信号的比力');
通用函数
Allnodes 合计树结点 appcoef 提取一维小波变更低频系数 appcoef2 提取二维小波分解低频系数 bestlevt 合计残缺最佳小波包树 besttree 合计最佳(优)树 * biorfilt 双正交样条小波滤波器组 biorwavf 双正交样条小波滤波器 * centfrq 求小波中间频率 cgauwavf Complex Gaussian小波 cmorwavf coiflets小波滤波器 cwt 一维不断小波变更 dbaux Daubechies小波滤波器合计 dbwavf Daubechies小波滤波器 dbwavf(W) W='dbN' N=1,2,3,...,50 ddencmp 取患上默认值阈值(软或者硬)熵尺度 depo2ind 将深度-位置结点方式转化成索引结点方式 detcoef 提取一维小波变更高频系数 detcoef2 提取二维小波分解高频系数 disp 展现文本或者矩阵 drawtree 画小波包分解树(GUI) dtree 妄想DTREE类 dwt 单尺度一维离散小波变更 dwt2 单尺度二维离散小波变更 dwtmode 离散小波变更拓展方式 * dyaddown 二元取样 * dyadup 二元插值 entrupd 更新小波包的熵值 fbspwavf B样条小波 gauswavf Gaussian小波 get 取患上工具属性值 idwt 单尺度一维离散小波逆变更 idwt2 单尺度二维离散小波逆变更 ind2depo 将索引结点方式转化成深度—位置结点方式 * intwave 积分小波数 isnode 分说结点是否存在 istnode 分说结点是否是开幕点并返回部署值 iswt 一维逆SWT(Stationary Wavelet Transform)变更 iswt2 二维逆SWT变更 leaves Determine terminal nodes mexihat 墨西哥帽小波 meyer Meyer小波 meyeraux Meyer小波辅助函数 morlet Morlet小波 nodease 合计上溯结点 nodedesc 合计下溯结点(子结点) nodejoin 重组结点 nodepar 追寻父结点 nodesplt 分割(分解)结点 noleaves Determine nonterminal nodes ntnode Number of terminal nodes ntree Constructor for the class NTREE * orthfilt 正交小波滤波器组 plot 绘制向量或者矩阵的图形 * qmf 镜像二次滤波器 rbiowavf Reverse biorthogonal spline wavelet filters read 读取二进制数据 readtree 读取小波包分解树 * scal2frq Scale to frequency set shanwavf Shannon wavelets swt 一维SWT(Stationary Wavelet Transform)变更 swt2 二维SWT变更 symaux Symlet wavelet filter computation. symwavf Symlets小波滤波器 thselect 信号消噪的阈值抉择 thodes References treedpth 求树的深度 treeord 求树妄想的叉数 upcoef 一维小波分解系数的直接重构 upcoef2 二维小波分解系数的直接重构 upwlev 单尺度一维小波分解的重构 upwlev2 单尺度二维小波分解的重构 wavedec 单尺度一维小波分解 wavedec2 多尺度二维小波分解 wavedemo 小波工具箱函数demo * wavefun 小波函数以及尺度函数 * wavefun2 二维小波函数以及尺度函数 wavemenu 小波工具箱函数menu图形界面调用函数 * wavemngr 小波规画函数 waverec 多尺度一维小波重构 waverec2 多尺度二维小波重构 wbmpen Penalized threshold for wavelet 1-D or 2-D de-noising wcodemat 对于矩阵妨碍量化编码 wdcbm Thresholds for wavelet 1-D using Birge-Massart strategy wdcbm2 Thresholds for wavelet 2-D using Birge-Massart strategy wden 用小波妨碍一维信号的消噪或者缩短 wdencmp De-noising or compression using wavelets wentropy 合计小波包的熵 wextend Extend a vector or a matrix * wfilters 小波滤波器 wkeep 提取向量或者矩阵中的一部份 * wmaxlev 合计小波分解的最大尺度 wnoise 发生含噪声的测试函数数据 wnoisest 估量一维小波的系数的尺度倾向 wp2wtree 从小波包树中提取小波树 wpcoef 合计小波包系数 wpcutree 剪切小波包分解树 wpdec 一维小波包的分解 wpdec2 二维小波包的分解 wpdencmp 用小波包妨碍信号的消噪或者缩短 wpfun 小波包函数 wpjoin 重组小波包 wprcoef 小波包分解系数的重构 wprec 一维小波包分解的重构 wprec2 二维小波包分解的重构 wpsplt 分割（分解）小波包 wpthcoef 妨碍小波包分解系数的阈值处置 wptree 展现小波包树妄想 wpviewcf Plot the colored wavelet packet coefficients. wrcoef 对于一维小波系数妨碍单支重构 wrcoef2 对于二维小波系数妨碍单支重构 wrev 向量逆序 write 向缓冲区内存写进数据 wtbo Constructor for the class WTBO wthcoef 一维信号的小波系数阈值处置 wthcoef2 二维信号的小波系数阈值处置 wthresh 妨碍软阈值或者硬阈值处置 wthrmngr 阈值配置规画 wtreemgr 规画树妄想
;)
