基于物体表面形貌的单相机视觉位姿测量方法

1 引言

在风洞实验中常常伴有气流扰动、震动等因素，视觉测量已经成为风洞内物体位姿测量的主要手段^[1]。视觉测量主要分为单相机视觉测量系统与多相机视觉测量系统。单相机测量系统结构简单，标定便利，有着更大的测量范围和更高的测量精度，广泛应用于空间交会、风洞非接触式检测等领域。

通常测量系统的目标可以分为合作目标与非合作目标两类，其中合作目标上具有固定数目、几何关系已知的光学标志点，例如刘巍^[2]设计了使用彩色光学标记点的单目测量系统来测量从飞行器上分离的外挂物的位姿变化。非合作目标存在于很多应用中，比如风洞测试的实心模型、传送带上的工件等。这类物体上无法设置光学标志，无法提供测量标志点。针对此类物体，周扬^[3]提出了基于匹配CAD模型离线库的单相机测量方法来测量工件的位姿，其测量精度在2°左右，其使用物体CAD模型的思路有着一定的参考价值。

为了测量不能设置光学标志点的飞行物体在恒流风场中的位姿变化，本文提出了一种基于物体表面形貌的单相机位姿测量方法。本文使用ORB(oriented FAST and rotation BRIEF)^[4]特征提取方法获取物体的特征角点，并结合物体形貌模型计算特征点的三维坐标，最后根据稳定多点透视算法^[5] (robust solution to the perspective-n-point problem，RPnP)求解物体的位姿变化。

2 位姿变化解算原理

运动中的物体与相机的位姿关系如图 1所示，以物体的初始位置为参考位置，需要测量的物体的位姿变化以旋转矩阵R_F，平移向量t_F表示。物体坐标系OXYZ在前后时刻与相机坐标系O_cX_cY_cZ_c的相对位姿关系以旋转矩阵R，R_r与平移向量t，t_r表示，两者有如下关系式：

\( \begin{array}{*{20}{c}} {{\mathit{\boldsymbol{O}}_{\rm{c}}} = \mathit{\boldsymbol{RO}} + \mathit{\boldsymbol{t}}}\\ {{\mathit{\boldsymbol{O}}_{\rm{c}}} = {\mathit{\boldsymbol{R}}_{\rm{r}}}{\mathit{\boldsymbol{O}}_{\rm{r}}} + {\mathit{\boldsymbol{t}}_{\rm{r}}}}\\ {\left[ {\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{R}}_{\rm{F}}}}\\ {{\mathit{\boldsymbol{t}}_{\rm{F}}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{R}}_{\rm{r}}^{ - 1}\mathit{\boldsymbol{R}}}\\ {\mathit{\boldsymbol{R}}_{\rm{r}}^{ - 1}\left( {\mathit{\boldsymbol{t}} - {\mathit{\boldsymbol{t}}_{\rm{r}}}} \right)} \end{array}} \right],} \end{array} \)(1)

式中：${\mathit{\boldsymbol{O}}_{\rm{c}}} = {\left[ {{X_{\rm{c}}}\; \; {Y_{\rm{c}}}\; \; {Z_{\rm{c}}}} \right]^{\rm{T}}}$为相机坐标，$\mathit{\boldsymbol{O}} = {\left[ {X\; \; Y\; \; Z} \right]^{\rm{T}}}$为运动前物体坐标系；${\mathit{\boldsymbol{O}}_1} = {\left[ {{X_1}\; \; {Y_1}\; \; {Z_1}} \right]^{\rm{T}}}$为运动后物体坐标系。

[Image omitted: See PDF]

根据多点透视原理(perspective-n-point problem，PnP)，相机与物体的相对位姿关系可以通过物体上一定数量的在物体坐标系下的坐标以及图像坐标都已知的特征点来计算。因无法设置光学标志，为了缩减一帧计算时间，本文使用物体的ORB特征角点对作为特征点。特征点在物体坐标系下的坐标通过物体的表面形貌模型求解。在之后的论述里，特征点在物体坐标系下的坐标简称为三维坐标。

3 获得特征点三维坐标 3.1 物体表面形貌模型信息的获取

在实际应用中，有多种方法获得待测目标表面形貌模型信息，比如CAD模型提取和多视图三维重建等。物体在实际加工与使用时会有缺陷以及磨损，多视图三维重建可以方便准确地获取物体表面的形貌信息^[6]。

该方法基于运动恢复结构(structure from motion，SFM)原理^[7]，基本重建过程如图 2所示。获取从多个角度拍摄的物体的图像序列，获取两两图像间特征角点匹配关系；根据匹配关系计算基础矩阵，并在已知重建所用相机内参矩阵的情况下，获得本征矩阵E；对本征矩阵进行SVD分解^[8]获得相机位置参数；使用三角法计算特征点对应的三维坐标，获得初始三维点云；使用光束法平差^[9]对三维点云与相机位置进行优化；最后以基于面元的多视图立体视觉(patch-based multiview stereo, PMVS) ^[10]算法进行稠密重建得到表面三维点云。

[Image omitted: See PDF]

如图 3所示，本文选择SFM多视图三维重建的方式，使用定位精度较高的SIFT特征点作为三维重建点，使用实验室的三维视觉重建系统获得待测飞行器模型的表面形貌三维点云。

[Image omitted: See PDF]

当待测物体表面属于可以解析的结构，可以通过点云拟合得到表面结构的几何表达式^[11]。为了方便之后的理论推导，设待测物体表面几何模型在物体坐标系下的拟合方程为

\( f\left( {X,Y,Z} \right) = 0。\)(2) 3.2 计算测量相机成像模型的参数

相机使用小孔成像模型，物体特征点的图像坐标(u，v)与三维坐标(X，Y，Z)的有如下关系式：

\( {Z_{\rm{c}}}\left[ {\begin{array}{*{20}{c}} u\\ v\\ 1 \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{f_x}}&0&{{u_0}}\\ 0&{{f_y}}&{{v_0}}\\ 0&0&1 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} \mathit{\boldsymbol{R}}&\mathit{\boldsymbol{t}}\\ {\vec 0}&1 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} X\\ Y\\ Z\\ 1 \end{array}} \right], \)(3)

简化为

\( {Z_{\rm{c}}}\left[ {\begin{array}{*{20}{c}} u\\ v\\ 1 \end{array}} \right] = {\mathit{\boldsymbol{M}}_1}{\mathit{\boldsymbol{M}}_2}{\left[ {\begin{array}{*{20}{c}} X&Y&Z&1 \end{array}} \right]^{\rm{T}}}, \)(4)

在成像模型中有两个参数需要计算：相机的内参矩阵M₁与相机坐标系与物体坐标系的位置关系M₂。

相机内参矩阵M₁通过张正友标定法^[11]使用棋盘格标定得到。相机与物体之间的位置关系M₂则使用相机拍摄的物体图像与多视图重建使用的图像集计算。

因三维重建使用SIFT特征点，所以将相机拍摄的物体图像与物体三维重建图像集进行SIFT^[12]特征角点匹配，因物体三维点云中的点与重建所使用的图像集的SIFT特征角点间有着一对多的对应关系，所以物体图像上得到的SIFT特征点根据此对应关系可以得到三维坐标。根据多点透视原理(PnP)，有4个以上确定三维坐标的的SIFT特征点即可使用RPnP算法计算得出相机坐标系与物体坐标系的位置关系M₂。

3.3 使用物体几何模型计算ORB特征点三维坐标

确定M₁与M₂之后，从式(3)推导出以下方程组：

\( \left\{ \begin{array}{l} {M_x}{Z_{\rm{c}}} - {B_x} = X\\ {M_y}{Z_{\rm{c}}} - {B_y} = Y\\ {M_z}{Z_{\rm{c}}} - {B_z} = Z \end{array} \right., \)(5)

其中：

\( \mathit{\boldsymbol{M}} = \left[ {\begin{array}{*{20}{c}} {{M_x}}\\ {{M_y}}\\ {{M_z}} \end{array}} \right] = {\mathit{\boldsymbol{R}}^{ - 1}}\left[ {\begin{array}{*{20}{c}} {\left( {u - {u_0}} \right)/{f_x}}\\ {\left( {v - {v_0}} \right)/{f_y}}\\ 1 \end{array}} \right], \) \( \mathit{\boldsymbol{B}} = \left[ {\begin{array}{*{20}{c}} {{B_x}}\\ {{B_y}}\\ {{B_z}} \end{array}} \right] = {\mathit{\boldsymbol{R}}^{ - 1}}\mathit{\boldsymbol{t}}。\)

式(5)为四元方程组，需要额外的约束条件以获得有限个解。式(2)为待测物体表面的拟合方程可以作为额外条件，则式(5)转化为

\( f\left( {{M_x}{Z_{\rm{c}}} - {B_x},{M_y}{Z_{\rm{c}}} - {B_y},{M_z}{Z_{\rm{c}}} - {B_z}} \right) = 0, \)(6)

问题转变成求式(6)的解集：S{S_i(X_i，Y_i，Z_i，Z_ci)| i=1, ..., n}。方程组(5)的物理意义为由相机坐标系原点O_c与特征点U(u，v)两点确定的一条空间直线。式(6)的解的物理意义为该直线与物体表面的交点S_i，其中Z_ci为交点S_i