基于超分辨率图像重建技术的运动估计研究
2015-05-29顾刚袁杰
顾刚+袁杰
摘 要： 阐述用于超分辨率重建对两帧连续图像运动估计的若干算法，对其中的块匹配算法以及相位相关法进行研究。对块匹配算法主要研究三步搜索策略以及快速搜索策略，实验说明块匹配算法当平移变换量在一定范围内，图像灰度变化连续稳定，但在处理压缩图像序列时易陷入局部最小值。而相位相关法使用的是傅里叶频域方法对有相对运动的两幅图像进行运动估计，实验表明相位相关法不仅可以用于平移变化，且对于旋转变化也可以通过坐标变换后再通过此法进行运动估计。
关键字： 运动估计； 超分辨率重建； 块匹配； 搜索策略
中图分类号： TN919?34； TP751.1 文献标识码： A 文章编号： 1004?373X（2015）10?0080?07
在数字图像的采集与处理过程中，有许多因素会导致图像分辨率的下降，其主要表现为模糊、噪声和变形。为了使获得的图像具有更高的分辨率、更好的质量，通常采用基于信号处理的软件方法对图像的空间分辨率进行提高，即超分辨率图像重建[1]。
超分辨率图像重建就是从一系列质量较差、分辨率较低的图像来重建图像质量更好、空间分辨率更高的图像的算法[2]。其核心思想就是用时间带宽（获取同一场景的多帧图像序列）换取空间分辨率，实现时间分辨率向空间分辨率的转换[3]。由于在成象的场景中一般有多个物体做不同的运动，如果直接按照不同类型的运动将图像分割成复杂的区域是比较困难的。因此通常采用基于块的物体运动表示法。
1 基于块的运动估计
1.1 块运动模型
假设图像运动可以用块运动来表征，块运动通常分为平移、旋转、仿射、透视等运动形式，一般情况下，块运动是这些运动的组合，称为变形运动。
1.2 块的运动估计算法
目前基于块的运动估计算法有如下几种：块匹配方法、相位相关法和傅里叶方法[4?5]，下面对其分别进行分析。
1.2.1 块匹配算法
块匹配算法的基本思想如图1所示。在第[k]帧中选择以[（x，y）]为中心、大小为[m×n]的块[W，]然后在第[k+1]帧中的一个较大的搜索窗口内寻找与块[W]尺寸相同的最佳匹配块的中心位移矢量[r=（Δx，Δy）]。搜索窗口一般是以第[k]帧中的块[W]为中心的一个对称窗口，其大小常常根据先验知识或经验来确定[6?8]。
图1 块匹配示意图
各种块匹配算法的差异主要体现在如下几个方面：匹配准则、搜索策略、块尺寸选择方法。
（1） 匹配准则
目前典型的匹配准则有：最大互相关准则、最小均方差准则、最小平均绝对值差、最小绝对误差和准则等[9?10]。由于最小绝对误差和准则（SAD）具有不需乘法运算、实现简单方便的优点，在本设计中选用此准则：
[SAD（Δx，Δy）=（x，y）∈WI（x，y，k）-I（x+Δx，y+Δy，k+1）] （1）
（2） 搜索策略
目前的搜索策略有全搜索策略、快速搜索策略、亚像素搜索，下面分别进行介绍：
全搜索策略：计算所有可能位移矢量对应的匹配误差，然后选择最小匹配误差对应的矢量就是最佳位移估计值。这种策略的最大优点是可以找到全局最优值，但十分浪费时间，因此，提出了各种快速搜索策略。
快速搜索策略：本设计中采用的是n步搜索或对数搜索快速搜索策略。其具体设计思路是：设窗口大小为15×15，当前像素值位于窗口中心，用“0”来标记，如图2（a）所示。第一步，选择标记为“0”和“1”的9个像素计算匹配准则函数，如果最佳匹配仍在“0”处，则无运动。第二步，以第一步最佳匹配对应的像素点为中心选择8个点（图中用标记“2”表示），计算这8个点的匹配准则函数值。第三步，以第二步最佳匹配对应的像素点为中心选择8个点（图中用标记“8”表示），计算这8个点的匹配准则函数值，最佳匹配值即为最后的最佳运动估计。由图2（a）可见，每进行一步，搜索距离减小一半，并且越来越接近精确解。
人们将上述搜索过程称为三步搜索[8?9]。进行超分辨率图像配准是可以继续在子像素级上进行搜索，以得到更精确的估计值，这样就需要大于三步的搜索，称之为[n]步搜索，由于搜索步数与窗口内像素个数是对数关系，因此，常将这种搜索称为对数搜索。另一种对数搜索策略是在每一步有4个搜索位置，它们以十字形或交叉形布置，如图2（b）所示。亚像素搜索：仅搜索到整像素精度的运动矢量远不能满足高质量的图像重建需要，经过三步搜索或全搜索得到整像素的偏移量后，再以最佳整点匹配宏块为中心进行半像素搜索。
图2 对数搜索示意图
算法步骤：
步骤1：设整像素精度的搜索结果为以点1为起点的宏块，设点1为（m，n），如图3所示。
图3 半像素搜索原理示意图
步骤2：用线性插值方法求得以点（m，n）周围8个半像素点为起点的宏块的各点灰度值（以点2，6为例）：
[G2=0.5×（I1（m，n-1）+I1（m，n））] （2）
[G6=0.25×（I1（m-1，n-1）+I1（m-1，n）+I1（m，n-1）+I1（m，n））] （3）
步骤3：仍然使用最小SAD准则，求得周围的8个半像素点以及点1为起点的宏块的匹配准则函数值。SAD值最小的就是半像素搜索结果。同理，可以再以搜索到的半像素点为中心，连同周围的8个[14]像素点进行搜索，使得运动估计的精度达到[14]像素。
1.2.2 相位相关法
（1） 仅有平移变换的两幅图像的配准
相位相关技术是一种非线性、基于傅氏功率谱的频域相关技术，经常被用来检测两幅图像之间的平移，根据频域信息，利用相关技术能够快速地找到最佳匹配位置[11?12]。
设[f1（x，y）]为参考图像，[f2（x，y）]为仅存在位移变化的实测图像，其位移值为[（x0，y0）]，即[f2（x，y）]=[f1（x-x0，y-y0）]，设[F2（ξ，η）]和[F1（ξ，η）]分别为[f2（x，y）]和[f1（x，y）]的傅氏变换，则有：
[F2（ξ，η）=e-j2π（ξx0+ηy0）*F1（ξ，η）] （4）
则功率谱定义为：
[F1*（ξ，η）×F2（ξ，η）F1*（ξ，η）×F2（ξ，η）=e-j2π（ξx0+ηy0）] （5）
式中[F1*]为[F1]的共轭[13?16]。
对式（4）再做傅氏反变换，则可检测到一δ函数，δ函数的位置即为其位移[（x0，y0）]。
（2） 存在平移和旋转的两幅图像的配准
假设[f1（x，y）]为参考图像，[f2（x，y）]为存在位移和旋转变化的实测图像，令：
[f2（x，y）=f1（xcosθ0+ysinθ0-x0，-xsinθ0+ycosθ0-y0）] （6）
根据傅里叶变换的平移和旋转理论，有：
[F2（ξ，η）=e-j2π（ξx0+ηy0）*F1（ξcosθ0+ηsinθ0，-ξsinθ0+ηcosθ0）] （7）
假设M1，M2为[F1]，[F2]的幅值，根据式（7），有：
[M2（ξ，η）=M1（ξcosθ0+ηsinθ0，-ξsinθ0+ηcosθ0）] （8）
从式（8）看出，[F1]，[F2]的幅值是相同的，只不过是相对旋转了[θ0]，式（8）用极坐标写成：
[M2（ρ，θ）=M1（ρ，θ-θ0）] （9）
由式（9）可见，当两幅图像同时又平移和旋转时，如果是用极坐标表示频谱的幅值，则可表示成只有旋转角度的平移，其幅度是不发生变化的，这样就可以用只发生平移变化的方法来检测旋转角度[17]。
（3） 坐标变换
由以上分析可知，要配准发生平移和旋转的两幅图像时，可先将两幅图像的频谱由迪卡尔坐标系转化为极坐标系，然后用相位相关法检测出两幅图像相差的角度，角度配准以后再配准平衡移量，可见从迪卡尔坐标系变到极坐标系也是用此方法进行图像配准比较关键的一步。因此，引入如下极坐标：
[ρ=ξ2+η2 ， θ=tan-1ηξ] （10）
则：
[ξ=ρcosθ，η=ρsinθ] （11）
在极坐标系中，横轴坐标为[θ]，对[θ]的采样间隔可取0.01，为了提高精度，[θ]的采样间隔可以更小些，但这样会降低计算速度。设采样间隔为[Δθ]，则变换后[ρ]?[θ]平面中图像的大小为[360Δθ][×][ρ]，设极坐标系中点[（ρ，θ）]的值为[G（ρ，θ）]，求出的[ξ]、[η]通常落在原频谱图像[F（ξ，η）]的4个像素之间，因此，要求[G（ρ，θ）]的值必须用插值的办法来计算。常用的插值方法有最邻近插值、双线性插值和高阶插值。由于双线性插值速度较快，而且插值效果也比较好，因此选用双线性插值来求[G（ρ，θ）]。在极坐标系中，[θ]的取值范围为0°～360°，横轴坐标[ρ]，纵轴坐标为[θ]，[θ]的采样间隔取0.01。
2 实验结果
2.1 块匹配算法
2.1.1 简单确定的降质模型
三步搜索法原理如图4所示。设图5中有以（i，j）为左上点的大小为L×L的块，要在图5（b）中搜索到与之匹配的块，则将图5（a）中的点（i，j）标记为点1，加上周围8个点，计算它们的SAD，SAD最小值对应的点，就是第一步搜索的结果。再以这一点为中心，加上周围8个点，此时步长减半，与第一步同样的方法进行第二步搜索，同理进行第三步搜索。到此为止，实现了精确到整像素的快速搜索。
图4 三步搜索法原理示意图
图5 “cir.bmp”原图和左移3像素后的图
以“cir.bmp”图为例，图5（a）为原图，图5（b）为原图向左平移3像素，根据上述原理，可以预想搜索结果为图中第7行第9列加粗的3号点。
用坐标图来表示每一步搜索中各点的SAD值，SAD值最小点便是这一步的搜索结果。由图6可以看到，第一步搜索中，3号点的SAD值最小，第二步搜索中，7号点的SAD值最小，第三步搜索中，3号点的SAD值最小，因此参考原理图4，可见的确搜索到了红色3号点。在块匹配运动估计的快速算法中，三步法等在希疏格点上的搜索算法，运算量低，但容易陷入局部极小点。因此对于图片有一定的要求，当两图之间变化不剧烈，差异较小，可以用快速匹配法来进行处理。而当两幅图像变化较为剧烈，或是像素的灰度值变化剧烈，容易因快速算法导致陷入局部最小值而得不到正确的结果。
图6 整像素搜索分步示意图
超分辨率的图像恢复是要求由一系列低分辨率图像中构造出一幅或多幅高分辨率的图像，接下来进行的是使精度达到半像素的搜索方法。有关半像素搜索的原理图见图3以及相关阐述。以大小为M×N的“cir.bmp”为原图，将其行和列分别缩小1倍，即每4个像素平均得到一个像素值，得到大小为[M2×N2]的图，如图7所示。
图7 图片变化示意图
将进行半像素搜索中的SAD值标出，如图8所示，可以看到，SAD最小值所在点为4号点。从代码运行中可以看出，整像素搜索部分的结果为x=0，y=-2。意思是（i，j）左移两个点为整像素搜索的结果，见图4，即为红色3号点。
半像素搜索的结果为4号点，4号点是1号点右边0.5个像素点，而1号点是左移2个像素点，综合可知，搜索结果为左移1.5像素，这正是制造出的图7中7.3和7.4的关系。
2.1.2 JPEG 2000压缩图像序列的验证
在实际应用中，常常遇到的是经过压缩的图像，而不是简单确定的模型，因此一个运动估计算法的好坏不能单凭简单模型的验证，而要通过更接近实际情况的压缩图像的验证来评判优劣。分别采用整像素搜索以及半像素搜索进行处理。
（1） 整像素搜索。模型如图9所示，两幅图都是JPEG 2000格式，相当于有噪声的降质图，图9（a）为原图，图9（b）为原图向右移3个像素，向上移1个像素所得。
图8 半像素搜索部分的SAD值
图9 整体像素搜索模型
三步搜索的 SAD值如图10所示：可以看到，三步的SAD最小值所在点分别是7号点，1号点，2号点，由图4可知，搜索结果为蓝色7号点左上的点，正是原（i，j）点向右3像素，向上1像素。可见，若是精确到整像素，快速搜索的块匹配方法可以对JPEG 2000压缩图像进行运动估计。
图10 JPEG 2000压缩模型的三步搜索SAD值
（2） 半像素搜索。 模型如图11所示：这两幅图是由图9（a）和图9（b）压缩而得到的，可以认为图11（b）是图11（a）向上移0.5像素，向右移1.5像素。
整像素搜索SAD值如图12所示。
可见，三步搜索的SAD最小值所在点分别为7号点，1号点，1号点，参考图4，[Δx1]=0，[Δy1]=4。在进行半像素搜索，SAD值如图13所示，可见SAD最小值所在点为6，即[Δx2=-0.5，Δy2=-0.5]。综上所述，有[Δx=Δx1+][Δx2=-0.5，Δy=Δy1+Δy2=3.5]。制造模型的[Δx=-0.5，Δy=1.5]，结果不符。注意到，整像素搜索的结果有较大偏差，而在图13中不难看出，点1虽然是最小值，SAD1=1 249，但点3的值与点1相近，如图4所示，SAD3=1 265。
图11 半像素搜索模型
图12 JPEG 2000再压缩图像三步搜索SAD值
的半像素搜索结果 ]
如果将点3作为整像素搜索第二步的结果，参考图3可知： [Δx1=0，Δy1=2]再进行半像素搜索，结果如图15所示。则：[Δx=Δx1+Δx2=Δy1+Δy2=][1.5-0.5，Δy=1.5]，与制造出的模型完全符合。
2.2 相位相关方法
2.2.1 仅存在平移变换的两幅图像的配准
（1） 简单确定的降质模型
图16（a）为lean原图，（b）为原因右移20个像素。对图16进行相位相关处理。具体算法见1.2.2节。
图17中，X=21，则[Δx]=20，和实际吻合。
半像素搜索结果 20个像素 ]
图17 仅存在平移变换的相位相关结果
（2） JPEG 2000压缩模型
对图9进行相位相关处理，处理后结果如图18所示，可以看到，x=4，y=104，则：
[Δx=4-1=3， Δy=104-1-104=-1]
即向右3像素，向上1像素，和模型的平移变换一致。
图18 JPEG 2000压缩图像平移变换相位相关
2.2.2 仅存在旋转的两幅图像的配准
在处理旋转变换时，最关键的就是进行坐标系的变换。把图像从迪卡尔坐标变换到极坐标系：
[ξ=ρcosθ， η=ρsinθ] （12）
将极坐标中的点映射到迪卡尔坐标系中，通常落在非整像素点上，再通过周围四个像素点根据权重进行插值，例如：[ξ]=14.3，[η]=15.2，见图19，则第14行占比重0.7，第15行占比重0.3，第15列占比重0.8，第16列占比重0.2，则点（14，15）的权重为0.7×0.8，点（14，16）的权重为0.7×0.2，点（15，15）的权重为0.3×0.8，点（15，16）的权重为0.3×0.2。其算法如图20所示。
对图21中存在旋转的两幅普通图像进行相位相关处理，结果如图22所示。
图19 坐标变换示意图
图20 存在平移和旋转变换时的相位相关算法
图21 lena原图和逆时针旋转30°图
图22 只有旋转变换的相位相关结果
从图22可以看到y=31，即[Δθ]=30°，结果完全正确。对图23中的JPEG 2000压缩原图（a）和其顺时针旋转[π20]的图（b）进行相位相关处理。处理结果如图24所示：可见，y=352，根据代码，[Δθ]的精度为1°，所以[Δθ]=352-1-360=-9°，即顺时针旋转9°，与模型中的顺时针旋转[π20]相一致。
实验验证结果说明了相位相关方法的特点：
（1） 尖锐的相关峰。当两幅图像确实相关时，由于检测结果为一个δ函数，有非常尖锐的相关峰，因此能实现图像的精确配准。
（2） 相位相关算法对噪声有较高的容忍度，检测结果与照度无关，同时可以大大减少几何失真对匹配性能的影响。因此对于简单确定的降质模型和高倍率压缩的模型都能达到精确的相关估计。
图23 压缩原图和顺时针旋转[π20]图
图24 JPEG 2000旋转变换的相位相关
3 结 语
运动估计是进行超分辨率图像重建时不可缺少的前提步骤，运动估计的精度对于重建的质量有着直接的影响，因此在进行运动估计时，力求算法稳定，精度高，对于噪声的容忍能力强。经过试验发现块匹配算法的应用有一定的局限性，如易陷入局部最小值、对于平移超过一定范围无法进行配准等；虽然基于平移的块模型的运动估计是简单的，但处理旋转、变形运动和不连续运动场时，效果不是很好。相位相关方法可以得到非常好的结果，对于平移变换，用文中的方法只能精确到整像素，对于旋转变换可以达到较高的精度，但是要在运算速度上做出牺牲。
注：本文通讯作者为袁杰。
参考文献
[1] 刘丁峰.超分辨率图像复原技术综述[J].软件导刊，2009（12）：183?185.
[2] WOOD S L， LAN Hsueh?Ban. Edge detection performance in super?resolution image reconstruction from camera arrays [J]. IEEE Transactions on Image Processing， 2006 （1）： 38?43.
[3] BAI Yun?fei， HU Jing， LUO Yu?pin. Self?adaptive blind super?resolution image reconstruction [J]. IEEE Transactions on Image Processing， 2010 （9）： 1208?1213.
[4] FENG Xiang?bin， CHEN Yong?hong. Digital image watermarking based on Super?Resolution Image Reconstruction [J]. IEEE Transactions on Image Processing， 2012 （10）： 1778?1782.
[5] 汤华莲，庄奕琪.运动估计的分层搜索算法及FPGA实现[J].现代电子技术，2004，27（4）： 84?86.
[6] YIN Bin， DUAN Hui?chuan. Image stabilization by combining gray?scale projection and block matching algorithm [J]. IEEE Transactions on Image Processing， 2009 （9）： 1262?1266
[7] LIU Yu?hong， TU Dan. A digital image stability Algorithm [J]. Computer Simulation， 2008， 48： 200?204.
[8] VELLA F， CASTORINA A， MANCUSO M， et al. Digital image stabilization by adaptive block motion vectors filtering [J]. IEEE Transactions on Consumer Electronics， 2002， 50： 796?801.
[9] ZAYNAB Ahmed， ABIR Jaafar Hussain， DHIYA Al?Jumeily. Edge detection for fast block?matching motion estimation to enhance mean predictive block matching algorithm [J]. IEEE Transactions on Consumer Electronics，2012， 52： 1?5.
[10] LIN Yih?Chum， TAI Shen?Chum. Fast full?search hock? atching algorithm for motion?compensated video compression [J]. IEEE Press，1996 （4）： 914?919.
[11] JIA Ruiming， ZHANG Hong， WANG Lei， et al. Digital image stabilization based on phase correlation [J]. IEEE Press， 2009 （3）： 485?489.
[12] YUAN Fei， ZHANG Hong， JIA Ruiming， digital image stabilization based on log?polar transform [C]// Proceedings of Fourth International Conference on Image and Graphics. [S.l.]： [s.n.]， 2007： 111?115.
[13] LANG Li?ying ZHANG Xiao?fang. Study on phase correlation algorithm in intelligent transportation systems [J]. IEEE Press， 2008 （2）： 275?279.
[14] LIU Ji?Lin， SONG Jia?Tao. Vehicle license plate recognition system with high performance [J]. Acta Automatica Sinica， 2003， 29（3）： 457?465.
[15] WEN C?Y， YU C?C， HUN Z?D. A 3D transformation to improve the legibility of license plate numbers [J]. Journalof Forensic Sciences， 2002， 47（3）： 578?585.
[16] ZITOVA B， F LUSSER J. Image registration methods： a survey [J]. Image and Vision Computing， 2003， 21（11）： 977?1000.
[17] ZHANG Jing， WANG Chang?shun， LIAO Wu?ling. An image mosaics algorithm based on improved phase correlation [J]. IEEE Press， 2009 （2）： 383?386.endprint
