可见光与热像融合的三维温度模型重建

张远辉， 徐栢锐， 朱俊江， 孙 坚

(中国计量大学 机电工程学院，浙江 杭州 310018)

1 引 言

热像仪通过检测场景中可见光频率之外的红外热辐射，可获取传统可见光相机无法获得的温度信息。近年来，红外监测技术飞速发展，由于其具有无接触、测温范围广、测温速度快等诸多优点，在各种工业场景中具有重要应用[1]。例如，热像仪获取的温度信息可用于快速确定模具是否需要冷却或者存在缺陷；在火灾现场快速定位火源，检测气体泄漏；在军事方面，热像仪作为侦察设备，用于了解敌方动态等[2,3]。但是，由于热像图只存在温度信息，缺少几何信息，因此，人眼可能难以通过温度分布情况分辨物体几何外观，进而对辨认物体形态以及定位造成困难。

为了解决这个问题，本文首先根据运动恢复结构(structure from motion, SFM)与多视图立体匹配(multiple view stereo, MVS)，完成相机位姿估计与深度图的重建，并通过基于互信息(mutual information, MI)的图像配准方法，完成了几何与温度信息的融合；最终，将二维平面的可见光图与热像图集合转化为三维热像温度模型，解决了在平面热像图上难以分辨物体表面几何形态的问题。

2 研究现状

目前，大部分针对热像图的三维重建的研究建立在可见光三维重建的基础上，Truong等人[4]与Maset等人[5]将用于可见光立体匹配中的SFM算法与MVS算法直接运用在热像图中获得深度图。然而，大部分可见光特征提取算法并不适用于热像图，并且由于传统立体视觉要求有一定数量的成对特征匹配，该方法要求目标场景中具有较明显的温度差异；否则，在热像图分辨率较低的情况下，热像图可提取用于对极匹配的特征点不足以实现精度较高的三维重建。

为了实现高精度的热像图三维重建，大多数研究通过可见光相机和红外热像仪之间的配准，利用可见光的三维重建结果，间接实现三维热像模型重建。

Yang等人[6]利用同时具有可见光与热像图输出的红外热像仪，以及额外两个可见光相机实现了热像图的三维重建，他们使用基于归一化互相关(normalized cross correlation, NCC)的配准方法，通过红外热像仪的可见光输出，间接实现了热像图与可见光图的配准，最后结合传统重建方法完成了温度信息的三维重建。不过，大部分红外热像仪不具有可见光的输出，导致该方法的适用范围较窄。

Yamaguchi等人[7]使用一台热像仪与一台可见光相机，采用直接视觉里程计(direct sparse odometry, DSO)进行相机位姿和深度图的估计，并利用预先处理好的校正信息，将热像图反向投影到可见光图像上，进而在可见光重建的基础上获得了热像的三维稀疏点云，该方法有效解决了直接利用热像图进行立体视觉重建效果不佳的问题，但由于DSO属于稀疏重建方法，因此最终三维热像点云模型较为稀疏，缺失了部分温度信息。为了解决这个问题，本文根据这个实验思路，结合基于几何和阴影的重建细化方法优化深度图的估计，实现了稠密的热像图三维重建。

此外，不少研究通过融合深度传感器和红外热像仪进行三维温度模型的重建。Matsumoto等人[8]结合Kinect深度相机和红外热像仪，在室内场景下，使用Kinect Fusion[9]的方法完成三维温度模型重建，实现了模型的在线重建。相比传统相机，深度传感器进行三维重建的速度更快，可以满足在线重建的要求，但由于RGB-D相机的深度检测对环境光照较为敏感，通常在光照干扰较多的场景下无法完成较好的重建。

3 总体流程

本文提出的方法使用固定在同一摄像架的一台可见光相机与一台热像仪进行三维热像模型的重建。首先通过将可见光相机与热像仪对同一标定物进行校正获取相机内外参数以及相对位移，随后，先利用增量式SFM估计各个视角的相机位姿，后使用基于阴影优化的MVS完成深度图的重建。由于单目视觉系统尺度因子未知，故无法利用现有深度图与校正信息进行配准。为了实现可见光图和热像图的配准，本文通过热像仪与可见光相机的深度图，利用互信息计算尺度因子完成坐标系归一化，实现了可见光相机与热像仪的配准；最后，结合配准结果和深度图，完成三维热像模型的重建。总体流程如图1所示。

图1 重建方法总体流程Fig.1 Reconstruction method overview

4 相机校正

本文基于张氏标定法[10]，对热像仪和可见光相机进行棋盘格标定计算可见光相机与热像仪的内参矩阵K和外参矩阵[R|t]。为了获取相对平移旋转变换，需要使用同一标定物同时对红外热像仪和可见光相机进行标定。然而，由于普通的光学标定板上的黑白格温度梯度过小，在红外图像上难以检测角点，不适合用于红外热像仪的标定。

为了解决这个问题，本文设计一种既适用于可见光相机，也可以用于热像仪的标定方法，我们使用具有一定隔热效果的白色PVC材料，设计了一块宽285 mm，高360 mm，厚1 mm，具有非对称圆孔的标定板，并将其置于黑色加热垫与相机之间。如图2所示，这种低成本的方法简单有效地在热像仪中产生了圆孔网格温差的同时，也在可见光相机中制造了明显的颜色差异，易实现圆心检测。

图2 圆孔标定板Fig.2 Grid-hole calibration board

此外，为了在下述尺度估计步骤中的简便计算，且因标定板各尺寸参数已知，校正过程中世界坐标系尺度定义为真实测量尺寸。

5 相机位姿估计

SFM是指从一系列图像集合中的匹配点恢复相机参数的一种算法[11]，主要可分为全局式SFM与增量式SFM。相较于全局式SFM，增量式SFM对错误的图像匹配有较高的鲁棒性，总体精度更高；但是，增量式SFM效率较低，且对大规模场景容易累积误差，造成位置漂移的现象。由于本文只对中、小场景中的物体进行三维重建，因此，选择增量式SFM作为复原相机参数的方法。

5.1 特征匹配

SFM算法首先需要检测可匹配的特征点，通常选择使用SIFT[12]、SURF[13]等完成特征匹配。由于SIFT、SURF具有旋转、尺度不变性，检测结果可以在邻近图像之间得到匹配。根据对极几何，空间中一个点与其在两个不同视图中的成像点具有共面约束[14]，故可依靠对极几何约束去除部分错误的匹配。由于特征点成对匹配的特性，匹配结果可以在邻近视图中不断扩展或者合并，以此建立特征点轨迹图，如图3所示，轨迹图用于在增量式SFM中选择下一个重建视图。

图3 轨迹图Fig.3 Track graph

5.2 初始位姿估计

根据增量式SFM，需要选择一对合适的视图用于初始位姿的估计，由初始位姿对特征点三角化的结果，不断从轨迹图中选择下一个视图进行重建。因此，一个合适的初始位姿对后续增量式重建的速度和准确度有较大的影响。本文通过8点法，根据本质矩阵估计初始位姿[15]，设其中一个视图位姿为P=[I3×3|O3×1]，根据本质矩阵的奇异值分解可估计另一个视图的相对位姿。设本质矩阵可分解为Udiag(1,1,0)V，则4个可能相对位姿如公式(1)所示，其中u3为U的最后一列。

P1=[UWVT|+u3]，P2=[UWVT|-u3]，P3=[UWTVT|+u3]，P4=[UWTVT|-u3]

(1)

其中

5.3 增量式重建

完成初始位姿估计后，在初始视图内可见的特征点均可完成三角化，即从二维像素坐标估计三维空间坐标，根据三角化结果，可选取在轨迹图中具有最多三角化特征的视图作为候选重建视图，进而完成增量式重建。

本文使用P3P算法[16]，根据三角化结果估计候选视图的姿态。如图4所示，为了优化位姿的估计，在位姿估计过程中进行光束法平差[17]最小化投影误差，即最小化目标函数(2)。

图4 光束法平差Fig.4 Bundle adjustment

(2)

其中

式中：Pi为相机光心；PiXj为三角化点Xj在视图i的投影二维位置；点Xj在视图i的原始像素位置为xij；d(x,y)2表示二维点x与y之间的欧式距离的平方。

6 深度图重建

为了实现稠密的三维重建，本文在SFM位姿估计的基础上进行了深度图的重建。由于传统几何立体匹配对弱纹理与无纹理区域深度估计的鲁棒性较差，所以，本文首先使用半全局立体匹配方法(semi-global matching, SGM)估计深度图，随后对初始深度图进行结合几何与阴影的深度优化。

6.1 深度图初始估计

SGM算法由Hirschmueller[18]首次提出，经过改进，Humenberger等人[19]提出基于Census变换代价计算的SGM算法。Census变换通过局部窗口计算视图相关性，对整体明暗不敏感，具有并行性高等优点，本文采用了基于Census变换的SGM算法进行深度图的初始化。

Census变换T(u,v)是将像素点p=(u,v)在m×n邻域的灰度强度值转换为比特串的一种变换，如式(3)所示。其中，I(u,v)为像素点p灰度强度值，⊕表示比特位的逐位连接计算。

(3)

其中

基于Census变换的像素代价值由计算像素点在两视图之间的汉明距离由式(4)得出：

C(p,d)=Hamming(Ti(u,v),Tj(Pij(p,d)))

(4)

式中：d为深度值；Pij(p,d)为像素从视图i到视图j的投影变换，由式(5)得出：

(5)

式中K与[R|t]分别表示相机内参和外参矩阵。

由于代价计算只考虑了局部相关性，为了获得较好的立体匹配效果，SGM算法采取全局能量优化策略，寻找每个像素的最佳深度使得全局能量函数(6)最小。

E(d)=E1(d)+E2(d)+E3(d)

(6)

其中

式中：P1和P2为用于平滑深度图的惩罚系数，分别对像素p在Np邻域内针对深度变化很小[|dp-dq|=1]以及深度变化较大[|dp-dq|>1]的情况进行惩罚。

为了高效地解决这个最优化问题，本文结合动态规划思想，从水平、垂直以及对角4个方向聚合得到路径代价值，选取代价最低值作为最优深度值，沿着某一路径r进行代价聚合，如式(7)所示。其中D为深度值取值范围，取决于SFM算法中的三角化结果。

(7)

其中

6.2 深度图优化

由于立体几何匹配算法对弱纹理或无纹理物体无法实现准确的深度估计，而基于阴影的重建方法对物体表面高频细节具有较好的重建效果[20]。因此，有不少学者提出了结合阴影与几何的重建算法，以几何结果作为输入，根据光照估计优化重建细节。本文采用Fabian等人[21]提出的基于Retinex理论[22]和图像梯度的优化算法。

根据Retinex理论，可以认为图像梯度小的区域主要由阴影明暗导致，而在图像梯度大的区域主要与物体表面反照率有关；故根据当前图像区域梯度，对梯度大的区域使用传统几何重建，而对梯度小的区域采用基于阴影的重建细化。因此，算法可转化为针对平衡几何数据项与阴影数据项的最优化问题。本文通过建立一个包含几何误差项与阴影明暗误差项的能量函数，使能量函数最小化从而实现深度图的优化。

6.2.1 几何误差

对于主视图中的像素点p与某个邻近视图i的几何误差，如公式(8)所示：

(8)

(9)

6.6.6 阴影误差

基于阴影的细化重建采用基于Lambertian反射的光照模型，利用三阶球谐函数Bh近似入射光强度，设a为反照率，l为光照系数，在像素点p的反射强度R(p)为：

(10)

(11)

(12)

(13)

6.6.6 能量函数

综上所述，结合几何误差与阴影误差，能量函数为：

E=EG+ES

(14)

其中

式中：NR为主视图与其邻近视图的集合；V为主视图与邻近视图均可视的像素点集合。由于Es在log域，故乘上系数α平衡尺度。

(15)

7 尺度估计

由于单目可见光相机估计的深度图尺度因子是未知的，无法根据校正信息使用式(5)进行投影转换，即由深度重建得到的深度值无法直接用于可见光图与热像图之间的配准，而需要一个额外的尺度因子α方可完成配准，如式(16)所示。

(16)

本文提出的方法通过互信息完成尺度因子α的估计。互信息用于描述两组随机变量之间相互依赖的程度，在医学图像配准中有广泛的应用[23]。由于RGB信息与温度信息非同一形式的变量，无法直接进行互信息比较，因此，本文选择使用深度值进行互信息的计算。

由于可见光相机与热像仪在每组视图中的相对位置都是固定的，因此可见光相机之间相对位姿等同于热像仪之间的相对位姿。故本文在可见光相机位姿估计的基础上，对热像图进行SGM算法的深度图估计，而因为热像图只存在温度信息，故不对热像图进行基于阴影明暗的深度优化。

为了估计最佳尺度，本文通过式(17)遍历尺度计算得分M，选取得分最高值对应的尺度作为当前视图的最佳尺度。

M(Ivis,Ithe)=n(H(Ivis)+H(Ithe)-H(Ivis,Ithe))

(17)

式中：Ithe表示当前热像图对应的深度图；Ivis表示当前尺度由热像图像素投影到可见光图中的区域；n表示此区域在整个可见光图中的占比。H(I)和H(I′)的定义如式(18)所示。

(18)

式中：pI(i)为在深度图I中深度值为i的概率；pII′(i,j)为在深度图I中深度值为i且在深度图I′中深度值为j的二维联合概率。

如图5所示，得分最高的尺度为视图当前最佳尺度。由于在校正过程中将世界坐标系定义为真实尺度，因此，本文通过测量目标物体的真实尺寸与对应重建结果中的尺寸得到尺度因子的真实值。

图5 尺度估计得分Fig.5 Scale estimation scores

然而，由于热像图深度估计的不准确，导致尺度估计结果存在不稳定性，当前得分最高的尺度并不一定为全局最佳尺度。为了选择一个全局最佳尺度，本文首先分别对每组视图进行尺度估计获取视图当前最佳尺度，根据得分权重，通过随机采样一致算法(random sample consensus, RANSAC)[24]，设定阈值并随机选择视图样本计算内点比例，以内点比例最高的尺度作为全局最佳尺度。

8 实验结果

本文实验装置如图6所示，可见光相机选型为Balser piA2400-17gc，分辨率为2 450×2 052，热像仪选型为Guide IPT640，分辨率为640×480。

图6 实验相机配置(左：可见光相机；右：热像仪)Fig.6 Camera setup(left: visible camera; right: thermal camera)

实验拍摄了3组可见光与热像图，如图7所示。

图7 实验图片输入Fig.7 Experiment figures

其中第一行图片由可见光相机拍摄，第二行图片由热像仪拍摄，3组输入在不同角度下的三维重建温度模型如图8所示。

图8 重建结果Fig.8 Reconstruction result

在尺度估计中，为了验证RANSAC算法的准确性，本文分别使用加权平均算法和RANSAC算法进行全局最佳尺度估计。如表1所示，相较于加权平均算法，使用RANSAC算法进行全局最佳尺度的估计的准确度更高。

表1 尺度估计算法对比Tab.1 Scale estimation algorithm comparison

9 结 论

综上所述，通过本文提出的三维温度模型重建方法，我们可以在三维模型上观测温度的分布情况，解决了二维平面热像图直观性不强的问题。此外，相较于Yang等人[6]的重建方法，本文使用的方法无需红外热像仪具有可见光输出功能，通过搭配高分辨率的可见光相机，能实现更高精度的三维重建，因此适用范围更高；通过结合几何与阴影进行深度图优化，本文提出的方法有效提高了热像图三维重建模型的稠密度，解决了Yamaguchi等人[7]的方法中三维模型过于稀疏的问题。

但是，由于本方法需要重建稠密的深度图，引入了全局最优化算法，导致算法耗时较长，因此无法支持实时在线的三维立体重建。若要实现在线重建，可根据稀疏三维点进行Delaunay三角剖分算法构建三维模型，省去深度图的重建，但这种做法产生的模型三角面数量可能较低，需要进一步优化提高精度。