朱霞,陈仁文,夏桦康,章飘艳
(南京航空航天大学 机械结构力学及控制国家重点实验室,江苏南京210016)
图像配准是指同一目标的两幅图像在空间位置的对准。图像配准技术过程称为图像匹配或者图像相关。用一幅浮动图像寻找一种空间几何变换,使它与另一幅参考图像在某些特征点上达到一致。它们的差异可能来源于不同的拍摄设备、不同拍摄位置、不同采集时间、不同传感器分辨率等,而这些将加大配准的难度[1-4]。通常情况下,图像配准由以下4 个基本的步骤组成:1)选择相似性测度算法;2)确定变换模型;3)确定图像插值方法;4)优化策略。本文将重点放在相似性测度算法和优化策略上面。相似性测度用来衡量两幅或多幅图像匹配的效果,运用一种准则使两幅或多幅图像在该准则下匹配效果达到最佳。互信息作为一种相似性测度,已被公认为是精度和鲁棒性最高的图像配准方法之一。互信息无需进行预处理,可以应用于单模或多模的配准,特别是当其中一幅图像有数据缺损时,也能得到很好的配准效果。基于互信息有好多种方法,如最大化互信息法、归一化互信息法等。由于互信息只统计两图像间的灰度相关信息,在噪声、分辨率等因素的影响下会导致两图像间的相关信息减少,很容易陷入局部极值,最终可能导致误配准。解决方法就是将互信息与其他判断条件相结合,如:梯度、形态学梯度等,从多方面改进互信息测度,提高配准的精确度。
图像配准中的优化策略可以归结为最优化问题,即如何在精度和速度之间取得最佳。目前求解优化问题的主要手段就是对优化问题的可行解空间进行搜索。搜索方法主要有:共轭梯度法、爬山法、遗传算法、模拟退火算法、粒子群优化(PSO)算法等。PSO 算法具有原理简单、易于实现、需要调整的参数较少的特点。但是传统的PSO 算法后期容易陷入早熟收敛状态,导致计算时间延长,计算量增大,运算效率较低,抗噪能力差,最终造成配准效果不好。魏本征等在原始PSO 算法的基础上提出改进型PSO(IPSO)算法,改进后的惯性权重因子由ω转变为(α+β×rand())ω,即对粒子群体而言,群体中惯性权重因子可以是域[αω,(α+β)ω]中的任意随机数,群体中既可存在较大的权重因子加强PSO算法的全局搜索能力,又可有较小的权重因子加强局部搜索能力,改进操作起到提高权衡全局搜索和局部搜索能力作用,提高了算法搜索能力和收敛速度,在应用中取得较好的效果[3]。魏本征等将其引入到医学图像的配准中,效果有了一定的改善,但是从速度、收敛效果、精确度来讲还可以提高。本文提出了将文化粒子群优化(CPSO)算法与最大化互信息相结合的改进方案,以解决图像配准中计算量过大、搜索速度慢等问题。
互信息(MI)的概念来源于信息论,是信息理论中的一个基本概念,互信息衡量的是两个或多个事件之间关系的紧密程度[2,5]。就一幅图像而言,可以看成是灰度的随机变量集,它含有信息量的多少可以用熵来表示。
式中:p(xi)为灰度值xi出现的概率。图像的互信息表示两幅图像所包含的共同信息,可用两幅图像的互信息作为图像配准的相似度函数,当互信息最大时两幅图像就实现了配准。两幅图的平均互信息可以由两幅图的熵和联合熵组成。
联合熵H(X,Y)的含义是图像X 和图像Y 中的事件同时发生的平均不确定性。定义为
对于图像X、图像Y 的互信息可表示为
这样,由(4)式可以看出,对于互信息的计算就转化为对两幅图像概率分布的计算。最大化互信息将信息论引入图像配准,假设有两幅图像,X(浮动图像)和Y(参考图像),当图像X 在某些特征点上与图像Y 达到一致时,互信息值最大,此时两幅图像的配准程度达到最佳。该方法不依赖人工干预,不需要对图像进行预处理,也不需要任何的假设性先验知识等,是实现配准自动化的一个非常好的选择[5]。因此,基于最大化互信息配准的目的就是寻找互信息最大时的配准参数:
式中:Fit(r1,r2,r3,…)为优化约束函数,配准参数r1,r2,r3,…为图像的缩放因子、旋转因子、平移因子等。
由于噪声的存在,互信息相似度函数会存在局部极大值或收敛缓慢等问题,本文针对这一问题探讨利用CPSO 算法进行图像配准的参数优化[6-7]。
PSO 算法最早是美国Kennedy 博士受鸟群觅食行为的启发,提出的一种生物进化算法。PSO 算法基于群体与适应度,它将每个个体视为无体积的粒子,在搜索空间以变速飞行,其飞行速度根据该粒子本身的历史经验以及同伴的历史经验进行动态调整。Reynolds 于1994年提出了文化算法,从微观和宏观不同层面模拟生物层面的进化和文化层面的进化,该算法的双层演化结构框架有利于算法的全局收敛性[8]。因此本文将文化算法引入到PSO 算法,提出了一种CPSO 算法,并将该算法用于求解图像配准的参数优化问题。
本文提出的算法包含两个进化空间:一个是由具体个体组成的群体空间,另一个是由在进化过程中获取的经验和知识组成的信念空间。其演化框架如图1所示,其中群体空间是算法进行问题求解的主空间[9-10]。通过演化操作和性能评价进行自身的迭代求解,形成个体经验,并通过accept()函数将个体经验传递到信念空间。信念空间将接收到的个体经验,不断地进行自身的性能评价与演化,形成群体经验。信念空间通过influence()函数对群体空间的个体行为规则进行修改,以使群体空间中得到更高的进化效率[11]。
图1 文化粒子群算法的基本框架Fig.1 The basic framework of cultural particle swarm optimization algorithm
由于PSO 算法是一种基于个体改进、种群协作与竞争的进化计算方法,具有理论简单、易于编码实现的特点,因此群体空间采用基本PSO 算法进行演化。粒子群可用下面5 元素形式描述:
式中:n 为群体规模;iter 为迭代次数;V 和X 分别表示所有粒子的速度空间和位置空间;Ffit为适应度,从位置空间映射到实数空间。在每一次迭代中,粒子通过跟踪个体极值pbest 与全局极值gbest 来更新自已的位置和速度,公式如下:
式中:v(t)是粒子的当前速度;x(t)是粒子的当前位置;rand()是[0,1]之间的随机函数;c1和c2是学习因子,通常都取值为2;ω 是加权系数,用来控制历史速度对当前速度的影响程度,一般在[0.1,0.9]之间取值。
本文针对权重问题提出一种改进的动态惯性权重调整策略,惯性权重调整公式如下:
信念空间定义为<S,N[n]>,其中:S 表示环境知识中存储的一组由历代种群所产生的最优个体集合;标准知识N 保存目标函数n 个主变量参数的变化区间。每个区间描述为<Ij,Lj,Uj>,j =1,2,…,n. 其中:Iq=[l,u]={l≤x≤u,x∈R},表示变量x 定义域边界的值。上下界u 和l 由给定的值域初始化,Lj表示参数j 在区间下限lj对应目标函数的适应值,Uj是参数j 在区间上限uj对应目标函数的适应值。
种群空间通过接受函数向信念空间提供一组最优子集,在最优化问题中一般是按一定的百分比取值,通常取20% ~25% ,或者根据问题环境的变化按一定的规则选取参数。
接受操作:在群体空间的粒子群演化过程中,每运行Acc 代时,用群体空间当前的全局最好值来替换信念空间中适应值差的粒子。
式中:Bnum和Dnum为文化算法中的两个参数,用于调节接受操作和影响操作的完成次数;kmax为粒子群预先设定的粒子群最大演化代数;k 为粒子群演化的当前代数。
标准知识N 的更新规则如下:
这里,第i 个个体影响参数j 的区间下边界。lkj表示参数j 第k 次迭代时的下边界;Lkj表示lkj对应目标函数的适应值。同理,第i 个个体也影响参数j 的区间上边界。ukj表示参数j 第k 次迭代时的上边界,表示ukj对应目标函数的适应值。
信念空间在形成群体经验后,通过影响函数对群体空间中个体的行为规则进行修改,以使个体空间得到更高的进化效率。影响函数决定哪种知识将被用于指导需要解决的问题。根据具体问题的不同,影响函数的设计也不同。不同的知识源对种群空间也有相应不同的影响函数。
影响操作为在群体空间的粒子群演化过程中,每运行Inf 代时,用信念空间当前的全局最好值来替换群体空间中适应值差的粒子。
结合互信息的计算方法,将CPSO 算法用于二幅图像的配准。算法流程为:
步骤1 确定算法的搜索空间,即浮动图像在水平、垂直方向的平移量,浮动图像绕原点旋转的角度。
步骤2 初始化设置。设置群体空间和信念空间的种群大小,群体规模n =40,惯性权因子ω =0.9,学习因子c1=c2=2,微粒子最小位置值xmin=0,最大位置值xmax=4.0,微粒子最小速度值vmin=-2.0,最大速度值vmax=2.0. 计算粒子的互信息值以及每一个粒子的最优值和全局最优值。
步骤3 根据群体空间粒子位置,计算粒子的实际输入和适应值函数(两幅图像的相似性测度——最大互信息)。
步骤4 更新信念空间粒子位置,计算粒子的实际输入和适应值函数。
步骤5 判断若k%Acc = =0 时,进行接受操作。
步骤6 更新群体空间粒子的位置和速度,根据最大化互信息判定检验粒子的速度和位置是否越界,若是,调整速度和位置。
步骤7 计算群体空间中各粒子的实际输入和适应值函数。
步骤8 更新信念空间粒子的位置和速度,根据最大化互信息判定检验粒子的速度和位置是否越界,若是,调整速度和位置。
步骤9 计算信念空间中粒子的实际输入和适应值函数。
步骤10 判断若k%Inf = =0 时,进行影响操作。
步骤11 当误差达到最初设定值或最大迭代次数时,输出gbest 对应的调度方案和目标值;否则,令k=k+1,返回步骤3.
本文在内存为2 GB、处理器为Intel Core(TM)2 DUO CPU T5250 @1.50 GHz 的PC 上使用MATLAB 7.0 编译软件进行图像配准实验。以图像配准为例,将PSO,IPSO,CPSO 算法独立运行100 次。各性能参数比较见表1. 图2(a)为浮动图像,图2(b)为参考图像,两幅图在水平方向、垂直方向和旋转角度上偏差。设计算法时将旋转角度范围初定为[-20°,+20°]. 图2(c)为使用PSO 配置后得到的图像,图2(d)为使用IPSO 配置后得到的图像,图2(e)为使用CPSO 配置后得到的图像。
表1 3 种算法中的参数比较Tab.1 The parameters of three algorithms
图2 基于不同粒子群优化算法的配准效果图Fig.2 The registration results based on different PSO algorithms
表1中给出3 种算法对图2仿真得到指标比较。从表1可以看出,本文提出的方法在标准差、熵和互信息这些指标上优于其他方法,只是灰度均值比其他方法略差一点。图3、图4和图5给出了实例中3 种算法执行后得到的误差收敛曲线。很容易看出,使用PSO 算法搜索时,经过2.5 s 之后就能使搜索到的位置达到最优;使用IPSO 算法,经过1.5 s之后到达最优;使用本文方法,经过1 s 就能使搜索达到标准。由此可知,CPSO 算法比PSO 算法和IPSO 算法提前收敛,运行时间更短。
图3 水平方向、垂直方向和角度误差(PSO 算法)Fig.3 The error of parameters in PSO algorithm
图4 水平方向、垂直方向和角度误差(IPSO 算法)Fig.4 The error of parameters in IPSO algorithm
图5 水平方向、垂直方向和角度误差(CPSO 算法)Fig.5 The error of parameters in CPSO algorithm
为进一步验证配准算法的可行性,利用归一化互相关通过转换一副图像来使两幅图像配准[12]。图6显示了两幅图像,其中左图为右图的子图。计算两幅图像之间的归一化互信息,如图7所示。归一化峰值最大的地方意味着这两幅图示最相关的地点,图像中总的变换量依赖于互相关中峰值的位置、图像的大小等。最后为了凸显配准效果,将全景图先处理为灰度图像,在全景图中显示彩色的子图,如图8所示,从中可以看出,两幅图像配准成功了。同样的算法,用于实际拍摄的图像图9(a),可以看出如果选取的子图像比较小,图9(b)的归一化相关图中的峰值就比较错乱。所以本文的算法比较适合子图在全景图中占百分比多一点的图像配准。
图6 待配准的图像(图像来源于MATLAB 素材库)Fig.6 Images under registration(original images from MATLAB)
图7 两幅图像的归一化相关Fig.7 Normalized correlation of two images
图8 最后配准的图像Fig.8 Registrated image
图像配准是将同一场景的两幅或多幅图像对齐的过程。图像的配准本质上是一个多维空间寻优的过程,针对具有平移、旋转和拍摄条件不同的图像,提出了一种基于文化粒子群与互信息的图像配准方法。大量实验结果表明,本算法对存在尺度以及旋转变化的图像能够实现准确的配准。用最大化互信息作为相似性测度可以解决图像间相关性的问题。CPSO 算法作为搜索方法可以解决传统PSO 算法搜索速度慢的问题。本文将此算法与传统PSO 算法和IPSO 算法的效果进行了比较,通过实验结果分析可知,本文提出的图像配准算法是一种比较稳健的配准算法,具有速度快、鲁棒性好等特点。
图9 配准后的图像(图像来源于实际拍摄)Fig.9 Registrated image (the images captured in the actual test environment)
References)
[1]丁南南,刘艳滢,张叶,等.基于SURF-DAISY 算法和随机kd 树的快速图像配准[J].光电子·激光,2012,23(7):1395 -1402.DING Nan-nan,LIU Yan-ying,ZHANG Ye,et al. Fast image registration based on SURF-DAISY algorithm and randomized kd trees [J]. Journal of Optoelectronics ·Laser,2012,23(7):1395 -1402.(in Chinese)
[2]张泽旭,李金宗,李冬冬. 基于光流场分析的红外图像自动配准方法研究[J].红外与毫米波学报,2003,22(4):307 -312.ZHANG Ze-xu,LI Jin-zong,LI Dong-dong .Research of automated image registration technique for infrared images based on optical flow field analysis [J]. Journal of Infrared and Milimeter Waves,2003,22(4):307 -312.(in Chinese)
[3]魏本征,赵志敏,宋一中. 基于IPSO 和综合信息的医学图像配准新方法[J].光电子·激光,2009,20(9):1271 -1274.WEI Ben-zheng,ZHAO Zhi-min,SONG Yi-zhong. A novel algorithm for multimodality medical image registration based on IPSO algorithm and hybrid information[J]. Journal of Optoelectronics·Laser,2009,20(9):1271 -1274.(in Chinese)
[4]柏连发,韩静,张毅,等.采用改进梯度互信息和粒子群优化算法的红外与可见光图像配准算法[J].红外与激光工程,2012,41(1):248 -254.BAI Lian-fa,HAN Jing,ZHANG Yi,et al. Registration algorithm of infrared and visible images based on improved gradient normalized mutual information and particle swarm optimization[J].Infrared and Laser Engineering,2012,41(1):248-254.(in Chinese)
[5]杨帆,张汗灵.蚁群算法和Powell 法结合的多分辨率三维图像配准[J].电子与信息学报,2007,29(3):622 -625.YANG Fan ,ZHANG Han-ling. Multiresolution 3D image registration using hybrid ant colony algorithm and Powell’s method[J].Journal of Electronics & Information Technology,2007,29(3):622 -625.(in Chinese)
[6]崔伟,刘圣霞,徐骞,等. 基于互信息和梯度的红外与可见光图像配准新方法[J].激光与红外,2011,41(2):224 -228.CUI Wei,LIU Sheng-xia,XU Qian. Novel infrared-visual image registration based on combined mutual and gradient information[J].Laser & Infrared,2011,41(2):224 -228.(in Chinese)
[7]尹丹艳,吴一全.灰色预测和混沌PSO 的红外小目标检测[J].智能系统学报,2011,6(2):126 -131.YIN Dan-yan,WU Yi-quan. The detection of a small infrared target based on gray prediction and chaotic PSO[J].CAAI Transaction on Intelligent System,2011,6(2):126 -131.(in Chinese)
[8]Reynolds R G. An introduction to cultural algorithms[C]∥Proceedings of the Third Annual Conference on Evolutionary Programming.River Edge,NJ:World Scientific Publishing,1994:131-139.
[9]朱霞. 一种求解作业车间调度的文化粒子群算法[J].计算机应用研究,2012,29(4):1234 -1240.ZHU Xia. Cultural particle swarm optimization algorithm for Jobshop scheduling problem[J].Application Research of Computers,2012,29(4):1234 -1240. (in Chinese)
[10]Reynolds R G,Peng B. Knowledge learning and social swarms in culture algorithms[J]. The Journal of Mathematic Sociology,2005,29(2):115 -132.
[11]Wang F,You H J,Fu X Y. Adapted anisotropic Gaussian SIFT matching strategy for sar registration[J]. IEEE Geoscience and Remote Sensing Letters,2015,12(1):160 -164.
[12]de Miranda P A V,Falcão A X,Udupa J K. Synergistic arc-weight estimation for interactive image segmentation using graphs[J].Computer Vision and Image Understanding,2010,114(1):85 -99.