基于NSST图像特征的多聚焦图像融合方法

张耀军，栗 磊，吴桂玲

(信阳农林学院计算机系，河南信阳464000)

多聚焦图像融合是指将原本针对同一场景的多幅聚焦于不同目标区域中的源图像中的清晰聚焦区域加以提取，并注入到一幅图像中的过程。该工作有助于提升图像的信息量，为后续的图像处理工作奠定良好基础。

目前，主流的多聚焦图像融合方法主要基于变换域理论而展开。离散小波变换(Discrete Wavelet Transform，DWT)较之前的金字塔变换在图像信息捕捉能力上已经有了较大的提升，但该理论仅限于捕捉奇异“点”信息，对二维的奇异“线”信息表达能力不足。为此，国内外专家学者针对DWT理论进行了探索和研究，先后提出了复小波变换(Complex Wavelet Transform，CWT)［1］、提升小波变换(Lifting Wavelet Transform，LWT)［2-3］、脊波变换(Ridgelet Transform，RT)［4］、Bandelet变换［5］、曲波变换［6］、轮廓波变换(Contourlet Transform，CT)［7-8］等理论，虽然这些理论在一定程度上改善了图像融合效果，但平移不变性的缺失导致了最终融合结果图像中出现了“振铃”效应，严重影响了视觉效果。为此，Cunha等人对CT理论进行进一步改进，摒弃了其中的下采样机制，构造了非下采样轮廓波变换(Non-Subsampled Contourlet Transform，NSCT)理论［9］。同CT理论相比，NSCT拥有更敏锐的图像细节和边缘信息捕捉能力，因此近年来，NSCT已被广泛用于图像融合领域［10-11］。相比以往经典的多尺度几何分析理论，NSCT的图像融合能力得以大幅提升，然而过高的计算复杂度大大制约了该理论的应用前景，尤其是诸如军事、医学等对实时性要求较高的实用领域。

2008年，Easley等人［12］构造并提出了非下采样剪切波变换(Non-Subsampled Shearlet Transform，NSST)理论，不同于NSCT理论的非下采样方向滤波器分解机制，NSST直接采取将剪切滤波器从伪极化网络系统映射到笛卡尔坐标系统，因此，NSST理论能够在保证具有良好图像处理性能的同时，大幅降低计算复杂度，在很大程度上克服了NSCT理论自身的缺陷，使其能被广泛应用于实时性领域成为可能。目前，该理论仍主要处于探索和发展阶段，涉及到图像融合领域的文献较为有限。文献［13-15］将NSST理论分别与其他领域内的模型相结合，取得了优良的图像融合效果。

在此背景下，本文拟仅依托NSST理论和待融合源图像自身的特征信息，提出一种基于NSST图像特征的多聚焦图像融合方法。由于低频子带图像包含了原始图像的主体信息，因此尽可能地将低频子带图像中的清晰区域加以提取，必将有助于提升后期最终融合图像的视觉效果，本文将选取能够描述图像清晰度水平的空间频率和局部能量作为依据，确定融合后的低频子带系数;另一方面，高频子带图像中主要保存了原始图像中的细节和边缘信息，因此可以直接利用边缘检测算子直接加以提取。仿真实验部分将本文方法与目前几种经典的方法进行了比较，结果表明本文方法可以在保证结果图像清晰度水平的前提下具备更低的计算复杂度。

1 NSST基本模型

当维数n=2时，具有合成膨胀的仿射系统表达式为

式中:ψ∈L2(R2)，A和 B均为 2×2的可逆矩阵，且。如果MAB(ψ)满足紧框架条件，则MAB(ψ)的元素称为合成小波。A称为各向异性膨胀矩阵，Ai与尺度变换相关联;B为剪切矩阵，Bj与保持面积不变的几何变换相关联。A和B的表达式分别为A=［a 0 0 a1/2］，B=［1 s0 1］。通常取 a=4，s=1，即 A=［4 0 0 2］，B=［1 1 0 1］。此时，该合成小波即为剪切波［12］。

NSST作为一种冗余变换，其分解过程分为多尺度分解和多方向分解两部分。其中，多尺度分解仍沿用非下采样金字塔分解方式加以实现，源图像每经过一级分解便可产生1幅低频子带图像和1幅对应的高频子带图像，此后每一级分解都将只针对低频子带图像进行。因此，在经过k级多尺度分解后，每一幅源图像共可衍生出1幅低频微观子图像和k幅高频微观子图像，这些子图像均与源图像具有相同尺寸。显然，NSST与NSCT在多尺度分解环节上是完全一致的。不同于NSCT，NSST的多方向分解是将标准剪切滤波器从伪极化网络系统映射到笛卡尔坐标系统，完全摒弃了下采样环节，从而实现了自身的平移不变性，有效地避免了类似于离散小波变换、轮廓波变换结果中的振铃效应。当NSST对某尺度子带图像进行l级非下采样剪切波方向分解时，可得到2l个与源图像具有相同尺寸大小的方向子带图像，从而实现频域中更为精确的方向分解。图1给出了图像Zoneplate的NSST分解示意图，其中，尺度分解数为2，每一级尺度对应的方向分解数均为4。

图1 Zoneplate图像的非下采样剪切波变换分解示意

显然，图1b反映了原始图像中的主体信息(中心部位的同心圆)，而图1c、图1d分别反映了不同分解尺度、不同方向下原始图像中的细节信息。

2 基于NSST图像特征的多聚焦图像融合框架

相比以往以NSCT为代表多尺度几何分析工具，NSST在具有平移不变性特征的同时具有更优越的图像信息捕捉能力。待融合源图像在经过NSST处理后，可以分别得到一幅低频子带图像和一系列高频子带图像。低频子带图像侧重反映原始图像中的主体和轮廓信息，含有绝大部分的特征，而高频子带图像主要包含了原始图像中的细节以及多个目标之间的边缘信息。因此，在图像融合过程中，针对低频子带图像应着眼于清晰轮廓框架信息的捕捉与提取，即尽可能地将其中的清晰区域提取出来，而针对高频子带图像则需要着重捕捉其中的边缘和细节信息，对低频提取信息加以有效补充。为此，本文选取了能够表达图像清晰度水平的空间频率和局部能量指标作为依据，针对源图像的若干幅低频子带图像系数加以选取。此外，选择Canny算子作为高频子带图像中边缘细节信息的提取工具。图2给出了本文的图像融合框架。

图2 基于NSST图像特征的多聚焦图像融合框架

其中，A，B，F分别表示两幅已配准源图像以及最终融合图像;X_low，X_high分别表示源图像经NSST分解后得到的低频、高频子带图像;low，high为融合后的低频、高频子带图像。

2．1 低频微观子带图像融合

空间频率［16］指数通过考查图像中各像素点(i，j)与相邻行像素点(i-1，j)以及相邻列像素点(i，j-1)之间灰度值的差异，从而刻画整幅图像的梯度特征，目前空间频率已被作为评价指标用于客观评价图像融合效果。然而，经典的空间频率指标只衡量了横、纵两个方向上的梯度差异，对对角线方向上的梯度差异并未加以刻画和描述。为此，文献［17］提出了一种新的空间频率指标，即

其中，RF，CF，MDF和 SDF分别代表水平、垂直、45°和135°对角线方向上的梯度水平。SF通过计算上述4个方向上的平方和水平的平方根，从而综合描述单个像素点的综合空间频率水平。

以尺寸为M×N的图像I(i，j)为例给出4个方向上的梯度水平表达式

其中，ωd表示距离权值，本文ωd取0.707。

显然，SF数值越大表示对应像素点的梯度水平越高，与周围像素点的差异越明显、越活泼突出，通常情况下也蕴含着更为丰富的信息。

此外，局部能量也从另一角度对像素点的特征进行了描述，局部能量的数学表达式为

式中:(i，j)为图像像素点;c表示子带图像系数;N×N为局部邻域尺寸。N的尺寸取决于图像的纹理水平，若纹理较为粗糙则适当减小N值，反之，若纹理较为平坦则适当增大N值。

低频子带图像系数的选择方法具体如下

2．2 高频微观子带图像融合

边缘检测算子可以用于对高频子带图像中的边缘细节信息进行捕捉。Canny算子基于最优化思想展开，该算子采用高斯函数对图像作平滑处理，具有较强的抑制噪声能力，此外，Canny算子采用双阈值算法检测和接连边缘，采用的多尺度检测和方向性搜索能力比其余算子具有明显优势。因此，本文采用Canny算子对高频子带图像中的细节边缘信息进行提取。

Canny算子求解边缘点的算法为:

1)高斯滤波器平滑图像处理

使用高斯函数

则

其中，p(x，y)为原始高频子带图像;q(x，y)为平滑处理后的子带图像。

2)使用一阶差分对平滑后的子带图像进行处理，求出偏导阵列

其中，M［x，y］和 θ［x，y］分别表示图像的边缘强度和边缘方向。M［x，y］取得最优值时对应的方向角 θ［x，y］即为边缘的走势方向。

3)用双阈值算法检测和连接边缘

对非极大值抑制图像施加双阈值 λ1和 λ2，且 λ1＜ λ2。将梯度值小于λ1的像素灰度值置为0，获取图像Ⅰ。同样地，将梯度值小于λ2的像素灰度值置为0，获取图像Ⅱ。紧接着，以图像Ⅱ为基础，使用图像Ⅰ对其进行边缘和细节部分的补充。

图3给出了一对多聚焦源图像的Canny算子处理结果图。

图3 高频子带图像Canny算子处理结果

显然，原始图像中的边缘细节信息在Canny算子的作用下得到了很好地保存。譬如右聚焦图像中右侧闹钟的数字信息较为清晰，图3c中对图3a对应的高频子带图像进行了Canny算子处理，实验结果表明源图像中的细节信息保留较为完好;左聚焦图像中左侧闹钟的商标信息较清晰，图3d中的结果来源于图3b对应的高频子带图像的Canny算子处理结果，如图3中的矩形框标识区域所示。

3 仿真实验

为了验证方法的实际性能，本节将采用MATLAB软件，仿真、比较与分析包括本文方法在内的几种经典融合方法的融合效果。

3．1 实验条件

实验仿真平台为一台 PC，配置情况为 Win 7系统，CPU2.1 GHz，2 Gbyte内存，500 Gbyte硬盘，软件环境为MATLAB2012。为了充分验证相关方法的融合性能，本节选取了三组多聚焦待融合源图像——Clock、Puma、Eyeball用于融合仿真实验，3组源图像为256级灰度、尺寸大小为512×512，图像均可以从www.imagefusion.org网址加以获取。

此外，还选取了5种方法与文中方法进行比较，分别为:基于剪切波变换的图像融合方法(方法1)［18］、基于NSCT与脉冲耦合神经网络的图像融合方法(方法2)［19］、基于非负矩阵分解和新轮廓波变换融合方法(方法3)［20］、基于NSCT与非负矩阵分解的图像融合方法(方法4)［21］以及基于NSCT的图像融合方法(方法5)［22］。为了保证融合效果的公正客观性，方法1～方法5中的参数仍沿用相关文献中的设置模式。在本文方法中，邻域大小为5×5，多尺度分解级数为3级，方向分解数从粗至细依次为8，16，16。

最终融合结果将采取主观视觉评价以及客观指标值评价两种方式。其中，将选取信息熵(IE)［17］、互信息量(MI)［17］以及空间频率(SF)［19］作为评价指标。IE数值表征图像中包含信息量的大小，数值越大，融合图像信息越丰富。MI数值反映了最终融合图像与源图像间的交互信息，MI数值越大，表征融合图像中涉及到源图像的原始信息越丰富。SF可以在一定程度上反映融合图像中各像素点的活跃程度，SF数值越大，代表融合图像包含的显著信息越丰富。IE，MI，SF的详细内容请参见文献［17，19］。

3．2 实验结果分析

3组仿真源图像如图4所示，其中，Clock图像中的聚焦区域分别集中于左、右两部分。图4a中，右侧闹钟处于聚焦区域，左侧闹钟处于离焦区域，相反的情况发生在图4b中;Puma图像中的聚焦区域集中于上、下两部分。图4c中，图像的下半区域尤其是两条腿部的区域较为清晰，但豹头以及豹身部分较模糊;相反情况出现在图4d中;Eyeball源图像的情况与Clock图像类似，图4e中左上角区域的一根血管非常清晰，但右侧区域较模糊，而图4f中的右侧血管区域较为清晰，但左上方的血管视觉效果不佳。

图4 仿真源图像

因此，尽可能地将源图像中的聚焦区域进行提取并注入到一幅图像中无疑将会大幅提升源图像的视觉效果。图5～图7给出了包括本文方法在内的6种方法的仿真结果图。

图5 Clock图像融合效果图

图6 Puma图像融合效果图

图7 Eyeball图像融合效果图

从视觉角度来看，相比其他5种方法，基于方法1的融合效果图像整体亮度水平较低，视觉效果水平不佳，例如图5a右侧闹钟的上方甚至出现了一定的“毛刺”虚假信息(如箭头所示);方法3～方法5的融合图像虽然亮度水平适中，但针对图像的细节信息捕捉能力较弱，譬如图5c中右侧闹钟的刻度较为模糊(如箭头所示)，图6c，6d，6e的豹头花斑也不清晰(如箭头所示)，图7c～7e中左上方血管信息几乎完全丢失(如箭头所示)。相比之下，方法2和本文方法对原始图像的主体和细节边缘信息进行了较好地提取和表达，图像整体亮度和对比度水平也更为适宜，二者的图像融合性能较为接近。然而，通过仔细比较不难发现，本文方法在图像清晰度水平上优势更为明显，譬如在图7中，基于本文方法的融合图像中左上方血管信息较方法2更为完整和清晰(如箭头所示)。

表1针对3组仿真实验列出了6种方法的比较结果。

从表1中不难发现，本文提出的方法较其他5种融合方法拥有更理想的客观指标结果。在第一组实验中，本文方法效果图对应的MI、SF指标值均为最优，IE数值也位于前列。第二组实验结果表明该组实验中基于本文方法的实验仿真图无论在信息量、互信息量和空间频率指标上均占有优势。第三组实验中，本文方法MI、SF指标值为最优，IE数值略低于方法1，这是因为基于方法1的融合图像中含有大量的噪声虚假信息，这些信息的存在客观上造成了相邻像素点间灰度值差异的增大，因此造成了IE数值也随之消极增长。很明显，表1中的客观数据与直观视觉效果在很大程度上是吻合的。

表1 6种融合方法的客观指标值比较

3．3 实验讨论

以上3组仿真实验分别从主观视觉效果和客观评价指标2个方面验证了本文方法的合理有效性。本节将对实验结果做进一步分析。一方面，由于方法2无论在直观视觉效果还是客观指标上均优于另外4种用于比较的方法，因此，本节选取方法2与本文方法做进一步深入比较与分析。为了更直观地观察方法2和本文方法在直观视觉效果层次的差异，以图7b和图7d为例对原结果图像中左上方血管部位作局部放大，得到图8的结果。

图8 Eyeball仿真结果局部放大图

显然，从图8所示的放大图像中可以发现，基于方法2的结果图像中血管清晰程度要弱于基于本文方法得出的结果图像，进一步验证了本文方法在图像信息捕捉能力方面的优势。

另一方面，计算代价也是衡量算法性能的重要指标之一。由于方法2与本文方法融合效果较为接近，且NSCT与NSST也同属近期极具代表性的非下采样理论多尺度几何分析体系，因此本节仅对这两种方法的算法复杂度加以研究。为了尽可能保证比较的客观公正性，本节选取第二组仿真实验为例，对每种方法均运行3次给出算法运行时间，如表2所示。

表2 算法运行时间比较s

由表2数据可知，本文方法在算法运行代价方面较方法2同样具有更为明显的优势，以期能够适用于一些实时性要求较高的应用领域。

不难看出，本文方法不仅具有更理想的图像融合效果，而且在算法运行效率上同样占有较为明显的优势。

4 结论

本文提出一种基于NSST图像特征的多聚焦图像融合方法。相比以往出现的相关方法，本文充分利用了NSST理论在图像信息捕捉以及计算复杂度方面的优势，提出依据图像自身特征解决NSST分解后的低频、高频子带图像的融合问题。仿真实验分别从主观视觉效果、客观评价指标值、仿真结果图像局部放大比较以及运行时间4个角度将本文方法与目前几种经典的方法进行了比较，结果表明本文方法可以在保证结果图像融合效果的前提下具备较低的计算复杂度，在一定程度上可以满足一些实时性要求较高的应用领域的实际需求。方法的进一步优化和性能提升是笔者未来努力的方向和研究重点。

［1］ SINGH R，KHARE A.Fusion of multimodal medical images using Daubechies complex wavelet transform-a multiresolution approach［J］.Information Fusion，2014，19(1):49-60.

［2］高颖，王阿敏，支朋飞，等.基于区域分割与提升小波变换的图像融合算法［J］.西北工业大学学报，2014，32(4):569-575.

［3］杨艳春，党建武，王阳萍.基于提升小波变换与自适应PCNN的医学图像融合方法［J］.计算机辅助设计与图形学学报，2012，24(4):494-499.

［4］赵子沂，郑永果.基于脊波的多光谱和全色图像融合方法研究［J］. 计算机工程与应用，2012，48(15):164-167.

［5］ 朱卫东，李全海，徐克科，等.基于二代Bandelet和主成分变换的高光谱遥感图像融合［J］.同济大学学报:自然科学版，2011，39(7):1068-1073.

［6］杨扬，戴明，周箩鱼.基于均匀离散曲波变换的多聚焦图像融合［J］. 红外与激光工程，2013，44(9):2547-2552.

［7］ DARWISH SM.Multi-level fuzzy contourlet-based image fusion for medical applications［J］.IET Image Processing，2013，7(7):694-700.

［8］ SAEEDI J，FAEZ K.A new pan-sharpening method using multiobjective particle swarm optimization and the shiftable contourlet trans-form［J］.Journal of Photogrammetry and Remote Sensing，2011，66(3):365-381.

［9］ CUNHA A L，ZHOU J P，DO M N.The nonsubsampled contourlet transform:theory，design，and application［J］.IEEE Trans.Image Processing，2006，15(10):3089-3101.

［10］ BHATNAGAR G，WU Q M J，ZHENG L.Directive contrast based multimodal medical image fusion in NSCT domain［J］.IEEE Trans.Multimedia，2013，15(5):1014-1024.

［11］ KONG W W，LEI Y J，LEI Y，et al.Image fusion technique based on non-subsampled contourlet transform and adaptive unitfast-linking pulse-coupled neural network［J］.IET Image Processing，2011，5(2):113-121.

［12］ GAOGR，XU L P，FENGD Z.Multi-focus image fusion based on non-subsampled shearlet transform［J］.IET Image Processing，2013，7(6):633-639.

［13］ KONG W W，LEI Y，ZHAO H X.Adaptive fusion method of visible light and infrared images based on non-subsampled shearlet transform and fast non-negative matrix factorization［J］.Infrared Physics＆ Technology，2014(67):161-172.

［14］ KONG W W，ZHANG L J，LEIY.Novel fusion method for visible light and infrared images based on NSST-SF-PCNN［J］.Infrared Physics＆ Technology，2014，65:103-112.

［15］ ESKICIOGLU A M，FISHER PS.Image quality measures and their performance［J］.IEEE Trans.Communications，1995，43(12):1959-2965.

［16］ ZHENG Y，ESSOCK E A，HANSEN B C，et al.A new metric based on extended spatial frequency and its application to DWT based fusion algorithms［J］.Information Fusion，2007，8(2):177-192.

［17］ MIAO Q G，SHI C，XU P F，et al.A novel algorithm of image fusion using shearlets［J］.Optics Communication，2011，284(6):1540-1547.

［18］ QU X B，YAN J W，XIAO H Z，et al.Image fusion algorithm based on spatial frequency-motivated pulse coupled neural networks in nonsubsampled contourlet transform domain［J］.Acta Automatica Sinica，2008，34(12):1508-1514.

［19］王斐，梁晓庚，崔彦凯，等.非负矩阵分解和新轮廓波变换的图像融合方法［J］.计算机工程与应用，2013，49(5):150-153.

［20］ KONG W W，LEI Y J，LEI Y，et al.Technique for image fusion based on non-subsampled contourlet transform domain improved NMF［J］.Science in China Series F-Information Sciences，2010，53(12):2429-2440.

［21］ Adu J H，GAN J H，WANG Y，et al.Image fusion based on nonsubsampled contourlet transform for infrared and visible light image［J］.Infrared Physics and Technology，2013，61(1):94-100.

［22］ LIU Z，BLASCH E，XUE Z Y，et al.Fusion algorithms for context enhancement in night vision:a comparative study［J］.IEEE Trans.Pattern Analysis and Machine Intelligence，2012，34(1):94-109.