摘 要：多模图像由于受到成像带宽的限制，往往在有些波段会受到较大噪声的干扰。对偶双边滤波是目前对多模图像去噪较为常规的方法，能同时对多幅图像进行联合去噪，但存在计算效率低下的问题。针对该问题提出一种基于交叉积分直方图的对偶双边滤波算法。该算法通过所构造的交叉积分直方图能够在常数时间里计算出滤波器的系数，并且算法的复杂度不受滤波半径的影响，因此有效提高了原始对偶双边滤波器的适用范围和计算效率。

关键词： 图像去噪；多模图像；联合双边滤波；积分直方图

0 引言

　　在图像处理与计算机视觉领域，图像去噪是一项重要的、基本的研究课题[1-2]。去除图像噪声能提高图像视觉质量，是对图像后续处理的前提，另外，图像去噪的方法在理论上与其他图像处理有着紧密的联系，研究图像去噪方法有一定的理论价值和实际意义。目前基于空间域的图像去噪分为局部方法和非局部方法[3]。

　　本文研究的是对多模图像的去噪，多模图像是指对同一目标在不同条件下获得的两幅或者两幅以上的图像。在对多模图像进行处理之前都必须进行图像预处理，如图像的配准。目前多模图像的柔性配准[4]是医学图像处理和遥感图像处理的一个热点，也是难点问题。在未来多模图像联合处理将会越来越多地应用于实际工程中。

　　目前的图像去噪方法中比较有代表性的一类是用高分辨率模态的图像来指导低分辨率的图像进行联合处理。TOMASI C等人[5-6]针对图像去噪问题提出的双边滤波算法，由于其简单易实现和无需迭代的特性在图像处理领域获得了广泛运用。PETSCHNIGG G等人[7]将双边滤波（Bilateral Filter，BF）扩展到了两张相关的图像，提出了联合双边滤波（Joint Bilateral Filter，JBF）的方法，用轮廓清晰的闪光图像来指导灰度信息保持较好的非闪光图像的去噪。Zhang Ke等人[8]又在此基础上提出了联合积分直方图方法，利用该方法能够在常数时间内执行联合双边滤波，有效降低了联合双边带滤波的计算复杂度。BENNETT E等人[9]针对包含了明显噪声的多模图像提出了对偶双边滤波（Dual Bilateral Filter，DBF），其几何权重使用欧式距离，灰度权值由所有分量共同计算得到，并采用统一的灰度权重。该算法能够在抑制噪声的同时对图像进行锐化。

　　本文主要考虑的是对存在明显噪声且已经配准好的多模图像进行交叉去噪。联合积分直方图（Joint Integral Histogram，JIH）[10]使得任意矩形区域的联合双边滤波结果能够在常数时间内计算出来。本文正是受此启发，对该算法进行了扩展，提出了交叉积分直方图（Cross Integral Histogram，CIH），并使用该数据结构对DBF算法进行了有效加速，使之能够在常数时间内得到去噪结果并且不受滤波块半径的影响。

1 相关工作

　　1.1 双边滤波

　　双边带滤波器由两个函数构成，一个函数gs（·）由几何空间距离决定滤波器系数，另一个函数gr（·）由像素差值决定滤波器系数。对于图像：

1.2 联合双边滤波

　　联合双边滤波（JBF）算法[7]是对双边滤波[5]算法的扩展。该算法在对一幅受噪声污染的图像I1进行滤波时所使用的滤波权重由另一幅不受噪声污染的图像I2来引导计算，公式如下：

其中，灰度滤波器的权重用无噪声图像的灰度进行计算。

1.3 对偶双边滤波

　　在多模图像处理领域，由于成像设备的限制不同模态的图像均可能被噪声污染。有些情况下在对一幅图像去噪时可能所参照的另一幅图像同样受到噪声污染。在这种情况下对偶双边滤波是较为广泛使用的一种方法。对偶滤波算法对于多模图像的每个分量采用统一的灰度测度权重，该权重由所有分量共同计算得到，公式如下：

　　其中，k=1，2代表对两个模态的图像分别进行去噪。

2 基于交叉积分直方图的双边滤波算法

　　本节对DBF算法使用交叉积分直方图进行加速，为了算法描述方便，首先考虑将式（3）中空间高斯函数代换为box函数的简单情况进行加速。后面将通过多个box空间滤波器的线性叠加来近似高斯滤波器。

　2.1 将空间高斯函数代换为box函数

　　为了简化计算先将空间距离的滤波函数由高斯函数近似为box函数。则对于图像中像k素点i的灰度值估计为：

　　对式（4）可以写成另外一种形式：

　　其中，s1、s2分别代表对图像I1和I2灰度量化后的灰度阶。相比于式（4）对于空间像素进行迭代累加，式（5）将迭代累加索引变换为图像灰度的量化阶。实际上，在分子中第二个累加和可以解释为对于给定的s1、s2，Lik δ（Ij1-s1）δ（Ij2-s2）代表对图像Ik的滤波块中，对满足特定条件的空间位置的灰度值进行累加。可以看到，该条件即在这些空间位置上I1和I2的灰度阶同时分别等于s1和s2。同样，分母中δ（Ij1-s1）δ（Ij2-s2）表示在空间滤波的块中，满足上述条件的空间位置的个数。

　　为了对式（5）加速，本文提出了交叉积分直方图，将图像从该点到图像原点的矩形区域中所有满足上述条件的灰度值累加结果保存在该数据结构中。式（5）分子部分的交叉积分直方图（CIH）定义为：

　　对算法的主要加速来自于对图像I1与图像I2构造交叉积分直方图，通过交叉积分直方图使得原始的乘积求和可以在常数时间得到并且算法的复杂度不受滤波半径的影响。

2.2 通过线性组合box函数来近似高斯函数

　　通常情况下空间滤波器设置为高斯函数，同参考文献[9]一样，本文使用二项式滤波器来近似高斯滤波器。而二项式滤波器可以被视为3个重叠的box滤波器的组合来进行近似。在本文中，3个box滤波器的半径分别为，相应的权重分别为={2，3，1}，因此式（10）可进一步写成：

　　其中，S（r，i）代表式（10）中的分子而W（r，i）代表分母，r为计算时的邻域半径。

3 实验与结果

　　在实验中，首先将所提出方法与对偶双边滤波方法（DBF）进行比较。如图1所示，第1列为从Brainweb上获得的模拟脑部核磁共振的多模图像“T1”和“T2”。第2列图像是对原始多模图像添加了高斯噪声后的图像，噪声图像的峰值信噪比为PSNR=20 dB。用DBF算法的去噪结果在第3列，所提出方法的去噪结果在第4列。第1行表示对“T1”的去噪结果，第2行表示对“T2”的去噪结果。通过图1可以看到，所提出方法与DBF方法所获得的去噪结果在视觉上几乎一样。在所提出方法中bins=32，bins的值代表对图像量化时将多少个像素值量化到一个值。

　　显然，量化参数s能够加快运算速度。对s采样越稀疏越能获得更快的速度。同样，也分析了由s的量化程度对去噪结果可能带来的影响。使用峰值信噪比（PSNR）来比较所提出方法的去噪结果与原始DBF方法的去噪结果。对于两幅图像I和I*，PSNR定义为10log，其中N为像素个数，图像灰度范围归一化到[0，1]。图2（a）为图1中测试结果的PSNR值。可以看到，所提出方法即使是在灰度被量化到bins很大时得到的结果与DBF方法也很近似。

　　图2（b）显示了在不同bins数目条件下的计算时间。从图中可以看到所提出方法的运行时间不受滤波半径影响，然而DBF算法的运行时间是随着滤波半径的增大而增大的。在合理的滤波半径下，可以明显看到所提出方法的加速性能。同时考虑到构建CIH的时间消耗，消耗的内存是随着图像分辨率以及灰度量化的程度而线性增加的。从图2（b）可以看到，适当减小bins的数量能得到理想的滤波效果同时内存的消耗也比较小。

4 结论

　　本文提出了一种基于交叉积分直方图的多模图像去噪方法。尽管这篇文章提出的方法与DBF方法很接近，但是实验结果表明该算法在计算效率上明显优于DBF方法，特别是针对需要较大滤波半径进行去噪的高分辨率噪声图像，该方法能更快速计算并取得与DBF相近的效果。该算法适合于各种图像去噪以及其他计算机视觉方面的应用。

参考文献

　　[1] 朱虹.数字图像处理基础[M].北京：科学出版社，2005.

　　[2] 章毓晋.图像处理和分析（第2版）[M].北京：清华大学出版社，2002.

　　[3] 成敏，马士友，董睿，等.基于非局部均值化的医学图像去噪[J].微型机与应用，2013，32（7）：34-35，39.

　　[4] SOTIRAS A， DAVATZIKOS C， PARAGIOS N. Deformable medical image registration： a survey[J]. IEEE Transactions on Medical Imaging， 2013， 32（7）：1153-1190.

　　[5] TOMASI C， MANDUCHI R. Bilateral filtering for gray and color images[C]. The 6th International Conference on Computer Vision， 1998： 839-846.

　　[6] ELAD M. On the origin of the bilateral flter and ways to improve it[J]. IEEE Transactions on Image Processing， 2002， 11（10）：1141-1151.

　　[7] PETSCHNIGG G， SZELISKI R， AGRAWALA M， et al.Digital photography with flash and no-flash image pairs[J].ACM Transactions on Graphics， 2004， 23（3）：661-669.

　　[8] Zhang Ke， GAUTHIER L， RUDY L， et al. Constant time joint bilateral filtering using joint integral histograms[J]. IEEE Transactions on Image Processing A Publication of the IEEE Signal Processing Society，2012，21（9）：4309-4314.

　　[9] BENNETT E P， MASON J L， MCMILLAN L. Multispectral bilateral video fusion[J]. IEEE Transactions on Image Processing， 2007， 16（5）：1185-1194.

　　[10] PORIKLI F. Integral histogram： A fast way to extract histograms in cartesian spaces[C]. IEEE Conference on Computer Vision and Pattern Recognition（CVPR 2005）， 2005： 829-836.