一种基于面特征的多平台点云匹配方法与流程

本发明公开了一种基于面特征的多平台点云匹配方法，属于移动测量点云匹配技术领域。

背景技术：

目前，车载、机载、单站lidar扫描仪受到搭载平台姿态、控制点精度、gps精度、惯导精度等影响，所采集到的多源数据坐标系之间存在微小的偏差。上述各平台采集方式均有一定的适用条件，在实际生产应用中，需要多平台点云数据进行匹配融合，以获得目标地物的多维、多时空点云数据。为了充分利用多源数据中的信息，需要匹配多源点云数据的坐标系并根据同名特征(点、线、面)完成点云的拼接。点云拼接的精度影响后续数据处理的精度，因此，点云拼接是数据处理的关键一步。

在体测量模式中，对同一表面的两次测量可能不存在公共点，因此需要研究基于模型或者基于特征的坐标转换方法；平面是一种广泛存在于自然界的规则模型，若能从点云数据中拟合出多个平面，并利用公共平面进行坐标转换，则能避免现有技术中，需要从海量的点云数据中寻找公共点的缺陷。

技术实现要素：

本发明公开了一种基于行驶方向结构特征约束的道路标识线提取方法，以解决现有技术中，需要从海量的点云数据中寻找公共点的问题。

一种基于面特征的多平台点云匹配方法，包括以下步骤：

s1.通过架设在不同平台上的数据采集设备采集目标区域的地物，获得不同平台下的原始点云数据并进行预处理；

s2.通过人机交互的方式，从预处理后的点云数据中均匀地选取三对或三对以上不共面的匹配平面，采用稳健ransac算法进行平面拟合并提取公共面的参数，进行平面点云的自身平整度评价；

s3.基于三对或三对以上对应平面的平面参数，将旋转参数化为非线性模型，进行线性化处理，平移参数不变，采用整体最小二乘法进行参数求解；通过基于平差结果的改正数，对待配准平面的相关性进行评价，对不满足要求的平面进行重新选取；

s4.基于坐标转换模型，完成点云的匹配。

步骤s1在数据处理阶段，将不同平台点云原始数据进行滤波去噪预处理，将预处理后点云坐标系转换到以重心为坐标原点、坐标轴方向不变的局部坐标系；局部坐标原点计算公式如下：其中，[x0，y0，z0]表示局部坐标系原点，[xi，yi，zi]表示第i个点云坐标；变换后点云坐标计算公式如下：[x′，y′，z′]＝[xi，yi，zi]-[x0，y0，z0]，其中，[x′，y′，z′]表示变换后点云坐标，m表示点云数目。

s2包括如下子步骤：

s2.1.设置相关参数：迭代次数n，有效点比例阈值p，有效阈值t；

s2.2.在点云集中随机选择不共线的3个点，计算对应的平面方程，公式如下：a1x+b1y+c1z＝d1，其中，(a1，b1，c1)表达该面的法向量，d1表示原点到平面的距离；

s2.3.计算所有观测值到平面的距离，公式如下：di＝(a1xj+b1yj+c1zj+d1)；如果di≤t，则认为该点为有效点或内点，否则为无效点；

s2.4.计算平面μ1中有效点数所占总点数的比例f，若f≥p，则认为平面μ1为有效平面，否则为无效平面，完成面点云的自身平整度评价；

s2.5.由设置的迭代次数n进行多次抽样，选出有效点占比最高的平面作为最终拟合平面；

s2.6.采用特征值法对s2.5拟合得到的平面进行参数计算，得到平面参数方程如下：ax+by+cz+d＝0，其中(a，b，c，d)为平面参数；

s2.7.在不同平台点云数据上，均匀选取k≥3对互相不平行的对于平面，规定法向量的正方向，计算其平面参数(a，b，c，d)，该平面参数被唯一确定。

步骤s3包括如下子步骤：

s3.1.构建坐标转换模型：

设同一平面在不同坐标系中的表示分别为(a1，b1，c1，d1)，(a2，b2，c2，d2)，则仿射转换模型为：其中，表示旋转矩阵参数，表示平移矩阵参数；

s3.2.旋转矩阵进行线性化，将9参数转换为3参数：

用坐标转换参数替换旋转矩阵和平移矩阵，公式如下：

其中，α为绕x轴的旋转角度，β为绕y轴的旋转角度，γ为绕z轴的旋转角度；

法向量(a2，b2，c2)计算公式如下：

a2＝cosγcosβa1+(cosγsinβsinα+sinγcosα)b1+(sinγsinα-cosγsinβcosα)c1；

b2＝-sinγcosβa1+(cosγcosα-sinγsinβsinα)b1+(sinγsinβcosα+cosγsinα)c1；

c2＝sinβa1-cosβsinαb1+cosβcosαc1；

s3.3.构建间接平差模型：

分别求α、β、γ关于a2、b2和c2的偏导数：

d2不需进行线性化处理，计算公式为：d2＝txa1+tyb1+tzc1+d1；

间接平差模型公式如下：

将上式改写成最小二乘误差方程的形式如下：

v＝ax-l；

其中，x＝(dαdβdγdtxdtydtz)′，

为(a2，b2，c2，d2)的估计值；

s3.4.整体最小二乘求解参数：

假设观测量l中观测误差为e，系数矩阵a中误差为ea，则最小二乘误差方程可写为：

(a+ea)x＝l+e，

即存在式中，e＝[eae]；

此时最优化限制约束条件为：

采用奇异值分解的方法进行求解，对增广矩阵[al]进行奇异值分解：

分解得到的m+1个特征向量组成的正交矩阵，

则参数x的估值可表示为：

s3.5.根据平差结果的改正数对待配准平面的相关性进行评价，改正数越大，则相关性越弱，平面对越差，则需对不满足要求的平面对进行重新选取。

步骤s4包括如下子步骤：

s4.1.根据计算的3个旋转参数，对不同平台的对应点云坐标系进行转换，转换模型为：

s4.2.将转换后的点云由局部坐标系变换到原始坐标系，完成基于面特征的多平台点云匹配，公式如下：[x，y，z]＝[x2，y2，z2]+[x0，y0，z0]，其中，[x，y，z]表示经匹配变换后的最终坐标。

与现有技术相比，本发明的有益效果是：

(1)通过建立局部坐标系，将不同平台的点云坐标转换至各自局部坐标系，避免不同平台点云坐标数值差值过大引起的匹配误差，提高了匹配精度；

(2)将平面转换模型的9参数旋转模型转换为3参数模型，考虑到观测向量和系数矩阵都可能存在误差，采用整体最小二乘法进行参数求解，有效提高了参数求解精度，增加了算法的鲁棒性，为点云匹配提供了有力支撑；

(3)由于本发明点云匹配方法是基于面特征，面的选取和拟合精度决定了匹配精度，因此设置了两种面的评价指标，以选出较优的平面，降低平面误差对匹配精度的影响。

附图说明

图1为基于面特征的多平台点云匹配方法的技术流程图。

具体实施方式

下面结合具体实施方式对本发明作进一步详细说明：

一种基于面特征的多平台点云匹配方法，技术流程图如图1所示，包括以下步骤：

s1.通过架设在不同平台上的数据采集设备采集目标区域的地物，获得不同平台下的原始点云数据并进行预处理；

s2.通过人机交互的方式，从预处理后的点云数据中均匀地选取三对或三对以上不共面的匹配平面，采用稳健ransac算法进行平面拟合并提取公共面的参数，进行平面点云的自身平整度评价；

s3.基于三对或三对以上对应平面的平面参数，将旋转参数化为非线性模型，进行线性化处理，平移参数不变，采用整体最小二乘法进行参数求解；通过基于平差结果的改正数，对待配准平面的相关性进行评价，对不满足要求的平面进行重新选取；

s4.基于坐标转换模型，完成点云的匹配。

步骤s1在数据处理阶段，将不同平台点云原始数据进行滤波去噪预处理，将预处理后点云坐标系转换到以重心为坐标原点、坐标轴方向不变的局部坐标系；局部坐标原点计算公式如下：其中，[x0，y0，z0]表示局部坐标系原点，[xi，yi，zi]表示第i个点云坐标；变换后点云坐标计算公式如下：[x′，y′，z′]＝[xi，yi，zi]-[x0，y0，z0]，其中，[x′，y′，z′]表示变换后点云坐标，m表示点云数目。

s2包括如下子步骤：

s2.1.设置相关参数：迭代次数n，有效点比例阈值p，有效阈值t；

s2.2.在点云集中随机选择不共线的3个点，计算对应的平面方程，公式如下：a1x+b1y+c1z＝d1，其中，(a1，b1，c1)表达该面的法向量，d1表示原点到平面的距离；

s2.3.计算所有观测值到平面的距离，公式如下：di＝(a1xj+b1yj+c1zj+d1)；如果di≤t，则认为该点为有效点或内点，否则为无效点；

s2.4.计算平面μ1中有效点数所占总点数的比例f，若f≥p，则认为平面μ1为有效平面，否则为无效平面，完成面点云的自身平整度评价；

s2.5.由设置的迭代次数n进行多次抽样，选出有效点占比最高的平面作为最终拟合平面；

s2.6.采用特征值法对s2.5拟合得到的平面进行参数计算，得到平面参数方程如下：ax+by+cz+d＝0，其中(a，b，c，d)为平面参数；

s2.7.在不同平台点云数据上，均匀选取k≥3对互相不平行的对于平面，规定法向量的正方向，计算其平面参数(a，b，c，d)，该平面参数被唯一确定。

步骤s3包括如下子步骤：

s3.1.构建坐标转换模型：

设同一平面在不同坐标系中的表示分别为(a1，b1，c1，d1)，(a2，b2，c2，d2)，则仿射转换模型为：其中，表示旋转矩阵参数，表示平移矩阵参数；

s3.2.旋转矩阵进行线性化，将9参数转换为3参数：

用坐标转换参数替换旋转矩阵和平移矩阵，公式如下：

其中，α为绕x轴的旋转角度，β为绕y轴的旋转角度，γ为绕z轴的旋转角度；

法向量(a2，b2，c2)计算公式如下：

a2＝cosγcosβa1+(cosγsinβsinα+sinγcosα)b1+(sinγsinα-cosγsinβcosα)c1；

b2＝-sinγcosβa1+(cosγcosα-sinγsinβsinα)b1+(sinγsinβcosα+cosγsinα)c1；

c2＝sinβa1-cosβsinαb1+cosβcosαc1；

s3.3.构建间接平差模型：

分别求α、β、γ关于a2、b2和c2的偏导数：

d2不需进行线性化处理，计算公式为：d2＝txa1+tyb1+tzc1+d1；

间接平差模型公式如下：

将上式改写成最小二乘误差方程的形式如下：

v＝ax-l；

其中，v＝(va2vb2vc2vd2)′，x＝(dαdβdγdtxdtydtz)′，

为(a2，b2，c2，d2)的估计值；

s3.4.整体最小二乘求解参数：

假设观测量l中观测误差为e，系数矩阵a中误差为ea，则最小二乘误差方程可写为：(a+ea)x＝l+e，

即存在式中，e＝[eae]；

此时最优化限制约束条件为：

采用奇异值分解的方法进行求解，对增广矩阵[al]进行奇异值分解：

分解得到的m+1个特征向量组成的正交矩阵，

则参数x的估值可表示为：

s3.5.根据平差结果的改正数对待配准平面的相关性进行评价，改正数越大，则相关性越弱，平面对越差，则需对不满足要求的平面对进行重新选取。

步骤s4包括如下子步骤：

s4.1.根据计算的3个旋转参数，对不同平台的对应点云坐标系进行转换，转换模型为：

s4.2.将转换后的点云由局部坐标系变换到原始坐标系，完成基于面特征的多平台点云匹配，公式如下：[x，y，z]＝[x2，y2，z2]+[x0，y0，z0]，其中，[x，y，z]表示经匹配变换后的最终坐标。

当然，上述说明并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的技术人员在本发明的实质范围内所做出的变化、改型、添加或替换，也应属于本发明的保护范围。