误导与重构：我们需要新的科学图像
日 14:36 来源：原载《文汇报》期
作者：江晓原
内容摘要：
作者简介：
　　科学与人文：“两种文化”的冲突及地位消长　　近年来，国内的科学史和科学哲学界人士也没有少谈“两种文化”，但我的感觉是，很长一段时间以来，这两种文化不仅没有在事实上相亲相爱，反而在观念上渐行渐远。而且有很多人已经明显感觉到，一种文化正在日益侵凌于另一种文化之上。　　C·P·斯诺著名的演讲《两种文化》已经出版了三个中译本。在斯诺演讲的年代（第一次是1959年），他觉得科学还处于被人文轻视的状况中，科学技术被认为只类似于工匠们摆弄的玩意儿。这倒很有点像中国古代的情形——工匠阶层是不能与士大夫们平起平坐的。斯诺要为科学争地位，争名份，要求让科学和人文平起平坐。他的这种主张自然得到科学界的热烈欢迎。　　从那时到现在已经过去将近半个世纪了，斯诺也去世27年（1980年）了。历史的钟摆早已摆向另一个端点。斯诺若生在今日中国，我想他恐怕要做另一个讲演了——他会为人文争地位，争名份，要求让人文能够和科学平起平坐。　　斯诺1959年的演讲是要抨击“对科学的傲慢与偏见”，然而有些人士却认为科学自身也充满着傲慢与偏见。在斯诺演讲之前7年，F·A·哈耶克早已经对此忧心忡忡了。他那本《科学的反革命——理性滥用之研究》（初版于1952年），从书名上就可以清楚感觉到他的立场和情绪。　　哈耶克的矛头并不是指向科学或科学家，而是指向那些认为科学可以解决一切问题的人。照他的意思，一个“唯科学主义”(scientism)者，很可能不是一个科学家——事实上，他认为这些唯科学主义者“几乎都不是显著丰富了我们的科学知识的人”。　　在哈耶克看来，计划经济的思想基础，就是唯科学主义——相信科学技术可以解决世间一切问题。计划经济幻想可以将人类的全部智慧集中起来，形成一个超级智慧，这个超级智慧知道人类的过去未来，知道历史发展的规律，可以为全人类指出前进的康庄大道。哈耶克反复指出：这样的超级智慧是不可能的；最终必然要求千百万人听命于一个人的头脑。而这样做的结果如何，世人早已经领教够了。　　这些关于科学本身的争论和思考，在以前闭关自守的年代，都是我们所听不到也想不到的。相反，在我们以往的宣传中，一直存在着对于科学的三大误导。
转载请注明来源：中国社会科学网
（责编：李秀伟）
用户昵称：
&(您填写的昵称将出现在评论列表中)
所有评论仅代表网友意见
最新发表的评论0条，总共0条
查看全部评论
中国社会科学网版权所有，未经书面授权禁止使用
Copyright (C)
by . all rights reserved查看: 8413|回复: 13|关注: 0
图像小波分解后重构问题
<h1 style="color:# 麦片财富积分
新手, 积分 5, 距离下一级还需 45 积分
关注者： 2
本帖最后由 yjccao123 于
15:46 编辑
我对一幅图进行三层小波分解，然后对其中的高频系数与低频系数进行了处理，我想利用新得到的这些系数矩阵重构出和原图像大小相同的图像，该怎么做？[C,S]=wavedec2(X,3,'db1');& &%对图像三级分解
ca1 = appcoef2(C,S,'db1',1);&&%提取尺度1的低频系数
ch1 = detcoef2('h',C,S,1);& & %提取尺度1的水平方向高频系数
cv1 = detcoef2('v',C,S,1);& & %提取尺度1的垂直方向高频系数
cd1= detcoef2('d',C,S,1);& &&&%提取尺度1的斜线方向高频系数
ca2 = appcoef2(C,S,'db1',2);
ch2 = detcoef2('h',C,S,2);
cv2 = detcoef2('v',C,S,2);
cd2= detcoef2('d',C,S,2);
ca3 = appcoef2(C,S,'db1',3);
ch3 = detcoef2('h',C,S,3);
cv3 = detcoef2('v',C,S,3);
cd3= detcoef2('d',C,S,3);
%%对上面12个系数矩阵进行处理
....复制代码本人刚学小波变换，希望各位大牛不吝赐教
<h1 style="color:# 麦片财富积分
sf, 表示不会
<h1 style="color:#0 麦片财富积分
关注者： 4
你已经得到了各个分量了,直接用这个函数
wrcoef2: 利用某个高频分量重建图像(与原图像图等大小)&&
<h1 style="color:# 麦片财富积分
关注者： 2
问题解决了，用的方法比较麻烦，希望能帮到一些同样困惑的朋友~
因为我小波分解都是用的matlab自带的工具箱
所以就基于小波分解函数 [C,S] = wavedec2(X,N,'wname') 来进行分析
这个函数的详细分析可以看这篇文章
（好吧，还不让发链接，哭了~~~~~~~）
百度搜索：wavedec2函数详解
了解了C和S具体的含义和内部结构之后，我们就可以进行还原了
我们都知道matlab自带的重构函数应该是这个 X = waverec2(C,S,'wname')
因此，如果想用这个函数对改动后的低频高频系数进行重构，必须改变C或者S的值
既然你了解了C与S表示什么意思，这个步骤就非常简单了
这个是我写的代码，供大家学习使用，里面有些不完善的地方，欢迎指正：%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 对小波分解后的低频系数，高频系数处理后的重构
%% 本例针对彩色图像
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all
clc
M1=double(imread('lena.jpg'))/255;
zt=3;& && && && &%小波分解级数
wtype = 'haar';&&%小波分解类型
[C,S] = wavedec2(M1, zt, wtype);%多尺度二维小波分解
%% 尺度1低频与高频系数提取
ca1 = appcoef2(C,S,wtype,1);&&%提取尺度1的低频系数
ch1 = detcoef2('h',C,S,1);& & %提取尺度1的水平方向高频系数
cv1 = detcoef2('v',C,S,1);& & %提取尺度1的垂直方向高频系数
cd1= detcoef2('d',C,S,1);& &&&%提取尺度1的斜线方向高频系数
%% 尺度2低频与高频系数提取
ca2 = appcoef2(C,S,wtype,2);&&%提取尺度2的低频系数
ch2 = detcoef2('h',C,S,2);& & %提取尺度2的水平方向高频系数
cv2 = detcoef2('v',C,S,2);& &%提取尺度2的垂直方向高频系数
cd2= detcoef2('d',C,S,2);& & %提取尺度2的斜线方向高频系数
%% 尺度3低频与高频系数提取
ca3 = appcoef2(C,S,wtype,3); %提取尺度3的低频系数
ch3 = detcoef2('h',C,S,3);& &%提取尺度3的水平方向高频系数
cv3 = detcoef2('v',C,S,3);& &%提取尺度3的垂直方向高频系数
cd3= detcoef2('d',C,S,3);& & %提取尺度3的斜线方向高频系数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 下面对这些系数矩阵进行处理，为了调试方便起见，我们对高频系数的进行去噪
%% 对低频系数不进行处理（另外由于C中只有尺度3的低频系数，因此处理低频系数要处理尺度3的）
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 对尺度1的高频系数去噪
d_ch1(:,:,1)=wiener2(ch1(:,:,1),[5 5]);
d_ch1(:,:,2)=wiener2(ch1(:,:,2),[5 5]);
d_ch1(:,:,3)=wiener2(ch1(:,:,3),[5 5]); %水平方向去噪
d_cv1(:,:,1)=wiener2(cv1(:,:,1),[5 5]);
d_cv1(:,:,2)=wiener2(cv1(:,:,2),[5 5]);
d_cv1(:,:,3)=wiener2(cv1(:,:,3),[5 5]); %垂直方向去噪
d_cd1(:,:,1)=wiener2(cd1(:,:,1),[5 5]);
d_cd1(:,:,2)=wiener2(cd1(:,:,2),[5 5]);
d_cd1(:,:,3)=wiener2(cd1(:,:,3),[5 5]); %斜线方向去噪
%% 对尺度2的高频系数去噪
d_ch2(:,:,1)=wiener2(ch2(:,:,1),[3 3]);
d_ch2(:,:,2)=wiener2(ch2(:,:,2),[3 3]);
d_ch2(:,:,3)=wiener2(ch2(:,:,3),[3 3]); %水平方向去噪
d_cv2(:,:,1)=wiener2(cv2(:,:,1),[3 3]);
d_cv2(:,:,2)=wiener2(cv2(:,:,2),[3 3]);
d_cv2(:,:,3)=wiener2(cv2(:,:,3),[3 3]); %垂直方向去噪
d_cd2(:,:,1)=wiener2(cd2(:,:,1),[3 3]);
d_cd2(:,:,2)=wiener2(cd2(:,:,2),[3 3]);
d_cd2(:,:,3)=wiener2(cd2(:,:,3),[3 3]); %斜线方向去噪
%% 对尺度3的高频系数去噪
d_ch3(:,:,1)=wiener2(ch3(:,:,1),[2 2]);
d_ch3(:,:,2)=wiener2(ch3(:,:,2),[2 2]);
d_ch3(:,:,3)=wiener2(ch3(:,:,3),[2 2]); %水平方向去噪
d_cv3(:,:,1)=wiener2(cv3(:,:,1),[2 2]);
d_cv3(:,:,2)=wiener2(cv3(:,:,2),[2 2]);
d_cv3(:,:,3)=wiener2(cv3(:,:,3),[2 2]); %垂直方向去噪
d_cd3(:,:,1)=wiener2(cd3(:,:,1),[2 2]);
d_cd3(:,:,2)=wiener2(cd3(:,:,2),[2 2]);
d_cd3(:,:,3)=wiener2(cd3(:,:,3),[2 2]); %斜线方向去噪
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 下一步进行重构了，如何把处理完的高频系数整合到C中呢，
%% 我们知道C的结构是[A(N)|H(N)|V(N)|D(N)|H(N-1)|V(N-1)|D(N-1)|H(N-2)|V(N-2)|D(N-2)|...|H(1)|V(1)|D(1)]
%% 在本程序中 C的结构是 [A(3)|H(3)|V(3)|D(3)|H(2)|V(2)|D(2)|H(1)|V(1)|D(1)]
%% S(1,:)是A(3)大小，S(2,:)是H(3)大小，S(3，：)是H(2)大小，S(4,:)是H(1)大小，S(5,:)是原图像大小
%% S(1,:)与S(2:)完全一致，即A(N)与H(N)大小是一样的
%% 可以容易算出C行向量的大小应该为
%% S(1,1)*S(1,2)*S(1,3)*4+S(3,1)*S(3,2)*S(3,3)*3+...S(N-1,1)*S(N-1,2)*S(N-1,3)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
len1=S(1,1)*S(1,2)*S(1,3);& & %A(3)，H(3),V(3),D(3)的大小
len2=S(3,1)*S(3,2)*S(3,3);& & %H(2),V(2),D(2)的大小
len3=S(4,1)*S(4,2)*S(4,3);& & %H(1),V(1),D(1)的大小
C(1:len1)=ca3(1:end);& &%%将A(3)放入相应位置
C(len1+1:2*len1)=d_ch3(1,end); %%将H(3),V(3),D(3)放入相应位置
C(2*len1+1:3*len1)=d_cv3(1,end);
C(3*len1+1:4*len1)=d_cd3(1,end);
C(4*len1+1:4*len1+len2)=d_ch2(1:end); %%将H(2),V(2),D(2)放入相应位置
C(4*len1+len2+1:4*len1+2*len2)=d_cv2(1:end);
C(4*len1+2*len2+1:4*len1+3*len2)=d_cd2(1:end);
C(4*len1+3*len2+1:4*len1+3*len2+len3)=d_ch1(1:end); %%将H(1),V(1),D(1)放入相应位置
C(4*len1+3*len2+len3+1:4*len1+3*len2+2*len3)=d_cv1(1:end);
C(4*len1+3*len2+2*len3+1:4*len1+3*len2+3*len3)=d_cd1(1:end);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 把一系列大小不同的矩阵依次放入一个大矩阵，
%% 肯定有更简单的方法，不用判断起始位置终止位置，
%% 本人能力有限，只能用这种笨方法
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 用函数waverec2进行重构
re_M1_1=waverec2(C,S,wtype);
figure(1);
subplot(1,2,1);imshow(M1);
subplot(1,2,2);imshow(re_M1_1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 因为考虑到只要知道C,S就可以进行重构，因此突发奇想自己构造C,S，可不可以呢
%% 我们以尺度2的低频高频系数为例，进行重构
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
S1=[S(3,1),S(3,2),S(3,3);S(3,1),S(3,2),S(3,3);S(4,1),S(4,2),S(4,3);S(5,1),S(5,2),S(5,3)];
leng1=S1(1,1)*S1(1,2)*S1(1,3);& & %A(2)，H(2),V(2),D(2)的大小
leng2=S1(3,1)*S1(3,2)*S1(3,3);& & %H(1),V(1),D(1)的大小
C1=zeros(size(leng1*4+leng2*3));&&
C1(1:leng1)=ca2(1:end);&&%%将A(2)放入相应位置
C1(leng1+1:2*leng1)=d_ch2(1:end); %%将H(2),V(2),D(2)放入相应位置
C1(2*leng1+1:3*leng1)=d_cv2(1:end);
C1(3*leng1+1:4*leng1)=d_cd2(1:end);
C1(4*leng1+1:4*leng1+leng2)=d_ch1(1:end); %%将H(1),V(1),D(1)放入相应位置
C1(4*leng1+leng2+1:4*leng1+2*leng2)=d_cv1(1:end);
C1(4*leng1+2*leng2+1:4*leng1+3*leng2)=d_cd1(1:end);
re_M1_2=waverec2(C1,S1,wtype);
figure(2);
subplot(1,2,1);imshow(M1);
subplot(1,2,2);imshow(re_M1_2);
复制代码
<h1 style="color:# 麦片财富积分
关注者： 2
libralibra 发表于
你已经得到了各个分量了,直接用这个函数
wrcoef2: 利用某个高频分量重建图像(与原图像图等大小)&&
呵呵，是啊 我在csdn也发了个，问题解决了 多谢
<h1 style="color:# 麦片财富积分
关注者： 1
你好，我最近在做这个，我是对分解得到的各层高频系数（水平，垂直，斜对角方向）进行软阈值去噪，完了之后也想把最低层的低频系数与处理过的各高层系数重构图像，但是不知道要怎么重构，可以讲解一下吗，十分谢谢
i=imread('C:\Documents and Settings\lena.bmp');
x=imnoise(i,'gaussian',0,0.01); %加入均值为0，方差为0.01的高斯噪声
%subplot(231);imshow(i);title('原始图像');
subplot(232); imshow(x);title('含噪图像');
s=double(x);
[a,b]=size(s);
%n=[1,2,3];%设置尺度向量
p=0.01*(sqrt(2*log(a*b)));%阈值
[C,S]=wavedec2(s,3,'db1');%采用db1小波对图像x1进行3层分解,c为各层分解系数s为各层分解系数长度
a3=appcoef2(C,S,'db1',3);
[h3 v3 d3]=detcoef2('all',C,S,3);
[h2 v2 d2]=detcoef2('all',C,S,2);
[h1 v1 d1]=detcoef2('all',C,S,1);
%对系数进行软阈值处理:把信号的绝对值与阈值进行比较
%小于或等于阈值的点变为零，大于阈值的点为该点值与阈值的差值。
H3=wthresh(h3,'s',p);
H2=wthresh(h2,'s',p);
H1=wthresh(h1,'s',p);
V3=wthresh(v3,'s',p);
V2=wthresh(v2,'s',p);
V1=wthresh(v1,'s',p);
D3=wthresh(d3,'s',p);
D2=wthresh(d2,'s',p);
D1=wthresh(d1,'s',p);
<h1 style="color:# 麦片财富积分
关注者： 1
你好，我最近在做这个，我是对分解得到的各层高频系数（水平，垂直，斜对角方向）进行软阈值去噪，完了之后也想把最低层的低频系数与处理过的各高层系数重构图像，但是不知道要怎么重构，可以讲解一下吗，十分谢谢
i=imread('C:\Documents and Settings\lena.bmp');
x=imnoise(i,'gaussian',0,0.01); %加入均值为0，方差为0.01的高斯噪声
%subplot(231);imshow(i);title('原始图像');
subplot(232); imshow(x);title('含噪图像');
s=double(x);
[a,b]=size(s);
%n=[1,2,3];%设置尺度向量
p=0.01*(sqrt(2*log(a*b)));%阈值
[C,S]=wavedec2(s,3,'db1');%采用db1小波对图像x1进行3层分解,c为各层分解系数s为各层分解系数长度
a3=appcoef2(C,S,'db1',3);
[h3 v3 d3]=detcoef2('all',C,S,3);
[h2 v2 d2]=detcoef2('all',C,S,2);
[h1 v1 d1]=detcoef2('all',C,S,1);
%对系数进行软阈值处理:把信号的绝对值与阈值进行比较
%小于或等于阈值的点变为零，大于阈值的点为该点值与阈值的差值。
H3=wthresh(h3,'s',p);
H2=wthresh(h2,'s',p);
H1=wthresh(h1,'s',p);
V3=wthresh(v3,'s',p);
V2=wthresh(v2,'s',p);
V1=wthresh(v1,'s',p);
D3=wthresh(d3,'s',p);
D2=wthresh(d2,'s',p);
D1=wthresh(d1,'s',p);
<h1 style="color:# 麦片财富积分
关注者： 1
你好，我最近在做这个，我是对分解得到的各层高频系数（水平，垂直，斜对角方向）进行软阈值去噪，完了之后也想把最低层的低频系数与处理过的各高层系数重构图像，但是不知道要怎么重构，可以讲解一下吗，十分谢谢
i=imread('C:\ena.bmp');
x=imnoise(i,'gaussian',0,0.01); %加入均值为0，方差为0.01的高斯噪声
%subplot(231);imshow(i);title('原始图像');
subplot(232); imshow(x);title('含噪图像');
s=double(x);
[a,b]=size(s);
%n=[1,2,3];%设置尺度向量
p=0.01*(sqrt(2*log(a*b)));%阈值
[C,S]=wavedec2(s,3,'db1');%采用db1小波对图像x1进行3层分解,c为各层分解系数s为各层分解系数长度
a3=appcoef2(C,S,'db1',3);
[h3 v3 d3]=detcoef2('all',C,S,3);
[h2 v2 d2]=detcoef2('all',C,S,2);
[h1 v1 d1]=detcoef2('all',C,S,1);
%对系数进行软阈值处理:把信号的绝对值与阈值进行比较
%小于或等于阈值的点变为零，大于阈值的点为该点值与阈值的差值。
H3=wthresh(h3,'s',p);
H2=wthresh(h2,'s',p);
H1=wthresh(h1,'s',p);
V3=wthresh(v3,'s',p);
V2=wthresh(v2,'s',p);
V1=wthresh(v1,'s',p);
D3=wthresh(d3,'s',p);
D2=wthresh(d2,'s',p);
D1=wthresh(d1,'s',p);
<h1 style="color:# 麦片财富积分
关注者： 2
lizhen213 发表于
你好，我最近在做这个，我是对分解得到的各层高频系数（水平，垂直，斜对角方向）进行软阈值去噪，完了之后 ...len1=S(1,1)*S(1,2)*S(1,3);& & %A(3)，H(3),V(3),D(3)的大小
len2=S(3,1)*S(3,2)*S(3,3);& & %H(2),V(2),D(2)的大小
len3=S(4,1)*S(4,2)*S(4,3);& & %H(1),V(1),D(1)的大小
C(1:len1)=a3(1:end);& &%%将A(3)放入相应位置
C(len1+1:2*len1)=H3(1,end); %%将H(3),V(3),D(3)放入相应位置
C(2*len1+1:3*len1)=V3(1,end);
C(3*len1+1:4*len1)=D3(1,end);
C(4*len1+1:4*len1+len2)=H2(1:end); %%将H(2),V(2),D(2)放入相应位置
C(4*len1+len2+1:4*len1+2*len2)=V2(1:end);
C(4*len1+2*len2+1:4*len1+3*len2)=D2(1:end);
C(4*len1+3*len2+1:4*len1+3*len2+len3)=H1(1:end); %%将H(1),V(1),D(1)放入相应位置
C(4*len1+3*len2+len3+1:4*len1+3*len2+2*len3)=V1(1:end);
C(4*len1+3*len2+2*len3+1:4*len1+3*len2+3*len3)=D1(1:end);
%% 用函数waverec2进行重构
re_M1_1=waverec2(C,S,'db1');复制代码试试行不行~
<h1 style="color:# 麦片财富积分
关注者： 1
yjccao123 发表于
试试行不行~
非常感谢啊
站长推荐 /3
Powered by(思远老师的好文章留着)
去年11月发布了一系列有关小波变换和图像处理的文章，把学习小波过程中的心得体会和编写的程序放在网上和大家共享交流。半年来，感谢大家的关注和帮助，在相互的讨论交流中，我不断地从大家提出的问题中拓展自己的知识面，对小波的理论及其应用有了更深入的了解和掌握。的确，分享，才是知识永恒之道！
接下来有新的课题，并且要准备考博，研学征途马上要展开新的探索。在此之前，得先整理下博客文章，根据和大家讨论交流中发现的问题，对博客中的程序进行修正。今天发布的文章，就是对小波图像分解和重构程序作出重要的修正，这个程序是了解小波、使用小波最根本的，希望对大家学习小波有所帮助。
首先指出我原来发布的小波图像分解重构的文章中存在的问题。相关的文章有：
1、自己动手编写小波信号分解与重构的Matlab程序
2、用自编的程序实现小波图像分解与重构
其存在的问题有：
（1）程序所用的小波函数只有非标准的Haar小波，其滤波器组为 Lo_D=[1/2 1/2], Hi_D=[-1/2 1/2]，是固化在 mydwt2.m 的程序中的，不能选择其他的小波函数；
（2）非标准的Haar小波，其分解出来的系数矩阵中，高频系数的细节内容（轮廓、边缘等特征）不明显；
（3）函数 mydwt2 中列变换的矩阵对象为输入矩阵，这是错误的，其矩阵对象应该是行变换后的缓存矩阵；
（4）函数 mydwt2 的输出用[LL,HL,LH,HH]表示，不是很规范，应改为[cA,cV,cH,cD]来表示，即一级小波变换输出的系数矩阵有4个部分：平均部分、垂直细节部分、水平细节部分和对角线细节部分。
（5）函数 mydwt2 的输出[LL,HL,LH,HH]中，HL 取的是小波行变换和列变换后所得矩阵的左下角部分，即垂直细节部分，但在最终的显示结果中 HL 是在右上角的位置的；类&#20284;地，LH 取行、列变换所得矩阵的右上角部分，即水平细节部分，但最终显示位置是在左下角，当时这样定义，是我从教材中看到小波一级分解后所得矩阵的平均、细节部分分布从上到下、从左到右依次为：平均部分、水平细节部分、垂直细节部分、对角线细节部分。但通过编程实现和仔细分析，我发现其分布次序应该是平均部分、垂直细节部分、水平细节部分、对角线细节部分才对，即水平和细节的位置对调，具体过程下文会给出图示，希望能和大家讨论下哪个次序才是正确的。
（5）函数 mywavedec2 的输出 y 是与输入矩阵 x 相同大小的矩阵，并且已将N级分解后所有的平均、细节系数组合成一体的。实际上，这种定义只对Haar小波有效。因为对于Haar小波来说，其滤波器长度仅为2，输入序列与滤波器系数卷积后，输出序列与输入序列长度是相同的（在Matlab中，若信号长度为n，滤波器长度为2N，则小波分解后cA、cD的长度均为：lnf
= floor((n-1)/2)&#43;N，若2N=2，则 lnf*2 = n），从而使小波系数矩阵与图像矩阵的大小相同。若使用其它类型的小波，由于 2N & 2，lnf*2 & n，这样函数 mydwt2 输出的系数矩阵将比输入矩阵大，那么 mywavedec2 的输出 y 将会显示出平均与细节部分边界重叠的问题。这说明了试图把平均和细节系数组合成一体，以显示小波分解的塔式结构，理论上是不可行的，至少对于Matlab来说是不可行的。当然，由公式
lnf = floor((n-1)/2)&#43;N ，我们可以对卷积后的序列进行截取，使卷积输出序列与输入序列长度相同，但通过编程实现发现，截取会找出下一级分解错误，所以函数 mywavedec2 的输出 y 的&#26684;式要作适当调整。
（6）原程序中要调用 modmat 函数对图像矩阵进行修剪，使之能被 2 的 N 次方整除，主要是为了生成塔式结构图像而设的，对上述问题修正后，这个&modmat 函数已不需使用了。
针对上述问题，我对程序作了以下几点修正：
（1）函数 mywavedec2 / mywaverec2 的输入参数增加了小波函数'wname'，可以选择任意类型的小波函数来对图像进行小波分解和重构。相应地，调用Matlab函数 wfilters 来求取小波函数的分解滤波器组系数[Lo_D,Hi_D]和重构滤波器组系数[Lo_R,Hi_R]。
（2）函数&mydwt2 的输出应改为[cA,cV,cH,cD]来表示，注意这4个部分的位置分布：cA 在系数矩阵的左上角，cV 在右上角，cH 在左下角，cD 在右下角。
（3）函数&mydwt2 的列变换，其处理对象由输入矩阵 x 改为 行变换后的缓存矩阵 x1，列变换后得到缓存矩阵 x2。输出 cA 取 x2 的左上角1/4部分，cV 取 x2 的右上角1/4部分，cH 取 x2 的左下角1/4部分，cD 取 x2 的右下角1/4部分。
（4）函数&mywavedec2 的输出 y &#26684;式改为与Matlab的 wavedec2 类&#20284;，但我的是以“细胞矩阵”的&#26684;式，而 wavedec2 是以“向量序列”的&#26684;式输出小波系数。小波分解的结果显示方式不再使用塔式结构方式，而是用 subplot 函数把各级平均、细节系数显示出来。
（5）重构函数 mywaverec2 根据分解函数 mywavedec2 的输出&#26684;式做了相应的修改，程序更简洁、易懂和合理。并且&mywaverec2 的重构级数可以与分解级数不同，能够实现不同级别的小波重构。
（6）函数&myidwt2 的行、列变换中，卷积后要对卷积输出序列进行裁剪，消除分解滤波器和重构滤波器对序列长度的影响，使重构矩阵的大小与上一分解级的系数矩阵大小一致。
以上是主要的修正点，更详细的修改信息请看下一篇文章给出的程序及其注释。
阅读(...) 评论()【转】小波图像分解与重构
11:15&4728人阅读&(6)&&
图像融合是将两幅或多幅图像融合在一起，以获取对同一场景的更为精确、更为全面、更为可靠的图像描述。融合算法应该充分利用各原图像的互补信息，使融合后的图像更适合人的视觉感受，适合进一步分析的需要；并且应该统一编码，压缩数据量，以便于传输。
图像融合可分为三个层次:
像素级融合
特征级融合
决策级融合
其中像素级融合是最低层次的融合，也是后两级的基础。它是将各原图像中对应的像素进行融合处理，保留了尽可能多的图像信息， 精度比较高，
因而倍受人们的重视。像素级的图像融合方法大致可分为三大类:
简单的图像融合方法
基于塔形分解（如Laplace塔形分解、比率塔等）的图像融合方法
基于小波变换的图像融合方法
小波变换是图像的多尺度、多分辨率分解，它可以聚焦到图像的任意细节，被称为数学上的显微镜。近年来，随着小波理论及其应用的发展，已将小波多分辨率分解用于像素级图像融合。小波变换的固有特性使其在图像处理中有如下优点：
完善的重构能力，保证信号在分解过程中没有信息损失和冗余信息；
把图像分解成平均图像和细节图像的组合，分别代表了图像的不同结构，因此容易提取原始图像的结构信息和细节信息；
具有快速算法，它在小波变换中的作用相当于FFT算法在傅立叶变换中的作用，为小波变换应用提供了必要的手段；
二维小波分析提供了与人类视觉系统方向相吻合的选择性图像。
——像素级图像融合的主要步骤
以两幅图像的融合为例。设A，B为两幅原始图像，F为融合后的图像。若对二维图像进行N层的小波分解，最终将有(3N+1)个不同频带，其中包含3N
个高频子图像和1个低频子图像。其融合处理的基本步骤如下：&
（1）对每一原图像分别进行小波变换，建立图像的小波塔型分解；&
（2）对各分解层分别进行融合处理。各分解层上的不同频率分量可采用不同的融合算子进行融合处理，最终得到融合后的小波金字塔；&
（3）对融合后所得小波金字塔进行小波重构，所得到的重构图像即为融合图像。
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
在图像融合过程中，小波基的种类和小波分解的层数对融合效果有很大的影响，对特定的图像来说，哪一种小波基的融合效果最好，分解到哪一层最合适，都是需要考虑的问题。为此可以通过引入融合效果的评价来构成一个闭环系统。如图2所示。
10:56&9744人阅读&(31)&&
P.S.：（）感谢网友‘李明杨艳’指出了本文程序中一维信号小波分解重构程序mydwt和myidwt存在的一个大Bug，现已修正，请参见今天发表的文章《一维信号的小波分解与重构程序》。
P.S.：（）去年11月发布了一系列有关小波变换和图像处理的文章，把学习小波过程中的心得体会和编写的程序放在网上和大家共享交流。半年来，感谢大家的关注和帮助，在相互的讨论交流中，我不断地从大家提出的问题中拓展自己的知识面，对小波的理论及其应用有了更深入的了解和掌握。根据和大家讨论交流中发现的问题，对博客中的程序进行修正。有关小波图像分解和重构的两篇文章中分享的程序，存在下列问题：
（1）程序所用的小波函数只有非标准的Haar小波，其滤波器组为 Lo_D=[1/2 1/2], Hi_D=[-1/2
1/2]，是固化在 mydwt2.m 的程序中的，不能选择其他的小波函数；
（2）非标准的Haar小波，其分解出来的系数矩阵中，高频系数的细节内容（轮廓、边缘等特征）不明显；
（3）函数 mydwt2
中列变换的矩阵对象为输入矩阵，这是错误的，其矩阵对象应该是行变换后的缓存矩阵；
（4）函数 mydwt2
的输出用[LL,HL,LH,HH]表示，不是很规范，应改为[cA,cV,cH,cD]来表示，即一级小波变换输出的系数矩阵有4个部分：平均部分、垂直细节部分、水平细节部分和对角线细节部分。
（5）函数 mywavedec2 的输出
y 是与输入矩阵 x
相同大小的矩阵，并且已将N级分解后所有的平均、细节系数组合成一体的。实际上，这种定义只对Haar小波有效。
（6）原程序中要调用 modmat
函数对图像矩阵进行修剪，使之能被 2 的 N 次方整除，主要是为了生成塔式结构图像而设的，对上述问题修正后，这个 modmat
函数已不需使用了。
针对上述问题，我对程序作了修正，发布在今天的3篇文章里，请大家点击查看。新修正的程序更为简洁易懂，功能也有所增强，可以用任意的小波函数进行小波分解，可根据小波分解系数矩阵重构出指定分解级的低频系数和原始图像。
1、《小波图像分解与重构程序存在的问题与解决办法》
2、《小波图像分解 Matlab 程序
- V2.0版》
3、《小波图像重构 Matlab 程序
- V2.0版》
（）Matlab小波分析工具箱丰富的函数和强大的仿真功能为我们学习小波、用好小波提供了方便、快捷的途径，但是，如果我们要深入掌握小波分析的原理，真正学好、用好小波，就应该尽量用自己编写的程序去实现小波变换和信号分析，尽量在自己的程序中少调用Matlab提供的函数，多用自己的理解去编写相关的小波函数，这样的过程是一个探索、求知的过程，更能让我们体会到小波的强大和学习的乐趣。下面，我把自己编写的小波一维、二维信号分解和重构Matlab程序共享出来，也希望有朋友共享自编的程序，共同学习，提高程序的效率和简洁性。
首先要说明的一点是，虽然是自己编写Matlab程序，但并不是说一点也不用Matlab的自带函数。我们要编写的是实现小波变换的主要功能函数，而绘图等基本功能还是要用到Matlab函数的。而且，根据小波变换的滤波器组原理，原始信号要通过低通、高通滤波器处理，这里就涉及到卷积这一运算步骤。卷积——FFT算法的实现，相信很多朋友都能用Matlab、C语言等来实现，不过与Matlab自带的用机器语言编写的FFT程序相比，运算速度一般会慢几倍、几十倍。所以，我的程序里边涉及卷积的就直接调用Matlab的conv()函数了。
我们知道，小波变换的一级分解过程是，原始信号分别进行低通、高通滤波，再分别进行二元下抽样，就得到低频、高频（也称为平均、细节）两部分系数；而多级分解则是对上一级分解得到的低频系数再进行小波分解，是一个递归过程。以下是一维小波分解的程序：
[cA,cD] = mydwt(x,lpd,hpd,dim);
[cA,cD]=MYDWT(X,LPD,HPD,DIM)
对输入序列x进行一维离散小波分解，输出分解序列[cA,cD]
输入参数：x——输入序列；
%&&&&&&&&&
lpd——低通滤波器；
%&&&&&&&&&
hpd——高通滤波器；
%&&&&&&&&&
dim——小波分解级数。
输出参数：cA——平均部分的小波分解系数；
%&&&&&&&&&&
cD——细节部分的小波分解系数。
cA=x;&&&&&&
% 初始化cA，cD
for i=1:dim
cvl=conv(cA,lpd);&& %
低通滤波，为了提高运行速度，调用MATLAB提供的卷积函数conv()
dnl=downspl(cvl);&& %
通过下抽样求出平均部分的分解系数
cvh=conv(cA,hpd);&& %
dnh=downspl(cvh);&& %
通过下抽样求出本层分解后的细节部分系数
cA=&&&&&&&&&&&&
% 下抽样后的平均部分系数进入下一层分解
cD=[cD,dnh];&&&&&&&
% 将本层分解所得的细节部分系数存入序列cD
y=downspl(x);
% 函数 Y=DOWMSPL(X)
对输入序列进行下抽样，输出序列 Y。
% 下抽样是对输入序列取其偶数位，舍弃奇数位。例如
x=[x1,x2,x3,x4,x5]，则 y=[x2,x4].
N=length(x);&&&&&&&
% 读取输入序列长度
M=floor(N/2);&&&&&&&
% 输出序列的长度是输入序列长度的一半（带小数时取整数部分）
y(i)=x(2*i);
而重构则是分解的逆过程，对低频系数、高频系数分别进行上抽样和低通、高通滤波处理。要注意重构时同一级的低频、高频系数的个数必须相等。
function y =
myidwt(cA,cD,lpr,hpr);
% 函数 MYIDWT()
对输入的小波分解系数进行逆离散小波变换，重构出信号序列 y
% 输入参数：cA ——
平均部分的小波分解系数；
%&&&&&&&&&&
cD —— 细节部分的小波分解系数；
%&&&&&&&&&&
lpr、hpr —— 重构所用的低通、高通滤波器。
lca=length(cA);&&&&&&&&&&&&
% 求出平均、细节部分分解系数的长度
lcd=length(cD);
(lcd)&=(lca)&&&&&&&&&
% 每一层重构中，cA 和 cD 的长度要相等，故每层重构后，
&&&&&&&&&&&&&&&&&&&&&&&&&&&
% 若lcd小于lca，则重构停止，这时的 cA 即为重构信号序列 y 。
upl=upspl(cA);&&&&&&&&&
% 对平均部分系数进行上抽样
cvl=conv(upl,lpr);&&&&&
% 低通卷积
cD_up=cD(lcd-lca+1:lcd);&&&
% 取出本层重构所需的细节部分系数，长度与本层平均部分系数的长度相等
uph=upspl(cD_up);&&&&&&
% 对细节部分系数进行上抽样
cvh=conv(uph,hpr);&&&&&
% 高通卷积
cA=cvl+&&&&&&&&&&&&
% 用本层重构的序列更新cA，以进行下一层重构
cD=cD(1:lcd-lca);&&&&&&
% 舍弃本层重构用到的细节部分系数，更新cD
lca=length(cA);&&&&&&&&
% 求出下一层重构所用的平均、细节部分系数的长度
lcd=length(cD);
end&&&&&&&&&&&&&&&&&&&&&&&&
% lcd & lca，重构完成，结束循环
y=cA;&&&&&&&&&&&&&&&&&&&&&&
% 输出的重构序列 y 等于重构完成后的平均部分系数序列 cA
y=upspl(x);
% 函数 Y = UPSPL(X)
对输入的一维序列x进行上抽样，即对序列x每个元素之间
% 插零，例如 x=[x1,x2,x3,x4],上抽样后为
y=[x1,0,x2,0,x3,0,x4];
N=length(x);&&&&&&&
% 读取输入序列长度
M=2*N-1;&&&&&&&&&&&
% 输出序列的长度是输入序列长度的2倍再减一
i=1:M&&&&&&&&&&
% 输出序列的偶数位为0，奇数位按次序等于相应位置的输入序列元素
if mod(i,2)
y(i)=x((i+1)/2);
&&&&我们知道，二维小波分解重构可以用一系列的一维小波分解重构来实现。以下程序是基于Haar小波的二维小波分解和重构过程：
[LL,HL,LH,HH]=mydwt2(x);
% 函数 MYDWT2() 对输入的r*c维矩阵 x
进行二维小波分解，输出四个分解系数子矩阵[LL,HL,LH,HH]
% 输入参数：x ——
输入矩阵，为r*c维矩阵。
% 输出参数：LL,HL,LH,HH ——
是分解系数矩阵的四个相等大小的子矩阵,大小均为 r/2 * c/2 维
%&&&&&&&&&&&&&&
LL：低频部分分解系数；&&&
HL：垂直方向分解系数；
%&&&&&&&&&&&&&&
LH：水平方向分解系数；&&&
HH：对角线方向分解系数。
lpd=[1/2 1/2];hpd=[-1/2
1/2];&&&&&&&&&&
% 默认的低通、高通滤波器
[row,col]=size(x);&&&&&&&&&&&&&&&&&&&&&
% 读取输入矩阵的大小
j=1:row&&&&&&&&&&&&&&&&&&&&&&&&&&&&
% 首先对输入矩阵的每一行序列进行一维离散小波分解
tmp1=x(j,:);
[ca1,cd1]=mydwt(tmp1,lpd,hpd,1);
x(j,:)=[ca1,cd1];&&&&&&&&&&&&&&&&&&
% 将分解系数序列再存入矩阵x中，得到[L|H]
k=1:col&&&&&&&&&&&&&&&&&&&&&&&&&&&&
% 再对输入矩阵的每一列序列进行一维离散小波分解
tmp2=x(:,k);
[ca2,cd2]=mydwt(tmp2,lpd,hpd,1);
x(:,k)=[ca2,cd2];&&&&&&&&&&&&&&&&&&
% 将分解所得系数存入矩阵x中，得到[LL,Hl;LH,HH]
LL=x(1:row/2,1:col/2);&&&&&&&&&&&&&&&&&
% LL是矩阵x的左上角部分
LH=x(row/2+1:row,1:col/2);&&&&&&&&&&&&&
% LH是矩阵x的左下角部分
HL=x(1:row/2,col/2+1:col);&&&&&&&&&&&&&
% HL是矩阵x的右上角部分
HH=x(row/2+1:row,col/2+1:col);&&&&&&&&&
% HH是矩阵x的右下角部分
y=myidwt2(LL,HL,LH,HH);
% 函数 MYIDWT2()
对输入的子矩阵序列进行逆小波变换，重构出矩阵 y
% 输入参数：LL,HL,LH,HH —— 是四个大小均为
r*c 维的矩阵
% 输出参数：y —— 是一个大小为 2r*2c
lpr=[1 1];hpr=[1
-1];&&&&&&&&&&
% 默认的低通、高通滤波器
tmp_mat=[LL,HL;LH,HH];&&&&&&&&&
% 将输入的四个矩阵组合为一个矩阵
[row,col]=size(tmp_mat);&&&&&&&
% 求出组合矩阵的行列数
k=1:col&&&&&&&&&&&&&&&&&&&&
% 首先对组合矩阵tmp_mat的每一列，分开成上下两半
ca1=tmp_mat(1:row/2,k);&&&&
% 分开的两部分分别作为平均系数序列ca1、细节系数序列cd1
cd1=tmp_mat(row/2+1:row,k);
tmp1=myidwt(ca1,cd1,lpr,hpr);&& %
yt(:,k)=tmp1;&&&&&&&&&&&&&&&
% 将重构序列存入待输出矩阵 yt 的相应列，此时 y=[L|H]
j=1:row&&&&&&&&&&&&&&&&&&&&
% 将输出矩阵 y 的每一行，分开成左右两半
ca2=yt(j,1:col/2);&&&&&&&&&&
% 分开的两部分分别作为平均系数序列ca2、细节系数序列cd2
cd2=yt(j,col/2+1:col);
tmp2=myidwt(ca2,cd2,lpr,hpr);&& %
yt(j,:)=tmp2;&&&&&&&&&&&&&&&
% 将重构序列存入待输出矩阵 yt 的相应行，得到最终的输出矩阵 y=yt
&&&&我编写程序时尽量实现模块化，方便不同功能程序之间的调用。以上程序实现了一维、二维小波分解与重构的功能，接下来我们就可以此为基础，探讨小波信号压缩和图像融合等技术了。
12:08&7634人阅读&(25)&&
去年11月发布了一系列有关小波变换和图像处理的文章，把学习小波过程中的心得体会和编写的程序放在网上和大家共享交流。半年来，感谢大家的关注和帮助，在相互的讨论交流中，我不断地从大家提出的问题中拓展自己的知识面，对小波的理论及其应用有了更深入的了解和掌握。根据和大家讨论交流中发现的问题，对博客中的程序进行修正。有关小波图像分解和重构的两篇文章中分享的程序，存在下列问题：
（1）程序所用的小波函数只有非标准的Haar小波，其滤波器组为 Lo_D=[1/2 1/2],
Hi_D=[-1/2 1/2]，是固化在 mydwt2.m 的程序中的，不能选择其他的小波函数；
（2）非标准的Haar小波，其分解出来的系数矩阵中，高频系数的细节内容（轮廓、边缘等特征）不明显；
（3）函数 mydwt2
中列变换的矩阵对象为输入矩阵，这是错误的，其矩阵对象应该是行变换后的缓存矩阵；
（4）函数 mydwt2
的输出用[LL,HL,LH,HH]表示，不是很规范，应改为[cA,cV,cH,cD]来表示，即一级小波变换输出的系数矩阵有4个部分：平均部分、垂直细节部分、水平细节部分和对角线细节部分。
（5）函数 mywavedec2 的输出 y 是与输入矩阵 x
相同大小的矩阵，并且已将N级分解后所有的平均、细节系数组合成一体的。实际上，这种定义只对Haar小波有效。
（6）原程序中要调用 modmat 函数对图像矩阵进行修剪，使之能被 2 的 N
次方整除，主要是为了生成塔式结构图像而设的，对上述问题修正后，这个 modmat 函数已不需使用了。
针对上述问题，我对程序作了修正，发布在今天的3篇文章里，请大家点击查看。新修正的程序更为简洁易懂，功能也有所增强，可以用任意的小波函数进行小波分解，可根据小波分解系数矩阵重构出指定分解级的低频系数和原始图像。
1、《小波图像分解与重构程序存在的问题与解决办法》
2、《小波图像分解 Matlab 程序 -
3、《小波图像重构 Matlab 程序 -
&&&&&&上一篇文章中我们实现了小波的一维、二维信号分解与重构，其中的二维信号分解与重构，只要稍作修改，就可以实现图像的分解和重构了。修改的工作，主要是对图像信号进行规范化处理、数据格式转换和绘图细节处理等。
&&&&&&&简单起见，我们从黑白（灰度）图像的分解、重构说起，因为彩色图像的处理要复杂一点。在本文中，我们使用著名的Lena图作为原始图像。
首先，为了实现图像的N层分解，对一幅m行n列的黑白图像，我们要对其进行规范化处理，使其能被2的N次方整除。以下的modmat()
函数实现此功能：
y=modmat(x,dim)
% 函数 MODMAT()
对输入矩阵x进行规范化，使其行列数均能被 2^dim 整除
% 输入参数：x —— r*c
%&&&&&&&&&&
dim —— 矩阵重构的维数
% 输出参数：y —— rt*ct
维矩阵，mod(rt,2^dim)=0，mod(ct,2^dim)=0
[row,col]=size(x);&&&&&&&&&
% 求出输入矩阵的行列数row,col
mod(row,2^dim);&&&
% 将row,col分别减去本身模 2^dim 得到的数
mod(col,2^dim);&&&
% 所得的差为rt、ct，均能被 2^dim 整除
y=x(1:rt,1:ct);&&&&&&&&&&&&
% 输出矩阵 y 为输入矩阵 x 的 rt*ct 维子矩阵
然后，将规范化后的图像的数据格式由适合显示图像的uint8格式转换为适合数值处理的double格式，再调用二维小波分解函数进行图像分解，最后为了清晰地显示分解图像的塔式结构，在图像的相应区域绘制若干分界线。具体程序如下：
y=mywavedec2(x,dim)
% 函数 MYWAVEDEC2() 对输入矩阵 x 进行
dim 层分解，得到相应的分解系数矩阵 y
% 输入参数：x ——
输入矩阵；
%&&&&&&&&&&
dim —— 分解层数。
% 输出参数：y ——
分解系数矩阵。
x=modmat(x,dim);&&&&&&&&&&&
% 首先规范化输入矩阵，使其行列数均能被 2^dim 整除
subplot(121);imshow(x);title('原始图像');&&
% 画出规范化后的源图像
[m,n]=size(x);&&&&&&&&&&&&&
% 求出规范化矩阵x的行列数
xd=double(x);&&&&&&&&&&&&&&
% 将矩阵x的数据格式转换为适合数值处理的double格式
for i=1:dim
xd=modmat(xd,1);
[dLL,dHL,dLH,dHH]=mydwt2(xd);&& %
矩阵小波分解
tmp=[dLL,dHL;dLH,dHH];&&&&&&&&&
% 将分解系数存入缓存矩阵
xd=dLL;&&&&&&&&&&&&&&&&&&&&&&&&
% 将缓存矩阵左上角部分的子矩阵作为下一层分解的源矩阵
[row,col]=size(tmp);&&&&&&&&&&&
% 求出缓存矩阵的行列数
y(1:row,1:col)=&&&&&&&&&&&&
% 将缓存矩阵存入输出矩阵的相应行列
yd=uint8(y);&&&&&&&&&&&
% 将输出矩阵的数据格式转换为适合显示图像的uint8格式
i=1:dim&&&&&&&&&&&&
% 对矩阵 yd 进行分界线处理，画出分解图像的分界线
m=m-mod(m,2);
n=n-mod(n,2);
yd(m/2,1:n)=255;
yd(1:m,n/2)=255;
m=m/2;n=n/2;
subplot(122);imshow(yd);title([
num2str(dim) ' 维小波分解图像']);
以下是程序的运行结果。
&&&&&&&上述的图像分解程序，其输出数据是double格式的，以便作为重构程序的输入。
接下来我们讨论图像的重构。我编写的重构程序中，为了比较分解图像和重构图像，首先绘出经过小波分解的图像，然后再进行重构。在绘制分解图像和重构图像的过程中，要注意数据格式的转换。
y=mywaverec2(x,dim)
%&函数&MYWAVEREC2()&对输入的分解系数矩阵x进行&dim&层重构，得到重构矩阵&y
%&输入参数：x&——分解系数矩阵；
%&&&&&&&&&&
dim&——重构层数；
%&输出参数：y&——重构矩阵。
%&绘制分解图像
xd=uint8(x);&&&&&&&&&&&
%&将输入矩阵的数据格式转换为适合显示图像的uint8格式
[m,n]=size(x);&&&&&&&&
&%&求出输入矩阵的行列数
i=1:dim&&&&&&&&&&&&
%&对转换矩阵xd进行分界线处理
m=m-mod(m,2);
n=n-mod(n,2);
xd(m/2,1:n)=255;
xd(1:m,n/2)=255;
m=m/2;n=n/2;
subplot(121);imshow(xd);title([ num2str(dim)
'&层小波分解图像']);&%&画出带有分界线的分解图像
%&重构图像
xr=double(x);&&&&&&&&&&
%&将输入矩阵的数据格式转换回适合数值处理的double格式
[row,col]=size(xr);&&&&
%&求出转换矩阵xr的行列数
i=dim:-1:1&&&&&&&&&
%&重构次序是从内层往外层进行，所以先抽取矩阵&xr&的最内层分解矩阵进行重构
&&&&tmp=xr(1:floor(row/2^(i-1)),1:floor(col/2^(i-1)));&&&&&&
%&重构的内层矩阵的行列数均为矩阵xr的2^(i-1)
[rt1,ct1]=size(tmp);&&&&&&&&&&&&&&&&&&&&&&&&
%&读取待重构矩阵&tmp&的行列数
rt=rt1-mod(rt1,2);ct=ct1-mod(ct1,2);
rLL=tmp(1:rt/2,1:ct/2);&&&&&&&&&&&&&&&&&&&
%&将待重构矩阵&tmp&分解为四个部分
rHL=tmp(1:rt/2,ct/2+1:ct);
rLH=tmp(rt/2+1:rt,1:ct/2);
rHH=tmp(rt/2+1:rt,ct/2+1:ct);
tmp(1:rt,1:ct)=myidwt2(rLL,rHL,rLH,rHH);&&&&&&&&&&&&&
%&将重构结果返回到矩阵&tmp
xr(1:rt1,1:ct1)=&&&&&&
%&把矩阵&tmp&的数据返回到矩阵&xr&的相应区域，准备下一个外层的重构
y=&&&&&&&&&&&&&&&&&&&
%&重构结束后得到的矩阵xr即为输出矩阵&y
yu=uint8(xr);&&&&&&&&&&&
%&将矩阵xr的数据格式转换为适合显示图像的uint8格式
subplot(122);imshow(yu);title('小波重构图像');
已投稿到：
以上网友发言只代表其个人观点，不代表新浪网的观点或立场。}