基于Matlab的图像融合研究设计

小咪咪 2022-09-09 02:06 221阅读 0赞

**基于MATLAB的图像融合系统仿真\[GUI界面，多方法\]**

**一、课题开发背景**

**1.1 图像融合的定义**

数字图像融合（Digital Image Fusion）是以图像为主要研究内容的数据融合技术，是把多个不同模式的图像传感器获得的同一场景的多幅图像或同一传感器在不同时刻获得的同一场景的多幅图像合成为一幅图像的过程。由于不同模式的图像传感器的成像机理不同，工作电磁波的波长不同，所以不同图像传感器获得的同一场景的多幅图像之间具有信息的冗余性和互补性，经图像融合技术得到的合成图像则可以更全面、更精确地描述所研究的对象。正是由于这一特点，图像融合技术现已广泛地应用于军事、遥感、计算机视觉、医学图像处理等领域中。

**1.2 手动配准与图象融合**

图像融合包含图像配准和无缝合成两个部分.由于成像时受到各种变形因素的影响，得到的各幅图像间存在着相对的几何差异。图像配准是通过数学模拟来对图像间存在着的几何差异进行校正，把相邻两幅图像合成到同一坐标系下，并使得相同景物在不同的局部图像中对应起来，以便于图像无缝合成。本文采用Matlab中的cpselect、cp2tform函数完成几何配准。cpselect函数显示图像界面，手动在两幅图像的重叠部分选取配准控制点，Matlab自动进行亚像素分析，由cp2tform函数值正重叠部分的几何差异。Matlab中的cp2tform函数能修正6种变形，分别是图(b)：线性相似；图(c)：仿射；图(d)：投影；图(e)：多顶式；图(o：分段线性；图(g)：局部加权平均。前4种为全局变换，后两种为局部变换，如图I所示，图(a)为修正结果；图(b)至图(g)为对应的变形。

![v2-95e1127a4bdf44a173f551a779e9a4f2\_b.jpg][v2-95e1127a4bdf44a173f551a779e9a4f2_b.jpg]

图1-1 cp2tform函数能修正的6种几何变形

图像配准之后，由于图像重叠区域之间差异的存在，如果将图像象素简单叠加，拼按处就会出现明显的拼接缝，因此需要一种技术修正待拼接图像拼接缝附近的颜色值，使之平滑过渡，实现无缝合成。传统的融合方法多是在时间域对图像进行算术运算，没有考虑处理图像时其相应频率域的变化。从数学上讲，拼接缝的消除相当于图像颜色或灰度曲面的光滑连接，但实际上图像的拼接与曲面的光滑不同，图像颜色或灰度曲面的光滑表现为对图像的模糊化，从而导致图像模糊不清。

**二、 图像融合方法**

迄今为止，数据融合方法主要是在像元级和特征级上进行的。常用的融合方法有HIS融合法、KL变换融合法、高通滤波融合法、小波变换融合法、金字塔变换融合法、样条变换融合法等。下面简要介绍其中的几种方法。

1． HIS融合法

HIS融合法在多传感器影象像元融合方面应用较广，例如：一低分辨率三波段图象与一高分辨率单波段图象进行 融合处理。这种方法将三个波段的低分辨率的数据通过HIS变换转换到HIS空间，同时将单波段高分辨率图象进行对比度拉伸以使其灰度的均值与方差和HIS空间中亮度分量图象一致，然后将拉伸过的高分辨率图象作为新的亮度分量代入HIS反变换到原始空间中。这样获得的高分辨率彩色图象既具有较高空间分辨率，同时又具有与影象相同的色调和饱和度，有利于目视解译和计算机识别。

2． KL变换融合法

KL变换融合法又称为主成分分析法。与HIS变换法类似，它将低分辨率的图象（三个波段或更多）作为输入分量进行主成分分析，而将高分辨率图象拉伸使其具有于第一主成分相同的均值和方差，然后用拉伸后的高分辨率影象代替主成分变换的第一分量进行逆变换。高空间分辨率数据与高光谱分辨率数据通过融合得到的新的数据包含了源图象的高分辨率和高光谱分辨率特征，保留了原图象的高频信息。这样，融合图象上目标细部特征更加清晰，光谱信息更加丰富。

3． 高通滤波融合法

高通滤波融合法将高分辨率图象中的边缘信息提取出来，加入到低分辨率高光谱图象中。首先，通过高通滤波器提取高分辨率图象中的高频分量，然后将高通滤波结果加入到高光谱分辨率的图象中，形成高频特征信息突出的融合影象。

4． 小波变换融合法

利用离散的小波变换，将N幅待融合的图象的每一幅分解成M幅子图象，然后在每一级上对来自N幅待融合图象的M幅子图象进行融合，得到该级的融合图象。在得到所有M级的融合图象后，实施逆变换得到融合结果。

**三、图像融合步骤**

目前国内外己有大量图像融合技术的研究报道，不论应用何种技术方法，必须遵守的基本原则是两张或多张图像上对应的每一点都应对位准确。由于研究对象、目的不同，图像融合方法亦可多种多样，其主要步骤归纳如下:

(1) 预处理:对获取的两种图像数据进行去噪、增强等处理，统一数据格式、图像大小和分辨率。对序列断层图像作三维重建和显示，根据目标特点建立数学模型;

图3-1 图像融合步骤示意图

(2) 分割目标和选择配准特征点:在二维或三维情况下，对目标物或兴趣区进行分割。选取的特点应是同一物理标记在两个图像上的对应点，该物理标记可以是人工标记，也可以是人体解剖特征点;

(3 )利用特征点进行图像配准:可视作两个数据集间的线性或非线性变换，使变换后的两个数据集的误差达到某种准则的最小值;

(4) 融合图像创建:配准后的两种模式的图像在同一坐标系下将各自的有用信息融合表达成二维或三维图像;

(5) 参数提取:从融合图像中提取和测量特征参数，定性、定量分析。

**四、各算法程序源码**

**4.1 一般方法**

图像融合是通过一个数学模型把来自不同传感器的多幅图像综合成一幅满足特定应用需求的图像的过程，从而可以有效地把不同图像传感器的优点结合起来，提高对图像信息分析和提取的能力。简单的图像融合方法不对参加融合的源图像进行任何变换和分解，而是直接对源图像中的各对应像素分别进行选择、平均或加权平均、多元回归或其它数学运算等处理后，最终合成一幅融合图像。

对于图像融合的对象，可以分为两大类，即多光谱图像(通常为RGB彩色图像)与灰度图像之间的融合，以及灰度图像之间的融合。灰度图像之间的融合，在大体上可分为三大类。一类是简单融合方法，包括将空间对准的两幅图像直接求加权平均值。另一类方法是基于金字塔形分解和重建算法的融合方法，主要包括梯度金字塔法、对比度和比率金字塔法以及拉普拉斯金字塔法等，它们首先构造输入图像的金字塔，再按一定的特征选择方法取值形成融合金字塔，通过对金字塔实施逆变换进行图像重建，最终生成融合图像，它们的融合效果要远优于第一类方法，然而它也有很多不尽如人意的地方。还有一类方法就是近几年兴起的基于小波变换的图像融合方法，它通常采用多分辨分析和Mallat快速算法，通过在各层的特征域上进行有针对性的融合，比较容易提取原始图像的结构信息和细节信息，所以融合效果要好于基于金字塔形变换的图像融合法。这是因为小波变换更为紧凑;小波表达式提供了方向信息，而金字塔表达式未将空间方向选择性引入分解过程;由于可以选择正交小波核，因此不同分辨率包含的信息是唯一的，而金字塔分解在两个不同的尺度之间含有冗余，另外金字塔不同级的数据相关，很难知道两级之间的相似性是由于冗余还是图像本身的性质引起的;金字塔的重构过程可能具有不稳定性，特别是两幅图像存在明显差异区域时，融合图像会出现斑块，而小波变换图像融合则没有类似的问题。此外，小波变换具有完善重建能力保证了信号在分解重构过程中没有信息损失和信息冗余产生。

**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;

**4.3 金字塔（Pyramid）算法程序**

金字塔图像融合法:用金字塔在空间上表示图像是一种简单方便的方法。概括地说金字塔图像融合法就是将参加融合的每幅源图像作金字塔表示，将所有图像的金字塔表示在各相应层上以一定的融合规则融合，可得到合成的金字塔。将合成的金字塔，用金字塔生成的逆过程重构图像，则可得到融合图像。金字塔可分为:Laplacian金字塔、Gaussian金字塔、梯度金字塔、数学形态金字塔等。

基于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;

**4.4 小波变换（DWT）算法程序**

在众多的图像融合技术中，基于小波变换的图像融合方法已成为现今研究的一个热点。这类算法主要是利用人眼对局部对比度的变化比较敏感这一事实，根据一定的融合规则，在多幅原图像中选择出最显著的特征，例如边缘、线段等，并将这些特征保留在最终的合成图像中。在一幅图像的小波变换中，绝对值较大的小波系数对应于边缘这些较为显著的特征，所以大部分基于小波变换的图像融合算法主要研究如何选择合成图像中的小波系数，也就是三个方向上的高频系数，从而达到保留图像边缘的目的。虽然小波系数(高频系数)的选择对于保留图像的边缘等特征具有非常主要的作用，但尺度系数(低频系数)决定了图像的轮廓，正确地选择尺度系数对提高合成图像的视觉效果具有举足轻重的作用。

基于SIDWT（Shift Invariance Discrete Wavelet Transform）小波变换的算法程序：

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;

**五、运行示意图**

下面将本文的算法用于多聚焦图像的融合。多聚焦图像指的是对相同的场景用不同的焦距进行拍摄，得到镜头聚焦目标不同的多个图像。经过图像融合技术后，就可以得到一个所有目标都聚焦清晰的图像。图5-1中左边的目标较为清晰，图5-2中右边的目标较为清晰。

![v2-ed10ca758f7039ace91aac1eb9c57319\_b.jpg][v2-ed10ca758f7039ace91aac1eb9c57319_b.jpg]

图5-1 聚焦在左边的图像

![v2-15eaec1bb3b3dccd990a7fc7eccce15d\_b.jpg][v2-15eaec1bb3b3dccd990a7fc7eccce15d_b.jpg]

图5-2 聚焦在右边的图像

我们分别利用基于PCA的算法、金字塔图像融合法和小波变换法的算法程序得到的的融合图像结果，如图5-3、图5-4、图5-5所示

![v2-fe423dce18b3529eda02ac6d516e4fd2\_b.jpg][v2-fe423dce18b3529eda02ac6d516e4fd2_b.jpg]

图5-3 基于PCA算法的融合图像

![v2-cbf5ba57adc359961edd98c3331e1eb0\_b.jpg][v2-cbf5ba57adc359961edd98c3331e1eb0_b.jpg]

图5-4 基于金字塔图像融合算法的融合图像

![v2-e781c1f3f4b59ccb7c1dfcbe57f52366\_b.jpg][v2-e781c1f3f4b59ccb7c1dfcbe57f52366_b.jpg]

图5-5 基于SIDWT小波变换的融合图像

从实验结果可以看出，三种方案都可以得到较满意的视觉效果，消除了原图像的聚焦差异，提高了图像的清晰度，在合成图像中左、右两边的目标都很清晰。但通过比较分析，我们可以看出基于小波变换的融合图像（图5-5）最为清晰，所表现的图像细节效果最好，重影现象消除得最干净。图5-3的清晰度不够，而图5-4的细节表现力较弱，只有图5-5的边缘最清晰，重影消除地最干净，细节得到了最好地保留。

[v2-95e1127a4bdf44a173f551a779e9a4f2_b.jpg]: /images/20220829/5dc8f531642c4d1285b733aa818409cd.png
[v2-ed10ca758f7039ace91aac1eb9c57319_b.jpg]: /images/20220829/7a78c3d22b62452e8164b205a560a072.png
[v2-15eaec1bb3b3dccd990a7fc7eccce15d_b.jpg]: /images/20220829/3050a70b0ee54f3db2544f2c3f860c2f.png
[v2-fe423dce18b3529eda02ac6d516e4fd2_b.jpg]: /images/20220829/85e40d09c7034ffea7ce8d5fdebec41a.png
[v2-cbf5ba57adc359961edd98c3331e1eb0_b.jpg]: /images/20220829/641f12ccf2554722a32c99261db8ef1a.png
[v2-e781c1f3f4b59ccb7c1dfcbe57f52366_b.jpg]: /images/20220829/27cc471920784e30b7b7f7c0f743c0b4.png