5月4日架接别人的,已调通,,,,,,,,用于学习小波算法,,,,,,,,
但是对于小波算法中滤 波器的选取仍然不太清楚,,,,,,呵
function coef=mywavedec2(x,N,wname)
clc;clear;
%button=questdlg('questdlg','确定'); %打开问题提示对话框,有三个按钮,分别为:yes,no和cancel,’questdlg’确定提示信息。
pic_in_wavelet = imread('123.jpg');%%读取图像
figure(1),imshow(pic_in_wavelet);
[m,n,p]=size(pic_in_wavelet);
for i=1:m
for j=1:n
gray_pic(i,j)=0.3*pic_in_wavelet(i,j,1)+0.59*pic_in_wavelet(i,j,2)+0.11*pic_in_wavelet(i,j,3);
end
end
figure(2),imshow(gray_pic);
x=gray_pic;
N=3;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 函数 MYWAVEDEC2() 对输入矩阵 x 进行 dim 层分解,得到相应的分解系数矩阵 y
% 输入参数:x —— 输入矩阵
% N —— 分解级数
% wname —— 分解所用的小波函数
% 输出参数:coef —— 分解系数矩阵,其结构如下:
% coef = {cA_N;cV_N;cH_N;cD_N;cV_N-1;cH_N-1;cD_N-1;……;cV_1;cH_1;cD_1}
% Copyright by Zou Yuhua ( chenyusiyuan ), original : 2007-11-10, modified: 2008-06-04
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 求出小波函数的滤波器组系数向量
[Lo_D,Hi_D] = wfilters('bior3.3','d');
% 画出原始图像
imshow(x);title('Original Image');
% 标明图像大小
[r,c]=size(x);
xlabel(['Size : ',num2str(r),'*',num2str(c)]);
% 将矩阵x的数据格式转换为适合数值处理的double格式
xd=double(x);
coef=[];
for i=1:N
[cA,cV,cH,cD]=mydwt2(xd,Lo_D,Hi_D);
% 第 i 级小波分解
xd=cA; % 将第 i 级分解得到的低频系数矩阵作为第 i+1 级分解的源矩阵
outmp={cV;cH;cD}; % 将第 i 级分解得到的高频系数矩阵cV,cH,cD存入细胞矩阵 outmp
% 注意细胞矩阵的赋值是用大括号“{}”的,而普通矩阵赋值是用方括号“[]”
% 细胞矩阵不要求其中的子矩阵的行列数都相同
coef=[outmp;coef]; % 将细胞矩阵 outmp 存入输出矩阵 coef,coef将由空矩阵变为细胞矩阵
% 注意这里的方括号不能用大括号取代
% 否则,使用大括号会将初始的coef空矩阵也作为细胞矩阵的子矩阵
% 而且,在迭代中 coef 将是一个不断嵌套的细胞矩阵,不便于后续处理和读取
% 上面这个语句是一种有效的在迭代过程中保存数据的方法
% 设待存数据为 data,可以是单个数、向量或矩阵
% 保存数据的矩阵为 mat,初始为空矩阵:mat=[]
% 则可按以下格式保存迭代过程产生的数据
% mat=[mat;data];
% 方括号内的分号“;”表示数据 data 是按“列”排序的方式存入矩阵 mat
% mat=[mat,data];
% 方括号内的逗号“,”表示数据 data 是按“行”排序的方式存入矩阵 mat
% data 也可以在 mat 前嵌入,即 mat=[data;mat] 或 mat=[data,mat]
end
% 迭代结束后,矩阵 coef 中保存的是各级分解中的高频系数矩阵
% 故需将迭代后得到的矩阵 cA,即第 dim 级低频矩阵存入矩阵 coef
coef=[cA;coef];
% {zh1},小波系数矩阵 coef 的结构如下
% coef = {cA_N;cV_N;cH_N;cD_N;cV_N-1;cH_N-1;cD_N-1;……;cV_1;cH_1;cD_1}
% 画出各级低频、高频系数矩阵
% 首先建立一个名为“Wavelet Decomposition -- Wavelet Type: , Levels: ”的图像窗口
%figure('Name',['Wavelet Decomposition -- Wavelet Type: ',wname,' , Levels: ',num2str(N)]);
% 图像的第1行显示低频系数,置中,左右两个subplot为空
subplot(N+1,3,2);
yt=uint8(coef{1});
[yrow,ycol]=size(yt);
imshow(yt);
title( ['Approximation A',num2str(N)]);
xlabel(['Size : ',num2str(yrow),'*',num2str(ycol)]);
% 第2-(N+1)行显示各级高频系数
titllist={['Vertical Detail V'];['Horizontal Detail H'];['Diagonal Detail D']};
pn=2; % pn 是子图的显示序号
for pr=1:N
for pc=1:3
subplot(N+1,3,pn+2);
yt=[]; % 为了使高频细节内容(轮廓、边缘)更清晰,将高频系数增加100灰度值
yt=uint8(coef{pn})+100;
[yrow,ycol]=size(yt);
imshow(yt);
title([ titllist{pc},num2str(N-pr+1)]);
xlabel(['Size : ',num2str(yrow),'*',num2str(ycol)]);
% 每行的{dy}个图像的Y轴,显示该行高频系数对应的分解级别
if mod(pn+2,3)==1
ylabel(['Level ',num2str(N-pr+1)]);
end
pn=pn+1;
end
end
function [cA,cV,cH,cD]=mydwt2(x,Lo_D,Hi_D)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 函数 MYDWT2() 对输入的r*c维矩阵 x 进行二维小波分解,输出四个分解系数子矩阵[LL,HL,LH,HH]
% 输入参数:x —— 输入矩阵,为r*c维矩阵。
% Lo_D,Hi_D —— 小波分解的滤波器组系数向量
% 输出参数:cA,cV,cH,cD —— 是小波分解系数矩阵的四个相等大小的子矩阵
% cA:低频部分分解系数; cV:垂直方向分解系数;
% cH:水平方向分解系数; cD:对角线方向分解系数。
% Copyright by Zou Yuhua ( chenyusiyuan ), original : 2007-11-10, modified: 2008-06-04
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[row,col]=size(x);
% 读取输入矩阵的大小
for j=1:row % 首先对输入矩阵的每一行序列进行一维离散小波分解
tmp1=x(j,:);
[ca1,cd1]=mydwt(tmp1,Lo_D,Hi_D,1);
% tmp1 长度为 row ,滤波器长度为 lnf ,则 [ca1,cd1] 的总长为 ( row + lnf -1 )
x1(j,:)=[ca1,cd1]; % 将分解系数序列存入缓存矩阵 x1 中
end
[row1,col1]=size(x1); % row1=row + lnf -1, col1=col+lnf-1
for k=1:col1 % 再对缓存矩阵 x1 的每一列序列进行一维离散小波分解
tmp2=x1(:,k);
[ca2,cd2]=mydwt(tmp2,Lo_D,Hi_D,1);
x2(:,k)=[ca2,cd2]' ; % 将分解所得系数存入缓存矩阵 x2 中
% 注意不要遗漏了上一行代码中的转置符号“ ’”。 Matlab 6.5 及以下较低的版本不支
% 持行、列向量的相互赋值,故要把行向量[ca2,cd2]转置为列向量,再存入 x2 的相应列
end
[row2,col2]=size(x2);
cA=x2(1:row2/2,1:col2/2); % cA是矩阵x2的左上角部分
cV=x2(1:row2/2,col2/2+1:col2); % cV是矩阵x2的右上角部分
cH=x2(row2/2+1:row2,1:col2/2); % cH是矩阵x2的左下角部分
cD=x2(row2/2+1:row2,col2/2+1:col2); % cD是矩阵x2的右下角部分
function [cA,cD] = mydwt(x,lpd,hpd,dim)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 函数 [cA,cD]=MYDWT(X,LPD,HPD,DIM) 对输入序列x进行一维离散小波分解,输出分解序列[cA,cD]
% 输入参数:x——输入序列;
% lpd——低通滤波器;
% hpd——高通滤波器;
% dim——小波分解层数。
% 输出参数:cA——平均部分的小波分解系数;
% cD——细节部分的小波分解系数。
% Copyright by Zou Yuhua ( chenyusiyuan ), original : 2007-11-10
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
cA=x; % 初始化cA,cD
cD=[];
for i=1:dim
cvl=conv(cA,lpd);
% 低通滤波,为了提高运行速度,调用MATLAB提供的卷积函数conv()
dnl=downspl(cvl); % 通过下抽样求出平均部分的分解系数
cvh=conv(cA,hpd); % 高通滤波
dnh=downspl(cvh); % 通过下抽样求出本层分解后的细节部分系数
cA=dnl; % 下抽样后的平均部分系数进入下一层分解
cD=[cD,dnh]; % 将本层分解所得的细节部分系数存入序列cD
end
function y=downspl(x)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 函数 Y=DOWMSPL(X) 对输入序列进行下抽样,输出序列 Y。
% 下抽样是对输入序列取其偶数位,舍弃奇数位。例如 x=[x1,x2,x3,x4,x5],则 y=[x2,x4].
% Copyright by Zou Yuhua ( chenyusiyuan ), original : 2007-11-10
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
N=length(x); % 读取输入序列长度
M=floor(N/2); % 输出序列的长度是输入序列长度的一半(带小数时取整数部分)
i=1:M;
y(i)=x(2*i);