CN106408581A - 一种快速的三维点云直线提取方法 - Google Patents

Abstract

本发明公开了一种快速的三维点云直线提取方法，包括以下步骤：S1、获取原始点云数据；S2、将原始点云分割成小平面集合；S3、提取带方向向量的小平面边缘点；S4、根据方向向量将边缘点聚类成直线。本发明能够快速准确的从大规模点云数据中提取出直线结构；另外本发明的算法具有较强的抗噪声能力和很强的抗孤立点能力。对于杂乱的现实场景的点云数据，本发明可以获得很高的直线检测率和很低的误检率。

Description

一种快速的三维点云直线提取方法

技术领域

本发明涉及三维点云处理领域，尤其涉及一种快速的三维点云直线提取方法。

背景技术

点云中常见的结构信息包括平面、直线、曲率线、骨架、圆柱、曲面等等。其中，直线作为人造场景中最普遍的结构之一，提供了重要的几何信息和拓扑结构。许多点云应用诸如建筑物重建、对称检测、相机标定、数据配准和定位等都需要直线提取作为其基础步骤。同时，直线的存储结构相当简单，每条直线只需要存储两个端点的坐标信息。因此使用直线来描述点云可以有效地保留点云的大部分信息，且只需要很小的存储空间。

现有的关于直线提取方面的研究主要针对于二维图像数据，相比之下，三维点云中的直线提取则较少受到关注，研究三维点云直线提取的难点主要体现在：(1)算法的效率问题；由于大部分的算法是针对小规模问题设计的，没有考虑到点云数据的海量性，算法的时间复杂度较高，因而不适用于大规模的场景。(2)算法的鲁棒性问题；由于现实场景的多样性和激光雷达设备的差异，不同场景或不同设备下采集的点云数据也有很大的差别。这些差别体现在精度、点云分布密度、拓扑结构等方面。现有的算法往往只能针对某一种类型的数据，如针对机载点云数据的算法很难直接应用于针对车载激光雷达产生的点云数据，反之，针对车载点云数据的算法同样很难应用于机载的数据。(3)直线描述能力的问题；传统的点云直线提取工作倾向于先提取点云中的平面结构，然后提取平面结构间的交线。而且，在提取平面结构时，传统的算法倾向于将点云进行下分割，即尽可能地将点云拟合成大的平面结构。这种做法丢失了点云的细节，使得提取到的直线不足以描述原始点云。

由上可知，直线结构是描述点云的一种有效手段，且可以推动许多后续的应用，但现有的直线结构提取算法仍然有巨大的改进空间。

发明内容

本发明的目的在于克服上述现有技术的不足，提供一种可快速从点云中提取直线结构的方法。

为实现上述目的，本发明采用以下技术方案：

一种快速的三维点云直线提取方法，包括以下步骤：

S1、获取原始点云数据；

S2、将原始点云分割成小平面集合；

S3、提取带方向向量的小平面边缘点；

S4、根据方向向量将边缘点聚类成直线。

其中S2将原始点云分割成小平面集合具体包括以下步骤：

S21、计算每个点的切平面以及平滑度；

S22、根据每个点的切平面和平滑度，采用区域增长的方式得到初始的小平面集合；

S23、在初始的小平面集合的基础上，通过K均值聚类的方法得到具有更好的边缘信息的小平面集合。

步骤S3中对于每个小平面使用α-shape来将其转化为多边形结构，从而提取出边缘点，小平面的边缘点定义为小平面的α-shape的顶点除去位于共面的小平面相交处的点。

步骤S4根据方向向量将边缘点聚类成直线具体如下：

首先通过区域增长来得到具有一致的边缘点的区域，然后将具有一致边缘点区域拟合成圆柱形结构，即应用最小中值二乘法先将点集拟合成直线，并将该直线作为圆柱形结构的中心线，然后取距离中心线最远的点到中心线的距离作为圆柱形结构的半径，从而得到包含这些边缘点的圆柱形结构；一旦获得了一个最优的圆柱形结构，通过验证它的NFA值来判断是否要保留该结构，如果该圆柱形结构被拒绝，则将其区域内的点都标记为未访问，这些点可以在下次的拟合中被用到；如果接受该圆柱形结构，则将该圆柱形结构的中心线作为输出线段。

对于包含n个边缘点的圆柱形结构c，圆柱形结构c的NFA定义为：

其中，N表示边缘点的个数，N²为边缘点中可能存在的圆柱形结构的数量，k为边缘点与圆柱形结构c一致的边缘点个数，p＝θ/π为边缘点的随机的方向向量与圆柱形结构c一致的概率；当NFA(c)>1，则拒绝该圆柱形结构c。

上述步骤中一致的定义为：当且仅当边缘点在圆柱形结构c的内部，且边缘点的方向向量与圆柱形结构c的中心线的夹角不超过θ时，则称边缘点与圆柱形结构c一致。

采用上述技术方案后，本发明与背景技术相比，具有如下优点：

本发明可以从大规模点云中快速地提取的直线结构，且提取的直线可以很好的描述原始点云，因此具有重要的应用价值。另外本发明的算法具有较强的抗噪声能力和很强的抗孤立点能力。对于杂乱的现实场景的点云数据，本发明可以获得很高的直线检测率和很低的误检率。

附图说明

图1为平面分割结果对比图，其中图1(a)为原始输入点云图(加入了光照渲染)；图1(b)为使用区域增长算法得到的分割结果；图1(c)为使用VCCS超像素分割的结果；图1(d)为使用本方法小平面分割的结果。

图2为边缘点提取的示意图，其中图2(a)为给定的点云和分割好的小平面结构(不同的深浅度代表不同的小平面)；图2(b)为提取出的边缘点。

图3为用于直线聚类的圆柱形结构，其中实心的点为符合“一致条件”的边缘点，空心的点为不符合“一致条件”的边缘点。

图4(a)～(f)为直线提取结果对比图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

实施例

本发明一种快速的三维点云直线提取方法，包括以下步骤：

S1、获取原始点云数据

S2、将原始点云分割成小平面集合

通过两阶段式的小平面结构提取算法先在第一阶段对点云进行粗分割，其目的在于快速地定位点云中的局部平面结构，然后在第二阶段细化粗分割的边缘，具体如下：

S21、计算每个点的切平面以及平滑度

先采用K近邻的方式来判断一个点是否与另一个点相邻，然后计算每个与其相邻的点所构成的切平面，相比于八叉树的方法，采用K近邻的方法不需要固定体素的分辨率，更适合分布密度不均匀和存在孤立点的点云数据。

对于原始点云x的每个输入点x_i，其切平面Tp(x_i)可以表示为由其中心点o_i和法向量组成的二元组，即：

三维空间内的任一点到Tp(x_i)的距离可以表示为：

记x_i的K近邻构成的集合为Nb_K(x_i)，通过求解下式可以得到最小二乘意义下的最佳拟合平面。

这里采用主成分分析法(Principal Component Analysis,简称PCA)来求解上式。具体的，定义Nb(x_i)的协方差矩阵C为：

其中

使用奇异值分解(Singular Value Decomposition)算法来计算C的三个特征值：λ₁≥λ₂≥λ₃，以及三个对应的特征向量则即为Tp(x_i)最小二乘意义下的平面法向量，如此，则得到最佳拟合平面Tp(x_i)。

然而单纯使用最小二乘法得到的平面容易受到噪声的影响，为了提高鲁棒性，采用了迭代重加权的最小二乘法来优化拟合的平面，记由PCA得到的平面为Tp⁰(x_i)，计算每个点p∈Nb(x_i)的一个权值w(p)如下：

其中，∈为用户给定的距离阈值，用于判断点是否在平面上。注意到在激光点云中∈仅与设备的精度有关，而不受点云分布密度不均匀的影响。也就是说使用同一设备采集的点云，稀疏区域和稠密区域的∈值都是一致的。上式的权重函数是一个截断函数，即忽略到平面距离超过∈的点，而对于剩下的点根据距离的大小赋予权重。解带权最小二乘方程(下式)即可得到优化后的拟合平面Tp¹(x_i)。

对平面Tp¹(x_i)重复上述过程，直到算法收敛，通常情况下算法都可以在2到3步内收敛。最终就得到了每个点的切平面Tp(x_i)。同时还可以得到每个点的平滑度s(x_i)。记最终构成切平面Tp(x_i)的协方差矩阵的三个特征值为λ₁，λ₂和λ₃，且满足λ₁≥λ₂≥λ₃。则点x_i的平滑度s(x_i)可以表示为：

S22、根据每个点的切平面和平滑度，采用区域增长的方式得到初始的小平面集合

在已知每个点的切平面和平滑度的条件下，采用迭代选取种子点，并从种子点开始进行区域增长的方式来计算小平面。对于种子点的选取采用平滑度高的点优先的原则，因为平滑度高的点更有可能扩展成平面。

我们将小平面f_i形式化的定义为一个由所属点集P_i、中心点o_i和法向量所构成的三元组，即：对每个种子点x_i，令其初始的小平面f_i的初始点集为{x_i}，中心点和法向量分别是Tp(x_i).o_i和然后，采用宽度优先的原则对f_i进行区域增长。特别地，对每个候选点x_j，若x_j满足

(1)x_j到x_i的距离小于R_seed；

(2)向量与的夹角小于θ；

(3)x_j到平面Tp(x_i)的距离小于∈；

则将x_j加入到f_i的点集f_i.P_i中。其中，R_seed是用于限制小平面大小的阈值，不同于VCCS中的R_seed，这里的R_seed给出的是小平面的上限，因此最后得到的小平面不是固定的大小，而是半径不超过R_seed、自适应大小的平面结构。

当f_i无法再扩展时，根据f_i的点集f_i.P_i使用最小二乘法拟合平面，并将更新为拟合平面的法向量。这样就得到了初始的小平面集合。

S23、在初始的小平面集合的基础上，通过K均值聚类的方法得到具有更好的边缘信息的小平面集合

在初始小平面的基础上采用局部K均值聚类来得到更好的边缘信息。局部K均值聚类迭代地将点赋值于小平面，并保证每个点到其所属的小平面的距离小于到其他小平面的距离。这里采用和VCCS一样的距离函数：

D(x,f)＝W_s*D_s+W_c*D_c+W_n*D_n

其中，D_s表示点到小平面所在平面的最短距离；不同于VCCS，本方法希望提取到的小平面可以更接近于平面，因此当D_s超过2_∈时，令D_s＝∞。这样可以保证提取到的小平面上的每一个点到小平面所在平面的距离都不超过2_∈。D_c表示颜色距离，由于只考虑点云的几何特性而不考虑颜色空间的距离，因此设置权值W_c＝0。D_n表示法向量距离，这里设置为点的切平面法向量和小平面法向量的点积。最后设置权值W_s＝1和W_n＝4。

算法依次地从每个小平面f_i的中心点o_i开始向其邻域扩展，在扩展过程中如果出现小平面f_j中的点x_j到小平面f_i的距离小于x_j到f_j的距离，则将x_j分配给小平面f_i。并将x_j标记使其用作下一轮扩展的点。当这一轮所有的点都扩展完毕时，更新每个小平面。具体做法是，对每个小平面f_i，将f_i的中心点o_i替换为当前属于f_i的点集f_i.P_i的中心点；将f_i的法向量更新为f_i.P_i中所有点的切平面的法向量的平均值。

小平面分割算法的优势是可以在不损失算法时间效率的同时保持更好的边缘信息。例如，在如图1(a)所示的原始输入点云图，经过上述方法处理，最后得到的小平面分割的结果如图1(d)所示，本发明的算法得到的是自适应大小的小平面结构，很好的保留了局部平面结构和边缘信息，如图中建筑物窗户的边界。为了证明算法的优越性，还给出了如图1(b)所示的传统的基于区域增长算法得到的(大)平面结构图，以及如图1(c)所示的VCCS超像素分割的结果图，传统的区域增长的算法倾向于将点云拟合成大平面，导致了细节上的丢失；而VCCS虽然得到了半径相近的分割结果，但是却不能很好的保留边缘信息。

S3、提取带方向向量的小平面边缘点

将点云分割成小平面之后，可以只考虑位于这些小平面边缘的点，对于每个小平面，可以使用α-shape来将其转化为多边形结构，从而提取出它的边缘点。但是对于两个小平面共面的情况，两个小平面之间的交点应该被忽略。如图2(a)所示为分割好的三个不同的小平面结构，图2(b)所示的框内的点即为两个共面的小平面之间的交点，边缘点为图2(b)中的深色点，定义边缘点为小平面的α-shape的顶点除去位于共面的小平面相交处的点。

在具体实现中，对每个小平面以及其相邻的小平面集合NF_i，定义F_i为{f_i}∪NF_i的一个子集，且满足：F_i中的每个小平面与小平面共面。令V_i为F_i的α-shape的顶点，则小平面f_i的边缘点即为V_i∩F_i。

另外，对每个边缘点p，定义边缘点p的方向向量为V_i中与p相邻的点在最小二乘意义下的主方向。这些方向向量在直线聚类中起到了重要的作用。

S4、根据方向向量将边缘点聚类成直线

直线聚类是本发明直线提取算法中的关键步骤，直线聚类算法的主要思想是对每个候选的直线，动态地维护一个圆柱形结构(如图3所示)。每个圆柱形结构可以由中心线l和半径r确定，是三维点云空间的一个几何区域。这些圆柱形结构最终都可以使用它们的中心线l来近似，可以通过这些圆柱形结构来滤除假阳性(false positive)的检测，这有助于降低杂乱场景中直线提取的错误率。

为了给出具体的算法描述，首先定义关于边缘点和圆柱形结构的“一致”的概念：当且仅当边缘点p在圆柱形结构c的内部，且边缘点p的方向向量与圆柱形结构c的中心线的夹角不超过θ时，则称边缘点p与圆柱形结构c一致。

“对齐”的概念为判断一个圆柱形结构能否近似成直线提供了一个直观方法。朴素的做法包括：验证圆柱形结构内对齐点的个数与总点数的比例是否超过某一个阈值。然而，考虑到来自不同设备的不同类型的点云数据，使用这种固定阈值的方法显然是不合适的。为了解决这一问题，本方法引入了NFA(Number of False Alarm)的概念。NFA是由Desolneux等人提出来的用于解决二维图像中像素对齐问题，并被进一步应用于二维直线检测中的概念。这里将NFA的概念扩展到三维，并将其作为验证圆柱形结构能否被近似成直线结构的工具。

用B来表示提取出的边缘点集合，假设另一个与B具有相同大小的边缘点集合B₀，B₀中的每个边缘点的方向向量是随机设置的，称B₀为一个随机模型。根据Helmholtz原理，希望从随机模型B₀中提取出的直线应该足够的少。若某种结构出现在B₀中是小概率事件，则根据反证法，可以在B中接受这一结构。特别地，可以估计圆柱形结构c在随机模型B₀中出现的期望，称这样的期望为NFA(Number of False Alarm)，其形式化定义如下：

定义NFA：对于给定包含n个边缘点的圆柱形结构c，其中与c一致的边缘点个数为k，圆柱形结构c的NFA定义为：

其中，N表示边缘点的个数，N²为边缘点中可能存在的圆柱形结构的数量，由于每两个边缘点就可以确定一个圆柱形结构，因此最多可能存在N²个圆柱形结构，p＝θ/π为边缘点的随机的方向向量与圆柱形结构一致的概率。

若圆柱形结构c的NFA足够小，即：NFA(c)≤ε，则称c为有意义的结构，将ε固定为1，即对于随机模型B₀，最多只允许出现一个被检测为有意义的圆柱形结构，这样就可以简单地通过拒绝掉NFA(c)>1的圆柱形结构来过滤假阳性检测。

在上述的圆柱形结构和NFA的基础上通过对边缘点进行聚类来提取直线。首先通过区域增长来得到具有一致的边缘点的区域，然后通过拟合的方法计算包含这些边缘点的圆柱形结构，在这里采用了最小中值二乘法(LMS)将具有一致的边缘点区域拟合成圆柱形结构，即应用最小中值二乘法(LMS)先将点集拟合成直线，并将该直线作为圆柱形结构的中心线；然后取距离中心线最远的点到中心线的距离作为圆柱形结构的半径。

一旦获得了一个最优的圆柱形结构c，通过验证它的NFA值来判断是否要保留该结构。如果圆柱形结构c被拒绝，则将其区域内的点R都标记为未访问的，这样这些点可以在下次的拟合中被用到。如果接受c，就将该圆柱形结构c的中心线作为输出线段。

本发明的最终处理结果如图4所示，其中，图4(a)为几何图形，图4(b)、图(d)和图4(e)为建筑图形，图4(c)为工业零件，图4(f)为道路风景图形。

以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。

Claims (2)

1.一种快速的三维点云直线提取方法，其特征在于，包括以下步骤：

S1、获取原始点云数据；

S2、将原始点云分割成小平面集合，具体包括以下步骤：

S21、先采用K近邻的方式来判断一个点是否与另一个点相邻，然后计算每个点与其相邻点所构成的切平面以及平滑度；

S22、根据每个点的切平面和平滑度，采用区域增长的方式得到初始的小平面集合，具体为：

在已知每个点的切平面和平滑度的条件下，依次选取种子点，并从种子点开始进行区域增长的方式来计算小平面，对于种子点的选取采用平滑度高的点优先的原则；

当小平面的区域无法再扩展时，使用最小二乘法拟合平面，即得到初始的小平面集合；

S23、在初始的小平面集合的基础上，通过K均值聚类的方法得到具有更好的边缘信息的小平面集合，具体为：

在初始小平面的基础上采用局部K均值聚类来得到更好的边缘信息，局部K均值聚类迭代地将点赋值于小平面，并保证每个点到其所属的小平面的距离小于到其他小平面的距离。这里采用和VCCS一样的距离函数：

D(x,f)＝W_s*D_s+W_c*D_c+W_n*D_n

其中，D_s表示点到小平面所在平面的最短距离；D_c表示颜色距离，设置权值W_c＝0。D_n表示法向量距离，这里设置为点的切平面法向量和小平面法向量的点积。最后设置权值W_s＝1和W_n＝4。

S3、提取带方向向量的小平面边缘点；

S4、根据方向向量将边缘点聚类成直线，具体步骤如下：

首先通过区域增长来得到具有一致的边缘点的区域，然后将符合“一致条件”的边缘点区域拟合成圆柱形结构，即应用最小中值二乘法先将点集拟合成直线，并将该直线作为圆柱形结构的中心线，然后取距离中心线最远的点到中心线的距离作为圆柱形结构的半径，从而得到包含这些边缘点的圆柱形结构；

对于包含n个边缘点的圆柱形结构c，圆柱形结构c的NFA定义为：

$NFA\left(c\right)=N^2\·\underset{i=k}{\overset{n}{\Σ}}\left(\begin{array}{c}i\\ n\end{array}\right)p^i{(1-p)}^{n-i}$

其中，N表示边缘点的个数，N²为边缘点中可能存在的圆柱形结构的数量，k为边缘点与圆柱形结构c一致的边缘点个数，p＝θ/π为边缘点的随机的方向向量与圆柱形结构c一致的概率；当NFA(c)>1，则拒绝该圆柱形结构c；

所述一致的定义为：当且仅当边缘点在圆柱形结构c的内部，且边缘点的方向向量与圆柱形结构c的中心线的夹角不超过θ(θ＝22.5。)时，则称边缘点与圆柱形结构c一致；

一旦获得了一个最优的圆柱形结构，通过验证它的NFA值来判断是否要保留该结构，如果该圆柱形结构被拒绝，则将其区域内的点都标记为未访问，这些点可以在下次的拟合中被用到；如果接受该圆柱形结构，

则将该圆柱形结构的中心线作为输出线段。

2.根据权利要求1所述的一种快速的三维点云直线提取方法，其特征在于：步骤S3中对于每个小平面使用α-shape来将其转化为多边形结构，从而提取出边缘点，小平面的边缘点定义为小平面的α-shape的顶点除去位于共面的小平面相交处的点。