基于单应性矩阵追踪的自动焦距标定算法

孙翌乔，王春勇，李振华，张宝辉

(1.南京理工大学 理学院，江苏 南京 210094；2.昆明物理研究所，云南 昆明 650223)

在旋转式全景成像系统中，由于视角在不断地拼接中发生了扩展，多张图像拼接时会出现比较严重的图像畸变现象，最突出的表现是拼接图像呈橄榄球状，部分区域甚至难以识别。为此，本文采用了柱面投影算法，将处于不同平面的图像投影到统一的曲面上，从而减少图像畸变，实现多张图像连续拼接。

柱面投影算法对焦距(或视角)这一参数要求相对较高，利用不精确焦距实现的柱面投影不仅无法解决前述的畸变问题，还会因为畸变大幅降低各种匹配算法获得的特征点数量和质量，甚至导致匹配失败。文献[1]和[2]直接采用了已知内参数的摄像机，从而绕过了焦距标定的过程，但是并非所有摄像机都可以提供足够精确的内参数。文献[3]通过摄像机旋转一周，根据图像张数估算视角来实现焦距的标定，但这一方法偏差明显较大。文献[4]在此基础上加入了图像重叠率这一参数，一定程度上提高了精度，但仍未达到理想的精度。文献[5]和[6]分别采用了光束平差法和基于纯旋转运动的列文伯格-马夸尔特(Levenberg-Marquardt，LM)迭代方法，它们的共同点是精度较高，但复杂度较高、耗时长，难以实现实时拼接。文献[7]对文献[6]进行了简化，提出了通过追踪单应性矩阵参数标定焦距的方法，但使用的折半查找法未能充分利用参数变化趋势，精度略有不足且不适合并行运算，有进一步优化的空间。

本文主要针对文献[7]存在的问题，使用二次函数求解法代替折半查找法，以进一步提高标定精度。

1 旋转式全景成像系统工作流程

如图1所示，本文所述的旋转式全景成像系统工作流程如下:首先对2张待拼接图像提取高速鲁棒特征点(Speed up robust features，SURF)并获取描述子，然后对2张图像中的特征点进行匹配获得匹配点集。如果焦距已知，可以但不推荐直接采用已知的焦距进行柱面投影，利用投影后的点集经随机抽样一致(Random sample consensus，RANSAC)算法提纯后获取单应性矩阵完成拼接。如果焦距未知、不准确或存在动态变焦，则使用本文所述的方法进行焦距标定后进行后续步骤。

图1 全景成像系统工作流程

2 柱面投影算法

2.1 像素焦距

在介绍柱面投影之前，有必要介绍一下像素焦距的概念，在下文中，除特别标明外，所有的“焦距”文字和f符号均代表像素焦距。

在实际生活中，成像设备往往会提供以mm为单位的物理焦距，如荣耀30S手机拍摄出的相片显示焦距6 mm；并给出等效焦距，如手机提供26 mm/35 mm等效焦距，这代表此设备与焦距为26 mm的镜头在35 mm标准胶片上的视角等效。

但在实际图像处理中需要把这一物理焦距转换到像素域中，形成以像素为单位的像素焦距以简化计算。转换的基础是物理焦距和像素焦距所对应的视角相等。

如图2所示，设W为图像宽度(单位为px)，θ为半视角，l为感光器件宽度，f像素和f实际分别是像素焦距和实际焦距，根据对顶角定理

整理式(1)可以得到焦距转换式(2)

图2 焦距转换示意图

2.2 柱面投影

如图3所示，w和h所构成的平面即为原始的像平面，以光心O为中心的柱体即为目标曲面，其半径为焦距f，在平面上任取一点P(x，y)，连接OP交柱体于P′(x′，y′)，θ是水平方向的半视角。沿h/2垂直于h，w平面切割柱体得到图3(b)，沿OP竖直切割柱体得到切面c。

首先确定投影后的图像宽度neww。由图3(b)可知，投影后图像宽度即为视角所对应的弧长

实践中neww往往不是整数，此处采用向下取整的方法。

图3 柱面投影示意图

接下来根据P(x，y)的坐标确定对应点P(x′，y′)的坐标。由于x′是一段圆弧，求解x′的问题可以转换成求解∠AOP的问题

因此

根据∠OPM～∠OP′M′可得式(6)

式(5)和(6)构成了柱面投影的正向投影公式。将式(5)和式(6)取反，即可获得逆向投影公式，如式(7)

柱面投影将原本矩形图像转换成中间高两边低且宽度有所收窄的新图像，随着焦距变大，收窄程度变大，最小高度变小。以荣耀30S手机拍摄的960×540照片为例，图4是同一照片经不同焦距的柱面投影所得图像。

图4 柱面投影效果图

3 基于单应性矩阵参数追踪的自标定算法

本文所述的焦距自标定算法分为粗标定和精标定2个步骤。粗标定的目的是获得一个接近实际焦距的焦距点，保证实际焦距落在其附近，以提高后续拟合的精度。精标定则是在此基础上将焦距修正到目标值。算法流程如图5所示。

图5 焦距标定算法流程

3.1 焦距粗标定

在进行焦距粗标定算法前，不妨关注一下假如不进行柱面投影，2张图像直接进行拼接的拼接效果。图6是某两张540×960图像直接进行拼接的结果。

在此情形下，拼接图像经透视变换被投影到被拼接图像同一平面上，形成了类似梯形的形状。现将这一场景表示在图7上。

图6 直接进行拼接的结果

图7 直接拼接示意图

由图7(a)可知，设拼接后最大宽度为newW，最大高度为newH，又OO′=f

在图7(b)中，OC=f，，则根据三角形相似

联立式(8)和式(9)，得到f的表达式

此算法对焦距的标定速度较快，只与图像本身相关，与摄像机旋转速度无关，也无需提供初始值，但精度相对较低，存在较大方差，这是由于忽略了可能存在的竖直方向位移、图像中心未必与摄像机中心位于同一水平面等原因导致的。但这个精度已经足以满足下一小节的焦距精标定要求，且二者共用一套匹配点集，不会产生更多的时间消耗。

取8张精标定焦距为824 px的图像两两进行本节的焦距粗标定算法，得到表1。

表1 焦距粗标定结果

3.2 焦距精标定算法

获得预测值后就要进行焦距的精标定。本文通过SURF算法获得特征点配对组，经过RANSAC算法筛选后计算获得单应性矩阵H作为透视变换的变换矩阵。文献[7]指出，h1、h2、h4、h5主要与旋转变换有关，h3、h6分别是横向和纵向位移，h7、h8则是与透视变换相关的主要参数。最后一项恒为1。

现选取两张相邻的图像P1、P2进行SURF特征点匹配，首先经RANSAC算法提纯后获得匹配点集Ω，在f粗附近选取100个f值记为f k，然后根据式(5)、式(6)得到Ω′k，其中Ω′k为Ω经以f k为焦距的变换得到的点集。分别计算单应性矩阵H′k，记录与透视变换相关的参数h7、h8的变化情况，结果如图8所示，其中近似水平的曲线代表h8的变化情况，近似二次函数的曲线代表h7的变化情况。

图8 h7、h8的变化情况

当h7=h8=0时，透视变换退化成仿射变换。由于h8始终在0附近且变化幅度较小，因此不能作为追踪参数。h7随f增加而增加且具有明显的平滑曲线特征，可以作为追踪参数。

现在f粗附近任取3个点，利用二次函数的三点式可得到唯一确定的二次函数，这个函数与h7实际变化的关系如图9

图9展示了在目标点附近实际变化曲线与拟合二次函数的关系，其中蓝色部分为h7的实际变化曲线，橙色曲线为拟合的二次函数，绿色垂线表示实际焦距位置。由图可见h7在实际焦距f±20%区域内的变化几乎完全符合二次函数的规律，计算可知其互相关系数高达0.999 9，因此可以通过求解此二次函数的零点来获得精确的焦距值。整个过程只需要3次求解单应性矩阵(且可以通过并行运算进一步压缩时间)即可获得足够精确的焦距。尽管在实际焦距附近这个变化几乎完全符合二次函数变化规律，但由图9其他部分可以明显看出，远离实际焦距时两条曲线将逐渐偏离，此时精标定算法将有失效的风险，甚至会输出明显错误的复数解。因此，在焦距精标定前先进行粗标定是必要的。

图9 二次函数与实际变化趋势关系

现通过以下方法估测本文算法精度:本文算法精度与粗标定精度有关，粗标定偏差越大则二次函数曲线与参数变化曲线在零点处偏差越大。以粗标定中偏差最大的结果第8组的粗标定值880 px为例，选取附近1 000个f值作为测试点，得到样本变化曲线，本文与文献[7]的精标定过程，本质上都是寻找这个曲线的零点。现选取这些点的首、中、末3个f值计算对应h7构建二次函数，并形成图10，寻找2条曲线零点所在区间，则最大误差σ＜829 px-828.1 px=0.9 px。相同条件下原文献算法标定值为830 px，误差为1.7 px。

此外，由于可以使用上一次精标定的值作为粗标定输入，此算法对h7变化曲线零点查找的精度随着运算次数的增加会迅速收敛至小于0.01 px，此后标定精度将不再受这一过程影响，而只与特征点提取和匹配过程中产生的误差有关。

图10 本文算法误差分析

由于不准确的焦距会导致图像高度发生变化，本文采用文献[7]提出的评价方法，此算法基于当焦距足够精确时，相邻图像之间只存在平移变换，因此同样大小的图像在拼接后大小不发生变化的原理，以最左端和最右端图像高度差作为评价方法，但文献[7]采用了单向连续拼接的方法，这一方法误差最大的部分位于图像一端，会造成误差的累积。本文采用的拼接方法通过交替使用左向拼接和右向拼接，平衡了误差的累计，误差最大的位置也不再位于两端，因此以拼接后图像左右分别去除半幅图像宽度后(排除柱面投影影响)最大高度与最小高度的差Dev=Hmax-Hmin作为拼接效果的评价方法，当这个差值大于由柱面投影导致的单幅图像内图像最大高度与最小高度的差dev=hmax-hmin时或拼接过程中因错误中断时判定拼接失败。不同算法的Dev值表现如表2。

表2 用差值表示的拼接效果

上述结论证明本文方法有效地实现了焦距的校准，并且进一步提升了焦距校准的精度，并能适应自动变焦的环境。图11展示了27张图像的拼接成果，拼接效果良好，畸变较小。

图11 本文连续图像拼接结果(27张)

3 结束语

本文首先通过直接进行图像拼接时的几何关系实现了焦距的快速粗标定，然后在此基础上对单应性矩阵参数的追踪拟合实现了旋转全景成像系统中的焦距精标定，为后续的图像拼接打下了坚实的基础。本文在文献[7]基础上进一步提高了拼接效果，减少了图像畸变。但是，这一算法依旧是建立在摄像机的旋转运动上的，具有一定的局限性，主要针对的仍是旋转式全景成像设备，对诸如多摄像头阵列等全景成像系统没有作用。接下来摄像机位置不处于旋转轴上的图像拼接将成为进一步研究的方向。