基于Matlab的图像融合研究设计
基于MATLAB的图像融合系统仿真[GUI界面，多方法]
一、课题开发背景
1.1 图像融合的定义
数字图像融合（Digital Image Fusion）是一种基于多源图像是主要研究对象的数据整合技术。它通过将来自不同模式或同一设备在同一场景下不同时间获取的所有图象集合并形成一个综合图象来实现信息整合过程。由于多种模式下的成像原理和使用的电磁波波长存在差异导致各来源获取到的信息存在冗余性与互补性经过特定算法处理后形成的合成图象能够更加全面且精确地反映研究对象的本质特征其显著的特点使其在多个领域得到了广泛应用
1.2 手动配准与图象融合
注
图1-1 cp2tform函数能修正的6种几何变形
基于图像配准的结果，在存在重叠区域之间差异的情况下，在直接叠加的方式下会导致明显的接缝问题。因此必须采取技术手段来修正待拼接区域的颜色值分布情况，并实现平滑过渡以达到无缝合成的目标。传统的融合方法主要在时域采用算术运算方法，并未考虑处理过程中的频域变化特性。从数学意义上讲，在理论上消除拼接缝等价于实现两个颜色或灰度曲面的平滑连接；然而实际上由于图片的空间模糊特性所导致的光滑化效果与预期相反。
二、 图像融合方法
迄今为止，在像素级和特征级上开展的数据融合工作主要集中在特定领域内。通常采用的方法包括HIS综合法、基于KL变换的数据整合技术、高通滤波数据结合策略、基于小波变换的数据混合模式以及多分辨率金字塔数据处理框架等多种形式。以下将简要介绍几种典型的方法。
1． HIS融合法
HIS综合法在多源图像像素级融合领域具有广泛应用，在实际应用中通常采用以下流程进行处理：首先选取一个低分辨率的三波段图像与一个高分辨率的单波段图像作为输入源图像进行融合处理。具体而言，在该方法中采用HIS变换将三个波段的低分辨率数据转换至HIS空间，并对对比度拉伸处理后的新图像具有与目标空间亮度分量一致的均值和方差特征；随后将经亮度分量替换后的新图像通过反向变换恢复至原始空间。
2． KL变换融合法
KL变换融合法也被认为是主成分分析法。类似于HIS变换法的方式，在该方法中利用低分辨度图像（包含三个或更多波段）作为输入分量展开主元分析；随后将高于分辨度图像经过调整使其均值及方差与第一主元分量相一致，并用调整后高于分辨度图像替代主元转换的第一分量来进行反向转换过程。通过融合处理后获得的新数据集能够整合源图像中的高空间和高光谱细节特征，并且成功保留了原始图像中的高频细节信息；这样一来，在融合后的图像中能够更加清晰地展现目标细节特征，并且光谱信息也变得更加丰富。
3． 高通滤波融合法
该方法通过高通滤波器提取高分辨率图像的边缘细节，并将其整合至低分辨率的高光谱图像中。随后，利用高通滤波器分离并获取高分辨率图像中的高频细节特征。这些高频分量被有效整合至低分辨率图象中后，在生成的融合图象中显著提升了高频特征信息的表现能力。
4． 小波变换融合法
基于离散小波变换的方法中，我们首先对包含N幅待融合的图像进行处理，在各个层级上对来自原始集合中N幅图像所分割出的M个子图像进行综合处理。随后，在获得所有M个层级上的融合结果后层之后,通过逆向小波变换步骤整合各层级处理结果,从而最终获得完整的图像融合输出
三、图像融合步骤
目前国内外已有大量关于图像融合技术的研究报道。无论采用何种技术手段，在多幅图像中对应点必须精确对齐。根据研究对象及目标的不同,图象融合的方法也相应地多种多样.其主要步骤归纳如下:
(1) 预处理工作：去除噪声并进行图像增强等预处理操作，并确保数据格式统一、尺寸一致以及分辨率相同。完成基于序列断层扫描的数据的三维重建与可视化展示，并根据目标特性建立相应的数学建模方案
图3-1 图像融合步骤示意图
(2) 划分目标区域并确定配准关键点:对于二维或三维场景中,需要将目标物体或兴趣区域进行划分.这些关键点应基于同一物理标记在两个图像中的对应位置来确定,其中这些关键点可为人工标注或基于人体解剖学的自然标志.
基于特征点的数据配准：可被视为两个数据集间的线性或非线性变换，并使经过该变换后的两组数据在满足特定准则下达到误差最小化的目标。
构建图像融合：两幅经过配准的图像在一个统一的空间坐标系中将各自的重要信息整合为二维图或三维模型；
(5) 参数提取:从融合图像中提取和测量特征参数，定性、定量分析。
四、各算法程序源码
4.1 一般方法
基于一个数学模型将来自不同传感器的多幅图象结合成满足特定应用需求的一幅综合图象的过程称为图象融合技术。这种基本的图象融合方法并未对参与融合的原始图象进行任何变换或分解处理工作；而是直接针对各对应像素分别实施选择操作、算术平均运算、加权平均计算或其他数学运算等基础操作之后；最终生成一张经过综合处理后的融合图象。
对于图像融合的对象而言可划分为两大类即多光谱图象通常表现为RGB彩色图象与灰度图象间的融合以及灰度图象间的融合灰度图象间的融合大致可分为三大类包括基本融合方法对比度和比率金字塔法拉普拉斯金字塔法以及小波变换下的多分辨率分析方法等它们均需首先构建输入图象的金字塔并依据特定特征选择参数生成相应的融合金字塔随后通过逆变换重建出最终的综合图象其中基于金字塔形分解算法的方法尽管具有较高的合成效果但也存在一些不足之处例如在处理具有较大差异区域时容易导致合成后的分块现象此外该算法在分解过程中未充分考虑空间方向选择性可能导致信息提取不够全面相比之下基于小波变换的方法则能更好地提取原始图象的空间定位信息及其细节特征这是因为小波变换所具有的正交性质能够有效避免因空间方向选择性带来的冗余问题同时其多分辨率特性使得不同层次的信息提取更为精确此外基于小波变换的方法还能够完全恢复信号能量避免因传统的金字塔分解过程而导致的信息丢失或引入人工干扰因素综上所述尽管基于金字塔形分解算法的方法在某些应用领域仍具一定的实用价值但其局限性已逐渐被新型的小波域融合方法所取代
4.2 PCA算法程序
基于主成分分析（PCA）的图像融合技术称作K-L变换法，在其几何意义上则将原始特征空间中的特征轴重新定向至与混合集群结构的主要方向平行的过程中形成了新的特征轴系统。在实际应用中，则是对原有各因素指标（其中一些指标间存在相关性）进行重新组合，在重组后的综合指标之间相互独立的情况下构建新的特征空间系统。仅需选取前几个主成分图即可充分提取原集群的关键信息，在冗余的信息得以去除的同时实现了对原数据本质的有效表征；该方法特别适用于处理包含多重相关因素的多源遥感数据场景，在信息融合效率方面表现出明显优势
基于PCA的融合算法程序：
function Y = fuse_pca(M1， M2)
%Y = fuse_pca(M1， M2) image fusion with PCA method
%
% M1 - input image #1
% M2 - input image #2
%
% Y - fused image
% (Oliver Rockinger 16.08.99)
% check inputs
[z1 s1] = size(M1);
[z2 s2] = size(M2);
if (z1 ~= z2) | (s1 ~= s2)
error('Input images are not of same size');
end;
% compute， select & normalize eigenvalues
[V， D] = eig(cov([M1(:) M2(:)]));
if (D(1，1) > D(2，2))
a = V(:，1)./sum(V(:，1));
else
a = V(:，2)./sum(V(:，2));
end;
% and fuse
Y = a(1)*M1+a(2)*M2;
AI写代码4.3 金字塔（Pyramid）算法程序
Pyramid-based image fusion technique: representing images in a pyramid structure is a straightforward and efficient approach. In simpler terms, the pyramid image fusion technique involves constructing multi-resolution representations for each input image. These representations are then merged across corresponding pyramid layers using specific fusion rules to produce a composite pyramid. This composite pyramid can subsequently be reconstructed into an image, yielding the fused result. The pyramid can be categorized into types such as Laplacian pyramids, Gaussian pyramids, gradient pyramids, and morphological pyramids.
基于FSD Pyramid的图像融合算法程序：
function Y = fuse_fsd(M1， M2， zt， ap， mp)
%Y = fuse_fsd(M1， M2， zt， ap， mp) image fusion with fsd pyramid
%
% M1 - input image A
% M2 - input image B
% zt - maximum decomposition level
% ap - coefficient selection highpass (see selc.m)
% mp - coefficient selection base image (see selb.m)
%
% Y - fused image
% (Oliver Rockinger 16.08.99)
% check inputs
[z1 s1] = size(M1);
[z2 s2] = size(M2);
if (z1 ~= z2) | (s1 ~= s2)
error('Input images are not of same size');
end;
% define filter
w = [1 4 6 4 1] / 16;
% cells for selected images
E = cell(1，zt);
% loop over decomposition depth -> analysis
for i1 = 1:zt
% calculate and store actual image size
[z s] = size(M1);
zl(i1) = z; sl(i1) = s;
% check if image expansion necessary
if (floor(z/2) ~= z/2)， ew(1) = 1; else， ew(1) = 0; end;
if (floor(s/2) ~= s/2)， ew(2) = 1; else， ew(2) = 0; end;
% perform expansion if necessary
if (any(ew))
M1 = adb(M1，ew);
M2 = adb(M2，ew);
end;
% perform filtering
G1 = conv2(conv2(es2(M1，2)， w， 'valid')，w'， 'valid');
G2 = conv2(conv2(es2(M2，2)， w， 'valid')，w'， 'valid');
% select coefficients and store them
E(i1) = {selc(M1-G1， M2-G2， ap)};
% decimate
M1 = dec2(G1);
M2 = dec2(G2);
end;
% select base coefficients of last decompostion stage
M1 = selb(M1，M2，mp);
% loop over decomposition depth -> synthesis
for i1 = zt:-1:1
% undecimate and interpolate
M1T = conv2(conv2(es2(undec2(M1)， 2)， 2*w， 'valid')， 2*w'， 'valid');
% add coefficients
M1 = M1T + E{i1};
% select valid image region
M1 = M1(1:zl(i1)，1:sl(i1));
end;
% copy image
Y = M1;
AI写代码4.4 小波变换（DWT）算法程序
在众多图像融合技术中
以小波变换理论为基础的算法的具体实现流程
function Y = fuse_sih(M1， M2， zt， ap， mp)
%Y = fuse_sih(M1， M2， zt， ap， mp) image fusion with SIDWT， Wavelet is Haar
%
% M1 - input image A
% M2 - input image B
% zt - maximum decomposition level
% ap - coefficient selection highpass (see selc.m)
% mp - coefficient selection base image (see selb.m)
%
% Y - fused image
% (Oliver Rockinger 16.08.99)
% check inputs
[z1 s1] = size(M1);
[z2 s2] = size(M2);
if (z1 ~= z2) | (s1 ~= s2)
error('Input images are not of same size');
end;
% cells for selected images
E = cell(3，zt);
% loop over decomposition depth -> analysis
for i1 = 1:zt
% calculate and store actual image size
[z s] = size(M1);
zl(i1) = z; sl(i1) = s;
% define actual filters (inserting zeros between coefficients)
h1 = [zeros(1，floor(2^(i1-2)))， 0.5， zeros(1，floor(2^(i1-1)-1))， 0.5， zeros(1，max([floor(2^(i1-2))，1]))];
g1 = [zeros(1，floor(2^(i1-2)))， 0.5， zeros(1，floor(2^(i1-1)-1))， -0.5， zeros(1，max([floor(2^(i1-2))，1]))];
fh = floor(length(h1)/2);
% image A
Z1 = conv2(es(M1， fh， 1)， g1， 'valid');
A1 = conv2(es(Z1， fh， 2)， g1'，'valid');
A2 = conv2(es(Z1， fh， 2)， h1'，'valid');
Z1 = conv2(es(M1， fh， 1)， h1， 'valid');
A3 = conv2(es(Z1， fh， 2)， g1'，'valid');
A4 = conv2(es(Z1， fh， 2)， h1'，'valid');
% image B
Z1 = conv2(es(M2， fh， 1)， g1， 'valid');
B1 = conv2(es(Z1， fh， 2)， g1'，'valid');
B2 = conv2(es(Z1， fh， 2)， h1'，'valid');
Z1 = conv2(es(M2， fh， 1)， h1， 'valid');
B3 = conv2(es(Z1， fh， 2)， g1'，'valid');
B4 = conv2(es(Z1， fh， 2)， h1'，'valid');
% select coefficients and store them
E(1，i1) = {selc(A1， B1， ap)};
E(2，i1) = {selc(A2， B2， ap)};
E(3，i1) = {selc(A3， B3， ap)};
% copy input image for next decomposition stage
M1 = A4;
M2 = B4;
end;
% select base coefficients of last decompostion stage
A4 = selb(A4，B4，mp);
% loop over decomposition depth -> synthesis
for i1 = zt:-1:1
% define actual filters (inserting zeros between coefficients)
h2 = fliplr([zeros(1，floor(2^(i1-2)))， 0.5， zeros(1，floor(2^(i1-1)-1))， 0.5， zeros(1，max([floor(2^(i1-2))，1]))]);
g2 = fliplr([zeros(1，floor(2^(i1-2)))， 0.5， zeros(1，floor(2^(i1-1)-1))， -0.5， zeros(1，max([floor(2^(i1-2))，1]))]);
fh = floor(length(h2)/2);
% filter (rows)
A4 = conv2(es(A4， fh， 2)， h2'， 'valid');
A3 = conv2(es(E{3，i1}， fh， 2)， g2'， 'valid');
A2 = conv2(es(E{2，i1}， fh， 2)， h2'， 'valid');
A1 = conv2(es(E{1，i1}， fh， 2)， g2'， 'valid');
% filter (columns)
A4 = conv2(es(A4+A3， fh， 1)， h2， 'valid');
A2 = conv2(es(A2+A1， fh， 1)， g2， 'valid');
% add images
A4 = A4 + A2;
end;
% copy image
Y = A4;
AI写代码五、运行示意图
下面将本文提出的算法应用于多焦点图像融合过程。这一系列多焦点图像均是在同一场景下通过不同焦距拍摄所得，并且在这些图像中各处都有一定的清晰度特性。通过应用融合技术处理这些输入后，则能够生成一个所有目标均达到较高清晰度水平的输出画面（见图5-1）。在图5-1所示的画面中左侧区域呈现较高的清晰度，在图5-2所示的画面中右侧区域呈现较高的清晰度
图5-1 聚焦在左边的图像
图5-2 聚焦在右边的图像
采用主成分分析为基础以及多分辨率金字塔图像融合技术和小波变换方法所开发出的算法程序能够获得的结果
图5-3 基于PCA算法的融合图像
图5-4 基于金字塔图像融合算法的融合图像
图5-5 基于SIDWT小波变换的融合图像
从实验结果可以看出