CN104036542B - 一种基于空间光线聚集性的像面特征点匹配方法 - Google Patents

Abstract

本发明提出了一种基于空间光线聚集性的像面特征点匹配方法，在唯一的三维空间，根据像面特征点对应重建光线的聚集性分析，实现特征点匹配，从而利用被测空间的唯一性，将二维的匹配问题回归到三维空间解决，包括如下步骤：(1)空间光线重建；(2)确定光线聚集阈值；(3)光线聚集性判断；(4)像面特征点匹配；(5)同名点合并。

Description

一种基于空间光线聚集性的像面特征点匹配方法

技术领域

本发明涉及计算机视觉领域，特别涉及一种基于空间光线聚集性的像面特征点匹配方法。

背景技术

特征点匹配是计算机视觉中的关键步骤，在三维重建、运动估计、图像检索、摄像机标定等领域中有着重要的应用。大尺寸数字摄影测量中，特征点匹配对于测量系统的精度、可靠性、自动化等方面都有重要影响。但由于拍摄时间、角度、环境的变化、多种传感器的使用和传感器本身的缺陷，使拍摄的图像不仅受噪声的影响，而且存在严重的灰度失真和几何畸变。在这种条件下，匹配算法如何达到精度高、匹配正确率高、速度快、鲁棒性和抗干扰性强成为人们追求的目标。

在摄影测量系统中，通常会采用两种带有自反射特性的特征点：一种是带有编码信息的编码特征点(图1(a))，编码有不同的设计方案，但在一次测量任务中，采用的编码点通常是唯一的；一种是没有任何编码信息的普通特征点(图1(b))，一般为圆形。

可以看出，对于编码特征点，只要正确识别了该特征点成像的编码信息，就可以对相同编码信息的编码特征点直接匹配。而对于普通特征点，由于其本身不携带任何信息，其成像后不能直接匹配，通常利用成像几何关系、空间分布位置、其他特征约束等条件进行匹配。因此，普通特征点的匹配相对于编码特征点而言，具有较高的难度。

对于目前已有的像面特征点匹配方法，根据其基本原理大致可以分为基于灰度相关的匹配、基于特征的匹配、基于模型的匹配、基于变换域的匹配以及基于几何关系的匹配等。“基于灰度相关的匹配”能够获得较高的定位精度，但是它计算量大，难以达到实时性要求，容易受到光 源的影响，而且存在受成像畸变的影响比较大、待匹配单元的大小不容易确定等问题。“基于特征的匹配”对于图像畸变、噪声、遮挡等具有一定的鲁棒性，由于特征仅仅是源图像数据的部分信息，所以后续的匹配过程处理的数据量会明显少于区域匹配，故特征点匹配相对比区域匹配的速度要快，并且由于特征不是直接依赖于灰度，所以对图像噪声具有抵抗性。“基于模型的匹配”在计算机视觉、摄影测量和模式识别等领域中的应用非常广泛，它可以分为刚体形状匹配和变形模板匹配两大类；频域匹配技术对噪声有较高的容忍程度，可处理图像之间的旋转和尺度变化。

从上述传统的像面特征点匹配方法来看，解决问题的思路多定位于二维像面空间的直接解决方式，即在像面上通过灰度、特征、几何关系等约束条件进行特征点匹配，但由于摄影测量本身以及任务条件的复杂性，使得传统方法难以具有普遍适应性。

发明内容

本发明提供一种基于空间光线聚集性的像面特征点匹配方法，包括如下步骤：

1).布置测量现场：根据被测目标对象，设置编码特征点、普通特征点、定向参照物、长度基准尺，其中分别布置编码特征点和普通特征点这两类特征点；

2).多次成像并进行图像处理和特征点识别：通过同一数字相机对被测空间内的目标对象进行多次成像，其中数字相机可移动，以从不同角度拍摄多幅图像，采用数字图像处理技术对像面进行处理，进行特征点中心定位，并识别其中的编码信息；

3).进行图像的空间定向：对拍摄的每幅图像，根据定向参照物的已知三维信息与对应像面信息，实现每幅图片的空间定向，即获取拍摄每一幅图像时相机拍摄的外方位参数；

4).匹配编码特征点：识别编码特征点并且利用编码信息自动匹配不同图像之间的编码特征点，然后利用匹配的这部分编码特征点，进行初步光束平差优化；

5).匹配普通特征点：利用获得的所有图片参数，对于每个普通特征点的像点坐标，根据空间成像几何关系，重构每个成像点对应的空间光线，然后对于任意一条空间光线，遍历其余未匹配的空间光线，根据预定阈值获得光线在空间的汇聚性，将找到的所有汇聚光线对应的像点记入匹配关系矩阵，并将其标记为已匹配点；

6).数值解算：在成功实现编码特征点和普通特征点的匹配、建立像面参数的基础上，利用光束平差优化算法，实现所有特征点空间坐标、相机参数的高精度解算，

其中在步骤5)中，采用如下两种阈值来判断光线在空间的汇聚性：

①以同一特征点的所有成像光线在该点附近形成的预定离散范围作为分散性阈值T₁；

②以空间光线之间的预定距离作为距离阈值T₂。

优选的，在所述步骤5)中，重构每个成像点对应的空间光线的步骤为：对于像面上的像点p₁(x₁,y₁)，其在相机坐标下的坐标是经过相机外方位参数构成的旋转、平移关系，将p₁′变换到空间坐标系下：

$\left[\begin{array}{c}{X}_1\\ Y_1\\ Z_1\end{array}\right]=\left[\begin{array}{ccc}{a}_1& a_2& a_3\\ b_1& b_2& b_3\\ c_1& c_2& c_3\end{array}\right]\left[\begin{array}{c}{x}_1\\ y_1\\ -f\end{array}\right]+\left[\begin{array}{c}{X}_{s1}\\ Y_{s1}\\ Z_{s1}\end{array}\right]...\left(1\right),$

其中，a_i,b_i,c_i(i＝1,2,3)是相机光轴方位在空间坐标系中的三个角度量ω，，κ的余弦组合，(X_s1,Y_s1,Z_s1)是相机投影中心在空间坐标系中的坐标，

如此，通过如下两点式直线方程重构经像点p₁和投影中心的光线：

$\frac{X_{s1}-X_1}{X_{s1}-X_1}=\frac{Y_{s1}-Y_1}{Y_{s1}-Y_1}=\frac{Z_{s1}-Z_1}{Z_{s1}-Z_1}...\left(2\right).$

优选的，在所述步骤5)中，确定分散性阈值T₁和空间光线距离阈值T₂的过程如下：

(1)取得所有编码点及其对应空间光线；

(2)选择一个编码点CP，计算其所有空间光线之间的距离和交点；

(3)为了保证统计数据的有效性，根据距离剔除其中的粗大误差；

(4)统计编码点CP的空间光线分散范围和光线距离的最大值；

(5)重复(2)，直到所有编码点统计结束；

(6)选择所有编码点统计结果中分散范围和光线距离的最大值分别作为分散性阈值T₁和空间光线距离阈值T₂。

优选的，在所述步骤5)中，判断光线聚集性的步骤包括：

(1)初始化：设置光线分散性阈值T₁、光线距离阈值T₂、候选匹配点集C、匹配关系矩阵M；

(2)从所有图片中，选择一幅带有未匹配特征点的图像，记为Img₁；

(3)从Img₁中选择任意一个未匹配特征点p₁及其对应重建空间光线L₁；

(4)除了Img₁以外的所有图片中，如果有未匹配的特征点，则该图片记为Img₂；

(5)遍历Img₂中所有未匹配特征点，并分别计算对应空间光线与p₁对应空间光线之间的距离d和异面直线公垂线中点P，如果d≤T，则将对应特征点计入p₁的候选匹配点集C，同时记录距离d和中点P；

(6)重复(4)～(5)，直到所有图片完成遍历；

(7)对集合C中的所有候选像点，根据其对应距离d和中点P，判断所有候选光线的聚集性；

(8)将该次找到的所有汇聚光线对应的像点记入匹配关系矩阵M，并将其标记为已匹配点；

(9)返回(2)，重复上述步骤，直到没有未匹配点，

其中，在上述步骤(7)中，判断所有候选匹配点对应光线汇聚性的步骤包括：

(1)对集合C中所有候选点对应的与L₁的公垂线中点，计算其相互之间的距离；

(2)对每个中点统计离其距离小于光线距离阈值T₂的点数；

(3)选择点数最多的中点P以及离其距离小于T的其余中点，构成点集C_m；

(4)计算点集C_m中所有点的中心P_m，即空间坐标的平均值；

(5)对集合C中所有候选点对应的公垂线中点，计算其与P_m之间的距离，如果小于分散性阈值T₁，则其对应光线确定为汇集光线。

优选的，在所述步骤5)中，还包括同名点合并步骤，所述步骤包括：

(1)初始化：设置同名点空间最小距离阈值T_d、所有特征点标记为 未分组、建立所有三维特征点关系矩阵M_g；

(2)根据矩阵M中已有的匹配关系，结合相机参数，利用前方交会法，计算所有具有匹配关系的特征点三维坐标，点数记为n；

(3)计算任意两个特征点p_i和p_j之间三维空间距离，如果距离超过阈值T_d，则矩阵M_g(i,j)和M_g(j,i)置0，否则置1；

(4)遍历所有三维特征点，如果该点未标记分组，则建立新的分组G，将该点计入G，并置已分组标记；

(5)根据关系矩阵M_g，将所有对应置1的点计入分组G，并置已分组标记；

(6)对新计入分组G的特征点，反复执行(5)，直到没有新的点计入分组G为止；

(7)重复(4)～(6)，直到没有未分组特征点为止。

应当理解，前述大体的描述和后续详尽的描述均为示例性说明和解释，并不应当用作对本发明所要求保护内容的限制。

附图说明

参考随附的附图，本发明更多的目的、功能和优点将通过本发明实施方式的如下描述得以阐明，其中：

图1(a)和图1(b)示出了两种常用的特征点样式。

图2示出了根据本发明的光学成像几何关系的线性模型。

图3示出了根据本发明的双目立体视觉系统模型图。

图4示出了根据本发明的摄影测量系统的流程图。

图5示出了根据本发明的普通特征点匹配的流程图。

图6示出了根据本发明的实验用大型室内三维控制场

具体实施方式

通过参考示范性实施例，本发明的目的和功能以及用于实现这些目的和功能的方法将得以阐明。然而，本发明并不受限于以下所公开的示范性实施例；可以通过不同形式来对其加以实现。说明书的实质仅仅是帮助相关领域技术人员综合理解本发明的具体细节。

在下文中，将参考附图描述本发明的实施例。在附图中，相同的附图标记代表相同或类似的部件，或者相同或类似的步骤。

本发明主要解决的是普通特征点的成像匹配问题。

本发明的构思依据于摄影测量系统的原理：像面信息的来源是三维被测空间，像面特征点与空间特征点之间满足共线条件。唯一的三维空间成像为多个二维像面空间，这是匹配问题的来源，也产生了特征点匹配关系的复杂性。因此，为了避免直接的复杂二维像面特征点匹配，可以利用被测空间的唯一性，通过将二维的匹配问题回归到三维空间解决，可以大大降低匹配问题的复杂性。

具体的，在理想成像条件下，每个像点都是由空间特征点反射光经透镜中心到达像面，反过来，同一个特征点的所有成像，从各个像面沿光线反方向投射到空间，将汇聚到同一个空间特征点。因此，利用每个像点重构的光线在三维空间的聚集性，可以判断像面特征点之间的对应关系，也就是汇聚到同一空间特征点的所有光线对应像点即为相互匹配特征点。

下面详述本发明的光学成像几何模型以及测量过程。

1.成像几何模型

1.1相机线性成像模型

相机模型是本发明的光学成像几何关系的最简单的模型，为线性模型，或称为针孔模型(如图2所示)。空间中任何一点P在图像上的成像位置可以用针孔模型近似表示，即任何点在图像上的投影位置p为光心O和P点的连线OP与图像平面的交点。这种关系也称为中心投影或者透视投影。

该模型用到了三个坐标系：物空间坐标系、像平面坐标系和相机坐标系。在物空间选择了一个基准坐标系，用来描述空间中任意物体的位置。可用来确定相机的位置，也可以描述空间点的位置，该坐标系即称为物空间坐标系。它由X_W,Y_W,Z_W轴组成。

在图2中，光轴Z_C与像平面垂直，相交于O_i点，O点(也即O_C点) 称为相机光心，O_CO_i为相机焦距。空间物体经过相机投影在一平面，该平面即为像平面，在像平面上选一坐标系用来描述像点的位置，以相机的光轴Z_C与像平面的交点O_i为像坐标系的原点。由点O_C与X_C,Y_C,Z_C轴组成的直角坐标系称为相机坐标系，X_C和Y_C轴与像平面坐标X_i轴和Y_i轴平行。在实际的测量系统中，相机采集的图像经过图像采集卡变成数字图像存储在计算机的外存中。图像由行和列的像素组成，这里定义一个直角坐标系u,v来描述像素在图像中的位置，如(u,v)表示像素的列数u，行数v，其中u,v都是非负整数。

上面提及的像坐标系是O_iX_iY_i是以毫米为单位描述像点的位置，且其X_i轴与u轴平行，Y_i和与v轴平行。所以，在描述像平面时，有两个坐标系，(u,v)表示像素单位得像平面坐标系，x,y表示毫米的像平面坐标系。并且假设一个像素在x轴与y轴方向上的长度分别为dx和dy毫米。

若O_i在u,v坐标系中的坐标为(u₀,v₀)，则图像中任意一个像素坐标在两个坐标系下的坐标转化关系如下：

$u=\frac{x}{\mathrm{dx}}+u_0...\left(1\right)$

$v=\frac{y}{\mathrm{dy}}+v_0$

也可以表示为：

$\left[\begin{array}{c}u\\ v\\ 1\end{array}\right]=\left[\begin{array}{ccc}\frac{1}{\mathrm{dx}}& 0& u_0\\ 0& \frac{1}{\mathrm{dy}}& v_0\\ 0& 0& 0\end{array}\right]\•\left[\begin{array}{c}x\\ y\\ 1\end{array}\right]...\left(2\right)$

逆关系可写成：

$\left[\begin{array}{c}x\\ y\\ 1\end{array}\right]=\left[\begin{array}{ccc}\mathrm{dx}& 0& -u_0\mathrm{dx}\\ 0& \mathrm{dy}& -\mathrm{vdy}\\ 0& 0& 1\end{array}\right]\•\left[\begin{array}{c}u\\ v\\ 1\end{array}\right]...\left(3\right)$

用旋转矩阵R与平移量t可描述相机坐标系与物空间坐标系之间的转换关系。假设点P在物空间坐标系的坐标是(X_w,Y_w,Z_w,1)^T，而在相机坐标系下的坐标是(X_c,Y_c,Z_c,1)^T，于是其转换关系如下：

$\left[\begin{array}{c}{X}_c\\ Y_c\\ Z_c\\ 1\end{array}\right]=\left[\begin{array}{cc}R& t\\ 0^T& 1\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]=M_1\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]...\left(4\right)$

其中，R为3×3正交单位矩阵；t为三维平移量；M₁为4×4矩阵。

空间任何一点P在图像上的投影位置P_i为光心O与P点连线OP与图像平面的交点。由比例关系有如下关系式：

$x=\frac{{\mathrm{fX}}_c}{Z_c}...\left(5\right)$

$y=\frac{{\mathrm{fY}}_c}{Z_c}$

其中，(x,y)为P_i点图像坐标；(X_c,Y_c,Z_c)为空间点P在相机坐标系下的坐标。上述关系表示为：

$Z_c\left[\begin{array}{c}x\\ y\\ 1\end{array}\right.==\left[\begin{array}{cccc}f& 0& 0& 0\\ 0& f& 0& 0\\ 0& 0& 1& 0\end{array}\right]\left[\begin{array}{c}{X}_c\\ Y_c\\ Z_c\\ 1\end{array}\right]...\left(6\right)$

将式(3)与式(4)代入上式，得到空间点P坐标(X_w,Y_w,Z_w,1)^T与其投影点P_i的坐标(u,v)的关系如下式所示：

$\begin{array}{c}{z}_c\left[\begin{array}{c}u\\ v\\ 1\end{array}\right]=\left[\begin{array}{ccc}\frac{1}{\mathrm{dx}}& 0& u_0\\ 0& \frac{1}{\mathrm{dy}}& 0\\ 0& 0& 1\end{array}\right]\left[\begin{array}{cccc}f& 0& 0& 0\\ 0& f& 0& 0\\ 0& 0& 1& 0\end{array}\right]\left[\begin{array}{cc}R& t\\ 0^T& 1\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]\\ =\left[\begin{array}{cccc}{a}_x& 0& u_0& 0\\ 0& a_y& v_0& 0\\ 0& 0& 1& 0\end{array}\right]\left[\begin{array}{cc}R& t\\ 0^T& 1\end{array}\right]\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\\ 1\end{array}\right]=M_1{M}_2{\stackrel{\‾}{X}}_w=M{\stackrel{\‾}{X}}_w\end{array}=...\left(7\right)$

式中，α_x＝fdx,α_y＝fdy，M称为投影矩阵，它是一个3×4矩阵。

从上式可知，如果已知空间某点P在某一相机的像平面上的像坐标P_i(u,v)，且投影矩阵也是已知的，也无法求解空间点的坐标(X_w,Y_w,Z_w)。因此，用一个像平面上的像点无法确定一个空间点的坐标，至少需要用两个像平面上的点才可唯一确定一个空间点。

因此，在本发明的测量系统中采用的是单相机多成像，它们可同时获得空间点(或物体)的多幅图像，即可获得一空间点的多个像点的二维坐标，然后由式(7)可以求解该物体的三维信息。

在讨论线性模型的时候，认为透镜是理想透镜，不带有任何的畸变。但是在实际的情况下，透镜是有畸变的，因此在精密测量时需要对由于畸变带来的误差进行补偿。

1.2立体视觉

由于至少需要用两个像平面上的点才可唯一确定一个空间点，因此根据本发明的立体视觉系统提供由多幅图像获取物体三维几何信息的方法。

在图3所示的双目立体视觉系统中，两个相机坐标系分别为O_c1X_c1Y_c1Z_c1和O_c2X_c2Y_c2Z_c2。空间点P通过光线O_c1P成像于左像面上一点p₁，通过光线O_c2P成像于右像面上一点p₂，p₁和p₂是一对同名点。平面PO_c1O_c2与左右两个像平面分别交于直线l₁和l₂。由于p₁的同名点p₂既位于右像平面上，又位于平面PO_c1O_c2上，因此p₂必位于PO_c1O_c2于右像平面的交线l₂上；同理，p₂的同名点p₁必位于交线l₁上。l₁称为右图上对应于p₁点的极线，l₂称为左图上对应于p₂点的极线。随着空间点P位置的变化，像点和对应的极线在图像上的位置和角度也发生变化，但是，由于所有的PO_c1O_c2平面都相交于直线O_c1O_c2，而O_c1O_c2交两个像平面于固定两点e₁和e₂，故左像平面上所有的极线相交于e₁，右像平面上所有极线相交于e₂。e₁是右相机光心O_c2在左像面的像点，叫做左极点；e₂是左相机光心O_c1在右像面的像点，叫做右极点。这就是极线几何约束条件，也是传统特征点匹配方法中常用的基本约束条件。

图3中，p₁′和p₂′是由于系统畸变等原因造成的实际成像点，可以看出，由于像点误差的存在，实际成像点偏离理想成像位置，使得同名点的空间光线在空间中无法交于一点，而是构成异面直线，光线之间存在一定的距离D。

根据式(7)，一个双目视觉系统的两个成像单元各自的透视投影方程为：

Z_c1u₁＝M₁X_P＝(M₁₁m₁)X_P…(8)

Z_c2u₂＝M₂X_P＝(M₂₁m₂)X_P…(9)其中，X_P是一空间点P在世界坐标系下的齐次坐标；u₁和u₂是分别是P在两个成像系统下的像p₁和p₂的齐次图像坐标；投影矩阵M被分为两部分，3维列向量m_i(i＝1,2)表示M的最后一列；M_j1(j＝1,2)表示投影矩阵左边3×3的矩阵。

将X_P＝(X_WP,Y_WP,Z_WP,1)^T记作X_P＝(X^T1)^T，其中X＝(X_WP,Y_WP,Z_WP)^T，则上面两式可展开为：

Z_c1u₁＝M₁₁X+m₁…(10)

Z_c2u₂＝M₂₁X+m₂…(11)

将上式消去X得：

Z_c2u₂-Z_c1M₂₁M₁₁ ^-1u₁＝m₂-M₂₁M₁₁ ^-1m₁…(12)

将上式等号右端的向量记作m，即：

m＝m₂-M₂₁M₁₁ ^-1m₁…(13)

将m的反对称矩阵记作[m]_×，并用它去乘式(12)的两端，由于[m]_××m＝0，得：

[m]_×(Z_c2u₂-Z_c1M₂₁M₁₁ ^-1u₁)＝0…(14)

将上式两端除以Z_c2，并且记得到：

[m]_×Z_cM₂₁M₁₁ ^-1u₁＝[m]_×u₂…(15)

上式等号右端的向量[m]_×u₂＝m×u₂，该向量与u₂正交，将u₂ ^T左乘上式两端，并将所得等式两边除以Z_c得到如下结果：

u₂ ^T[m]_×M₂₁M₁₁ ^-1u₁＝0…(16)

式(16)给出了对应物空间同一点P的同名像点u₁和u₂之间必须满足的关系。可以看出，在给定u₁的情况下，式(16)是一个关于u₁横、纵坐标关系的线性方程，即对应于u₁在像平面I₂上的极线；反之，在给出u₂的情况下，式(16)是一个关于u₁横、纵坐标关系的线性方程，即对应于u₂在像平面I₁上的极线。同时，式(16)还表明了对于已标定的双目立体系统，极线方程仅与投影矩阵M₁和M₂有关。

令F＝[m]_×M₂₁M₁₁ ^-1，则F给出了双目之间的极线约束关系，将式 (16)写作u₂ ^TFu₁＝0，F是立体视觉中的基本矩阵。

2.摄影测量流程

根据上述摄影测量的成像几何模型和立体几何基本原理，设计摄影测量系统的主要流程，如图4所示，包括如下步骤：

1).布置测量现场：根据被测目标对象，设置编码特征点、普通特征点、定向参照物、长度基准尺等，其中分别布置编码特征点和普通特征点这两类特征点；

2).进行图像处理和特征点识别：通过数字相机对被测空间内的目标对象进行多次成像，其中数字相机可移动，以从不同角度拍摄多幅图像。采用数字图像处理技术对像面进行处理，进行特征点中心定位，并识别其中的编码信息；

3).进行图像的空间定向：对拍摄的每幅图像，根据定向参照物的已知三维信息与对应像面信息，利用后方交会算法(见下文)，可以实现每幅图片的空间定向，即获取拍摄每一幅图像时相机拍摄的外方位参数；

4).匹配编码特征点：识别编码特征点并且利用编码信息自动匹配不同图像之间的编码特征点，然后利用匹配的这部分编码特征点，进行初步光束平差优化，以提高相机参数的精度，减少由于相机镜头畸变、图像处理带来的误差对后续匹配的影响；

5).匹配普通特征点：利用获得的所有图片参数，结合每个普通特征点的像点坐标，构建成像几何关系，通过极线等约束条件搜索可能的匹配点；

6).数值解算：在成功实现编码特征点和普通特征点的匹配、建立像面参数的基础上，利用光束平差优化算法，实现所有特征点空间坐标、相机参数等的高精度解算，这之前的计算一般是部分点的，在这一步骤中，将所有点都参与运算，精度可以提高。另外，还可以在这一步上加入标准距离约束，进一步提高精度。

3.普通特征点匹配的流程

对上述流程中第5步骤中的普通特征点匹配的流程具体设计如图5所示，该匹配方法需要解决两个主要问题：空间光线的重建和空间光线聚集性判断。

3.1空间光线的重建

3.1.1方位参数初始化

1)后方交会

根据成像几何关系，重建每个成像点对应空间光线的前提条件是已知相机内、外方位参数和像点坐标。其中，相机内参数可以通过实验条件独立标定，在测量过程中作为初始参数；像点坐标通过数字图像处理和亚像素中心定位技术获得，普遍可以达到1/20～1/50像素精度；外方位参数是与相片拍摄时刻相机站位姿态相关，因此，需要根据现场参照物及其对应像面信息解算获得，常用的方法是后方交会法。

2)初步光束平差

受现场参照物限制，后方交会法标定的相机外参数精度不高，会对系统的后续匹配、参数优化等过程产生较大影响。测量过程中，编码点的作用除了用于大视场拼接外，还可以利用其已知匹配对应关系，进行初步光束平差优化，即相机内、外参数和编码点空间坐标等根据成像模型获得最优解的过程，以进一步提高相机参数和坐标的精度。光束平差法的本质是最优迭代问题，在进行光束平差之前，待优化变量中相机内、外方位参数的初始化已经建立，而编码点空间点坐标可以根据相机参数和其多幅图像上的像点坐标，采用多条光线前方交会法计算获得。

3.1.2空间光线方程重构

在空间重构光线关系中，匹配的准确率受空间光线重构精度影响，而空间光线重构精度取决于像面特征点定位精度和相机参数，此外，镜头畸变也是造成光线重构误差的重要因素。因此，为了提高匹配的准确 率，在初步优化相机参数的基础上，需要利用畸变模型参数对特征点成像进行畸变校正，进一步提高光线重构精度。由于在后续匹配时需要多次遍历所有未匹配像点，所以，为了提高匹配速度，降低运算量，在匹配之前对所有像点进行集中前期校正，并依此作为后续匹配处理用像点。

如图3所示，两个相机坐标系分别为O_c1X_c1Y_c1Z_c1和O_c2X_c2Y_c2Z_c2。空间点P通过光线O_c1P成像于左像面上一点p₁，通过光线O_c2P成像于右像面上一点p₂，p₁和p₂是一对同名点。平面PO_c1O_c2交左右两个像平面分别于直线l₁和l₂。由于p₁的同名点p₂既位于右像平面上，又位于平面PO_c1O_c2上，因此p₂必位于PO_c1O_c2于右像平面的交线l₂上；同理，p₂的同名点p₁必位于交线l₂上。l₂称为右图上应于p₁点的极线，l₂称为左图上对应于p₂点的极线。随着空间点P位置的变化，像点和对应的极线在图像上的位置和角度也发生变化，但是，由于所有的PO_c1O_c2平面都相交于直线O_c1O_c2，而O_c1O_c2交两个像平面于固定两点e₁和e₂，故左像平面上所有的极线相交于e₁，右像平面上所有极线相交于e₂。e₁是右相机光心O_c1在左像面的像点，叫做左极点；e₂是左相机光心O_c2在右像面的像点，叫做右极点。这就是极线几何约束条件，也是传统特征点匹配方法中常用的基本约束条件。

通过像面I和II中的像点p₁和p₂来解算空间点P的空间坐标。以像面I为例，要得到每个像点重构光线在空间中的方程，可以通过投影中心和像点在世界坐标系下的坐标，由两点法来完成。以像面上的像点p₁(x₁,y₁)为例，在相机坐标下的坐标是经过相机外方位参数构成的旋转、平移关系，可以将p₁′变换到物空间坐标系下：

$\left[\begin{array}{c}{X}_1\\ Y_1\\ Z_1\end{array}\right]=\left[\begin{array}{ccc}{a}_1& a_2& a_3\\ b_1& b_2& b_3\\ c_1& c_2& c_3\end{array}\right]\left[\begin{array}{c}{x}_1\\ y_1\\ -f\end{array}\right]+\left[\begin{array}{c}{X}_{s1}\\ Y_{s1}\\ Z_{s1}\end{array}\right]...\left(17\right)$

其中，a_i,b_i,c_i(i＝1,2,3)是相机光轴方位在物空间坐标系中的三个角度量ω，，κ的余弦组合，(X_s1,Y_s1,Z_s1)是相机投影中心在物空间坐标系中的坐标。所以，经像点p1和投影中心的光线通过两点式直线方程可以建立：

$\frac{X_{s1}-X_1}{X_{s1}-X_1}=\frac{Y_{s1}-Y_1}{Y_{s1}-Y_1}=\frac{Z_{s1}-Z_1}{Z_{s1}-Z_1}...\left(18\right)$

同样，为了减少运算量，提高处理速度，避免重复计算光线方程， 在匹配之前集中建立所有空间光线方程，并存储。

3.2、空间光线聚集性判断

3.2.1、光线聚集阈值

光线在三维空间的汇聚性在理想条件下表现为共点，即所有匹配光线应交会于同一空间点。但实际情况下，由于镜头畸变、相机标定误差等因素的影响，空间光线的重建存在一定的误差，同名特征点的所有匹配光线并不能完全交于其对应空间特征点，而且相互之间也往往不相交，存在一定的距离，构成空间异面直线。

判断光线在空间中是否汇聚的方法有两种：

一种是判断同名光线之间的汇聚分散性阈值，即同一特征点的所有成像光线应该在该点附近形成汇聚，即带有一定离散范围的汇聚。因此，设置该离散范围的阈值，称之为分散性阈值，记为T₁。

另一种是判断空间光线之间的距离阈值，即在进行匹配时先挑选满足一定条件的特征点作为候选匹配点，并选用空间光线之间的距离作为候选条件，因此，设置空间光线距离阈值，记为T₂。

为了实现测量系统的自动化，要使T₁和T₂能够根据不同的测量任务需求自动确定。

在前面的相机标定过程中，为了提高标定初始参数的精度，利用编码点的已知匹配关系，进行了初步光束平差优化。同样，自动确定分散性阈值和空间光线距离阈值同样可以利用编码点。在一次测量任务中，编码点的布设往往是全视场内分布的，因此，通过统计编码点所有重建光线在三维空间中的汇聚范围和光线距离，可以代表同样测量任务中其余普通特征点重建光线的分散性和距离，从而实现阈值的自动确定。

确定分散性阈值T₁和空间光线距离阈值T₂的过程如下：

(1)取得所有编码点及其对应空间光线；

(2)选择一个编码点CP，计算其所有空间光线之间的距离和交点；

(3)为了保证统计数据的有效性，根据距离剔除其中的粗大误 差；

(4)统计编码点CP的空间光线分散范围和光线距离的最大值；

(5)重复(2)，直到所有编码点统计结束；

(6)选择所有编码点统计结果中分散范围和光线距离的最大值分别作为分散性阈值T₁和空间光线距离阈值T₂

(7)结束。

3.2.2、光线聚集性

所有空间重建光线的聚集性判断思路是，对于任意一条空间光线，遍历其余未匹配空间光线与其距离，对所有满足空间光线距离阈值的光线，再对其光线公垂线中点之间的距离，利用光线分散性阈值T₁获得中点的聚集性，从而表现为对应空间光线的聚集性，也即代表了其像面特征点的匹配关系。

光线聚集性判断的具体流程：

(1)初始化：设置光线分散性阈值T₁、光线距离阈值T₂、候选匹配点集C、匹配关系矩阵M；

(2)从所有图片中，选择一幅带有未匹配特征点的图像，记为Img₁；

(3)从Img₁中选择任意一个未匹配特征点p₁及其对应重建空间光线L₁；

(4)除了Img₁以外的所有图片中，如果有未匹配的特征点，则该图片记为Img₂；

(5)遍历Img₂中所有未匹配特征点，并分别计算未匹配特征点所对应的空间光线与p₁对应的空间光线之间的距离d和异面直线公垂线中点P，如果d≤T，则将对应特征点计入p₁的候选匹配点集C，同时记录距离d和中点P；

(6)重复(4)～(5)，直到所有图片完成遍历；

(7)对集合C中的所有候选像点，根据其对应距离d和中点P，判断所有候选光线的聚集性(具体方法见下)；

(8)将该次找到的所有汇聚光线对应的像点记入匹配关系矩阵 M，并将其标记为已匹配点；

(9)返回(2)，重复上述步骤，直到没有未匹配点；

(10)结束。

所有候选匹配点对应光线汇聚性的判断方法：

(1)对集合C中所有候选点对应的与L₁的公垂线中点，计算其相互之间的距离；

(2)对每个中点统计离其距离小于光线距离阈值T₂的点数；

(3)选择点数最多的中点P以及离其距离小于T的其余中点，构成点集C_m；

(4)计算点集C_m中所有点的中心P_m，即空间坐标的平均值；

(5)对集合C中所有候选点对应的公垂线中点，计算其与P_m之间的距离，如果小于分散性阈值T₁，则其对应光线确定为汇集光线；

(6)结束。

3.2.3、同名点合并

通过空间光线聚集性确定的特征点匹配关系矩阵M中，基本完成了特征点匹配要求，但还有一些特殊情况需要考虑和处理，例如，在匹配过程中，受阈值限制、像面处理精度、拍摄质量等因素的影响，某些特征点的空间光线离线范围较大，超出了设定阈值的范围。按照上述匹配过程，在这种情况下，同名特征点有可能被匹配成多个点，每个点由部分全部空间光线汇集而成。针对该情况，需要对同名特征点合并，为此，需要建立同名点空间最小距离阈值T_d，该阈值表示空间非同名点的最小距离，即小于该距离阈值的空间点即视为同名点，其值通常可以由手工或经验确定。

同名点合并过程如下：

(1)初始化：设置同名点空间最小距离阈值T_d、所有特征点标记为未分组、建立所有三维特征点关系矩阵M_g；

(2)根据矩阵M中已有的匹配关系，结合相机参数，利用前方交 会法，计算所有具有匹配关系的特征点三维坐标，点数记为n；

(3)计算任意两个特征点p_i和p_j之间三维空间距离，如果距离超过阈值T_d，则矩阵M_g(i,j)和M_g(j,i)置0，否则置1；

(4)遍历所有三维特征点，如果该点未标记分组，则建立新的分组G，将该点计入G，并置已分组标记；

(5)根据关系矩阵M_g，将所有对应置1的点计入分组G，并置已分组标记；

(6)对新计入分组G的特征点，反复执行(5)，直到没有新的点计入分组G为止；

(7)重复(4)～(6)，直到没有未分组特征点为止；

(8)结束。

4.实验与结果分析

为了验证基于空间聚集性的特征点匹配方法的有效性，将其应用于大尺寸静态近景摄影测量系统中，进行特征点匹配实验，见图6。

实验对象是大型室内三维控制场，该控制场尺寸5m×2.5m×1.5m，主要由框架、三节臂、测点转换附件等组成，底部设计隔震层以保证其稳定性。三节臂长度可调，顶端安装测点转换附件，可以实现摄影测量系统与激光跟踪仪之间的测点转换，便于进行测量实验、精度评价等。

对于摄影测量任务而言，该控制场的结构不同于一般的连续被测物体表面，属于比较复杂的测量对象，在不同位置和角度下的特征点成像变化很大，甚至会产生不同点之间成像遮挡、重合等极端情况，大大增加了匹配的难度。因此，用该控制场进行匹配方法的有效性检验是合适的。

控制场中满视场布置了普通特征点和编码点，并辅以定向参照物、长度基准尺等附件。采用Nikon D2X非量测数字相机进行图像采集，镜头主距20mm，像面分辨率1200万。为了充分验证该匹配方法的有效性，围绕控制场进行10组测量成像，每组特征点布置方案不同，采集大约50～70幅图片，并且由不同人员完成。

通过本发明的基于空间光线聚集性的特征点匹配方法进行普通特征点匹配，实验结果如表1所示：

表1匹配实验结果

表1中的特征点指的是普通特征点，其中，特征点识别率表示全局特征点正确识别率，像面特征点识别正确率表示所有像面上特征点匹配准确情况，而漏匹配率则表示像面上特征点被剔除的情况，漏匹配率的相关统计信息见表2。

表2漏匹配率统计信息表

结合表1和表2的结果可以看出，本发明的基于空间光线聚集性的像面特征点匹配方法很好地解决了全局特征点识别和匹配关系准确性的问题，而漏掉的匹配点，占比不超过总数的0.5％，实验表明，这些漏掉的匹配点不会对全局光束平差优化结果产生显著影响。因此，该匹配方法可以满足摄影测量像面特征点匹配自动化需求。

利用相机标定参数和像面特征点信息重建空间光线，并根据光线在三维空间的聚集性判断相应特征点的匹配关系，可以很好地解决多图像之间的像面特征点匹配问题。其优点在于，将多二维空间问题转变至单一三维空间内处理，也就是将图与图之间的匹配转变为所有光线的集中 聚集性判断，有效降低了问题的复杂性，获得良好的匹配效果。

结合这里披露的本发明的说明和实践，本发明的其他实施例对于本领域技术人员都是易于想到和理解的。说明和实施例仅被认为是示例性的，本发明的真正范围和主旨均由权利要求所限定。

Claims (5)

1.一种基于空间光线聚集性的像面特征点匹配方法，包括如下步骤：

1).布置测量现场：根据被测目标对象，设置编码特征点、普通特征点、定向参照物、长度基准尺，其中分别布置编码特征点和普通特征点这两类特征点；

2).多次成像并进行图像处理和特征点识别：通过同一数字相机对被测空间内的目标对象进行多次成像，其中数字相机可移动，以从不同角度拍摄多幅图像，采用数字图像处理技术对像面进行处理，进行特征点中心定位，并识别其中的编码信息；

3).进行图像的空间定向：对拍摄的每幅图像，根据定向参照物的已知三维信息与对应像面信息，实现每幅图片的空间定向，即获取拍摄每一幅图像时相机拍摄的外方位参数；

4).匹配编码特征点：识别编码特征点并且利用编码信息自动匹配不同图像之间的编码特征点，然后利用匹配的这部分编码特征点，进行初步光束平差优化；

5).匹配普通特征点：利用获得的所有图片参数，对于每个普通特征点的像点坐标，根据空间成像几何关系，重构每个成像点对应的空间光线，然后对于任意一条空间光线，遍历其余未匹配的空间光线，对所有满足空间光线距离阈值的光线，再对其光线公垂线中点之间的距离，利用光线分散性阈值获得中点的聚集性，从而表现为对应空间光线的聚集性，将找到的所有汇聚光线对应的像点记入匹配关系矩阵，并将其标记为已匹配点；

6).数值解算：在成功实现编码特征点和普通特征点的匹配、建立像面参数的基础上，利用光束平差优化算法，实现所有特征点空间坐标、相机参数的高精度解算，

其中在步骤5)中，采用如下两种阈值来判断光线在空间的聚集性：

①以同一特征点的所有成像光线在该点附近形成的预定离散范围作为分散性阈值T₁；

②以空间光线之间的预定距离作为距离阈值T₂。

2.如权利要求1所述的方法，其中，在所述步骤5)中，重构每个成像点对应的空间光线的步骤为：对于像面上的像点p₁(x₁,y₁)，其在相机坐标下的坐标是p′₁(x₁,y₁,-f)，经过相机外方位参数构成的旋转、平移关系，将p′₁变换到空间坐标系下：

$\left[\begin{array}{c}{X}_1\\ Y_1\\ Z_1\end{array}\right]=\left[\begin{array}{ccc}{a}_1& a_2& a_3\\ b_1& b_2& b_3\\ c_1& c_2& c_3\end{array}\right]\left[\begin{array}{c}{x}_1\\ y_1\\ -f\end{array}\right]+\left[\begin{array}{c}{X}_{s1}\\ Y_{s1}\\ Z_{s1}\end{array}\right]...\left(1\right),$

其中，a_i,b_i,c_i(i＝1,2,3)是相机光轴方位在空间坐标系中的三个角度量ω，κ的余弦组合，(X_s1,Y_s1,Z_s1)是相机投影中心在空间坐标系中的坐标，

如此，通过如下两点式直线方程重构经像点p₁和投影中心的光线：

$\frac{X_{s1}-X_1}{X_{s1}-X_1}=\frac{Y_{s1}-Y_1}{Y_{s1}-Y_1}=\frac{Z_{s1}-Z_1}{Z_{s1}-Z_1}...\left(2\right).$

3.根据权利要求1所述的方法，其中在所述步骤5)中，确定分散性阈值T₁和空间光线距离阈值T₂的过程如下：

(1)取得所有编码点及其对应空间光线；

(2)选择一个编码点CP，计算其所有空间光线之间的距离和交点；

(3)为了保证统计数据的有效性，根据距离剔除其中的粗大误差；

(4)统计编码点CP的空间光线分散范围和光线距离的最大值；

(5)重复(2)，直到所有编码点统计结束；

(6)选择所有编码点统计结果中分散范围和光线距离的最大值分别作为分散性阈值T₁和空间光线距离阈值T₂。

4.根据权利要求1所述的方法，其中在所述步骤5)中，判断光线聚集性的步骤包括：

(1)初始化：设置光线分散性阈值T₁、光线距离阈值T₂、候选匹配点集C、匹配关系矩阵M；

(2)从所有图片中，选择一幅带有未匹配特征点的图像，记为Img₁；

(3)从Img₁中选择任意一个未匹配特征点p₁及其对应重建空间光线L₁；

(4)除了Img₁以外的所有图片中，如果有未匹配的特征点，则该图片记为Img₂；

(5)遍历Img₂中所有未匹配特征点，并分别计算对应空间光线与p₁对应空间光线之间的距离d和异面直线公垂线中点P，如果d≤T₂，则将对应特征点计入p₁的候选匹配点集C，同时记录距离d和中点P；

(6)重复(4)～(5)，直到所有图片完成遍历；

(7)对集合C中的所有候选像点，根据其对应距离d和中点P，判断所有候选光线的聚集性；

(8)将该次找到的所有汇聚光线对应的像点记入匹配关系矩阵M，并将其标记为已匹配点；

(9)返回(2)，重复上述步骤，直到没有未匹配点，

其中，在上述步骤(7)中，判断所有候选光线聚集性的步骤包括：

(1)对集合C中所有候选点对应的与L₁的公垂线中点，计算其相互之间的距离；

(2)对每个中点统计离其距离小于光线距离阈值T₂的点数；

(3)选择点数最多的中点P以及离其距离小于T的其余中点，构成点集C_m；

(4)计算点集C_m中所有点的中心P_m，即空间坐标的平均值；

(5)对集合C中所有候选点对应的公垂线中点，计算其与P_m之间的距离，如果小于分散性阈值T₁，则其对应光线确定为汇聚光线。

5.根据权利要求4所述的匹配方法，其中在所述步骤5)中，还包括同名点合并步骤，所述步骤包括：

(1)初始化：设置同名点空间最小距离阈值T_d、所有特征点标记为未分组、建立所有三维特征点关系矩阵M_g；

(2)根据矩阵M中已有的匹配关系，结合相机参数，利用前方交会法，计算所有具有匹配关系的特征点三维坐标，点数记为n；

(3)计算任意两个特征点p_i和p_j之间三维空间距离，如果距离超过阈值T_d，则矩阵M_g(i,j)和M_g(j,i)置0，否则置1；

(4)遍历所有三维特征点，如果该点未标记分组，则建立新的分组G，将该点计入G，并置已分组标记；

(5)根据关系矩阵M_g，将所有对应置1的点计入分组G，并置已分组标记；

(6)对新计入分组G的特征点，反复执行(5)，直到没有新的点计入分组G为止；

(7)重复(4)～(6)，直到没有未分组特征点为止。