本发明涉及一种激光点云数据重建实物三维模型的方法，尤其是涉及一种基于动态网格k邻域搜索的激光点云精简方法。
背景技术
现有的三维激光扫描仪获取到的点云数据含有丰富的细节特征，数据量庞大，其中包含大量的冗余数据点，若不进行必要的数据精简预处理，则会严重影响模型重建过程的效率，过密的数据点也会影响目标曲面重构的光滑性，甚至会由于存在噪声点而无法实现模型重建。因此在保留模型特征的前提下对点云数据进行精简处理是十分重要且具有实际意义的工作。
近年来，对散乱点云的精简方法主要分为两类：基于三角网格模型的点云精简方法和直接基于数据点的点云精简方法。
基于三角网格模型的点云精简方法需要先对点云数据进行三角剖分处理，建立其相应的三角网格拓扑结构，然后再对该三角网格进行处理，从而达到点云数据精简的目的。该方法过程比较复杂，需要消耗大量的计算机系统资源，并且该方法的抗噪能力较弱，对含有噪声的点云数据，构造的三角网格可能会出现变形等情况，精简后的点云模型与原始模型可能大不相同。
直接基于点云数据的精简方法依据点云数据点之间的空间位置关系来建立点云的拓扑连接关系，并根据建立的拓扑连接关系计算点云数据中每个数据点的几何特征信息，最后根据这些特征信息来对点云数据进行点云精简处理。相比基于三角网格的点云精简方法，直接基于点云数据点的精简方法由于无需计算和存储复杂的三角网格结构，使得其精简的效率相对较高，同时避免了三角网格可能会出现变形的情况，成为点云精简方法的主流。
直接基于点云数据的精简方法中，建立空间点云数据中点与点之间的拓扑连接关系，需要构建数据点的k邻域。建立k邻域最直接简单的方法是查找点云数据中指定点除外的每个点并计算欧氏距离，取距离最近的k个点，显然这是一种效率较低的方法。目前散乱点云的k近邻搜索常用的方法有空间栅格法、八叉树法和k-dtree法。空间栅格法算法简单运行速度快，但不适用于非均匀点云数据。八叉树法一般采用递归的数据结构，将非叶子节点再次划分为八个子节点，大量消耗了计算机的内存空间。而改进的自适应性八叉树虽然通过一次精简减少了存储消耗，但依然存在相邻节点之间需要编号使得过程繁琐等问题。
技术实现要素：
本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种基于动态网格k邻域搜索的激光点云精简方法。
本发明的目的可以通过以下技术方案来实现：
一种基于动态网格k邻域搜索的激光点云精简方法，该方法包括以下步骤：
(一)构建动态约束网格，获取点云数据点的k邻域点；
101)扩展一定边长的立方体范围，并能根据实际情况改变边长的大小，以此获取一个动态约束网格；
102)在扩展的立方体范围内，求出此范围内点集的点数m；
103)当点数m的值为αk≤m≤βk时，在该点集内搜索k邻域点，当m≤αk时，将边长的取值增大，扩大网格的约束范围，回到步骤102)，当m≥βk时，将边长的取值减小，缩小网格的约束范围，回到步骤102)，其中，α和β为调节系数。
(二)利用pca方法对点云数据点的k邻域点进行协方差分析，估算点云数据点的法向量，在法向量估算的基础上进行点云数据点的曲率估算。
(三)计算点云数据点与邻域点法向夹角平均值及平均距离；
点与邻域点法向夹角平均值θ(gi)的计算公式为：
式中，gi与gj分别为数据点及其邻域点，
式中，
点与邻域点的平均距离d(gi)的计算公式为：
(四)根据获取的曲率、法向夹角平均值及平均距离，定义特征判别参数及判别阈值：
特征判别参数w(gi)定义为：
式中，λh，λθ，λd分别为曲率系数、夹角系数和距离系数，hi为点云数据点的曲率；
特征判别阈值δ定义为：
式中，η为特征点数量控制系数，n为点云数据点个数。
判断特征判别参数是否大于判别阈值，若是，将该点云数据点作为特征点保留，输出点云数据，否则，将该点云数据点构成非特征点集，执行下一步。
优选地，曲率系数λh选取200，邻域点个数k选取15。
(五)采用包围盒精简法对非特征点进行二次精简，获取最终点云数据；
501)将判定的非特征点作为待处理点云集读入，求出点云数据点在x、y、z方向上的最大值和最小值，获取大包围盒的边长；
502)根据精简率的需要确定小包围盒的边长，将大包围盒分割成均匀大小的小包围盒；
503)将点云集中所有点按其三维坐标划分到不同的小包围盒内，对每一个小包围盒，获取其内部所有点到其中心的距离并进行排序，保留距离中心最近的点；
504)将步骤503)中保留的点与步骤4)判定的特征点数据进行拼接整合，获取二次精简后的点云数据。
(六)采用最大误差和平均误差对精简后的点云数据与原始点云进行评估，最大误差的表达式为：
平均误差的表达式为：
式中，d(g，s*)为原始曲面s上采样点g到二次精简点云曲面s*上投影点g*的欧式距离，假设采样点g的法向量为np，则d(g，s*)＝np*(g*-g)。
与现有技术相比，本发明具有以下优点：
(1)本发明方法以点集中某点为中心构建可扩展可收缩的动态网格，通过在动态网格内获取每一点到中心点的距离来获取点云数据点的k邻域点，与传统的包围盒搜索k邻域方法相比，避免了在划分空间网格时，由于网格内点数的不确定性而存在的繁琐的编号程序，本发明方法可根据点的密度，随意扩大或缩小网格，从而快速求得k邻域点，适用于非均匀点云数据和均匀点云数据；
(2)本发明方法对点云数据进行了两次精简，采用了基于包围盒精简法对点云数据的非特征点进行二次精简，将精简后的点云与特征点拼接，最终得到精简后的点云数据，能高精度的保留模型几何特征，有效弥补了一次精简造成的空白区域出现的问题；
(3)本发明方法采用基于回归方法的pca法矢估计方法来计算点云数据点的曲率，通过构造协方差矩阵估计曲面变化及法向量，算法简单且速度快，可提高方法的效率。
附图说明
图1为本发明方法的流程图；
图2为数据模型中特征点与曲率关系图；
图3为本发明实施例中初步特征点的提取图像，其中，图3(a)为兔子图像的原始点云，图3(b)为兔子图像的特征点，图3(c)为椅子图像的原始点云，3(d)为椅子图像的特征点；
图4为本发明实施例中采用二次精简的兔子图像的处理结果；
图5为本发明实施例中采用包围盒法的精简效果图，其中，图5(a)、图5(b)、图5(c)分别为精简率为15％、68.5％、91.3％的兔子图像的精简效果图，图5(d)、图5(e)、图5(f)分别为精简率为19.51％、66.08％、90.94％的椅子图像的精简效果图；
图6为本发明实施例中采用随机采样法的精简效果图，其中，图6(a)、图6(b)、图6(c)分别为精简率为15％、68％、91％的兔子图像的精简效果图，图6(d)、图6(e)、图6(f)分别为精简率为19％、66％、91％的椅子图像的精简效果图；
图7为本发明实施例中采用本发明二次精简方法的精简效果图，其中，图7(a)、图7(b)、图7(c)分别为精简率为14.34％、69.23％、92.17％的兔子图像的精简效果图，图7(d)、图7(e)、图7(f)分别为精简率为18.65％、64.89％、90.91％的椅子图像的精简效果图。
具体实施方式
下面结合附图和具体实施例对本发明进行详细说明。
实施例
如图1所示，本发明涉及一种基于动态网格k邻域搜索的激光点云精简方法，该方法应用于激光点云数据重建实物三维模型，该方法包括以下步骤：
s1：构建动态约束网格，获取点云数据点的k邻域点：
101)以点集中的某一数据点为中心，向四周沿x、y、z轴正反方向扩展l/2的距离，形成边长为l的立方体；
102)在扩展的立方体范围内，求出此范围内的点数m；
103)如果m≥αk，转到步骤104)，否则，转到步骤105)；
104)如果m≤βk，转到步骤106)，否则，l＝l-δl，缩小立方体的范围，转步骤102)，其中，α和β为调节系数，当测量点分布均匀时，α和β的值取小，否则值取大；
105)令l＝l+δl，扩大立方体的范围，转到步骤102)；
106)求出范围内每一点到中心点的距离，并按距离升序排列，得到数据点的k邻域点；
107)结束。
s2：利用pca方法通过对点云数据点的k个邻域点进行协方差分析，估算点云数据点的法向量，在法向量估算的基础上进行点云数据点的曲率估算；
如图2所示，实心点是模型点，水平线为切平面，空心点表示模型点的邻域点。
g1点位于曲率较大的曲线上，g1点的邻域点与切平面间的距离相对较远，g1点为特征点；g2点位于曲率较小的曲线上时，g2点的邻域点与切平面的距离较近，此点为非特征点，g2点为非特征点。
采用pca方法估计曲率及法向量：设点云数据集为g＝{gi(xi，yi，zi)|i＝1，2，…，n}，点gi的邻域点集为{gij(xij，yij，zij)j＝1，2，…，k}，邻域重心可表示为：
则点gi的协方差矩阵ti定义为：
协方差矩阵ti是一个半正定的对称矩阵，定义了局部曲面的几何信息。使用jacobi方法求解协方差矩阵ti，解出矩阵ti的三个特征值λ1、λ2、λ3以及特征值对应的特征向量n1、n2、n3。若λ3≥λ2≥λ1，则最小特征值λ1对应的特征向量n1就是数据点gi局部曲面的法向量，λ1描述了曲面沿法向量的变化，而λ2和λ3表示数据点在切平面上的分布情况。因此数据点gi在k邻域内的曲面变分定义为：
点云模型在数据点gi的曲率hi可近似为在该点的曲面变分τi，即hi≈τi。
s3：计算点与邻域点法向夹角平均值及平均距离；
点与邻域点法向夹角平均值计算：
设数据点gi是点云模型g中的任意点，gj是gi的邻域点，它们的法向分别为
其中，
对数据点与所有邻域点的法向夹角取平均值得到法向夹角参数：
法向夹角参数综合考虑了数据点gi的所有邻域点对数据点弯曲程度的影响，若θ(gi)值越大，数据点gi及邻域的曲面弯曲程度越大，数据点gi邻域所在区域是特征区域的可能性越大；θ(gi)值越小，数据点gi邻域的曲面弯曲程度越小，模型表面越平滑，数据点gi邻域所在区域是特征区域的可能性越小。
点与邻域点平均距离计算：
当数据点附近点云分布比较密集时，点与周围邻域点的平均距离较小，该区域为模型特征区域的可能性较大；反之，当数据点附近点云分布比较稀疏时，点与周围邻域点的平均距离较大，该区域是模型平坦区域的可能性较大。
对任一点gi与邻域点的平均距离为：
s4：根据数据点的曲率及步骤s3获取的法向夹角平均值和平均距离两个参数，定义特征判别参数及判别阈值；
由于在非均匀点云模型中特征点所在区域曲面变化较明显、点云密度较高，曲率hi和点与邻域点的法向夹角参数θ(gi)越大、点与邻域点的平均距离d(gi)越小，数据点gi为特征点的可能性越大。
因此曲率hi和点与邻域点的法向夹角参数θ(gi)与特征判别参数w(gi)成正比，平均距离d(gi)与特征判别参数w(gi)成反比。
定义特征判别参数w(gi)为：
其中，λh，λθ，λd分别为调节系数，λh为曲率系数，λθ为夹角系数，λd为距离系数。
经过对不同数据实验结果的分析，曲率系数λh对计算结果影响较大。邻域点个数的取值取决于点云数据模型的密度和分布的均匀性，当点云模型密度大，并且数据密集时该值可取小些，一般该值取10-30之间是一个较好的选择。当点云模型中包含噪声或数据分布过于密集时，需要对距离系数λd进行调整。
为了避免特征判别参数中的调节系数会根据模型不同而差异过大难以设置，本方法中，定义特征判别阈值δ为：
其中，η为特征点数量控制系数，n为点云数据点个数。
由此定义特征点的标识：
若某个数据点的特征判别参数w(gi)大于阈值δ，ω(gi)＝1，那么这个点gi被判定为特征点；若特征判别参数w(gi)小于阈值δ，ω(gi)＝0，gi为非特征点。
s5：采用包围盒精简法对非特征点进行二次精简，获取最终点云数据；
501)将步骤s4标识为0的非特征点作为待处理点云集读入，求出点云数据点在x、y、z方向上的最大值和最小值，得到大包围盒的边长；
502)根据精简率的需要确定小包围盒的边长，把大包围盒分割成均匀大小的小包围盒；
503)将点云集中所有点按其三维坐标划分到不同的小包围盒内，对每一个小包围盒，求出其内部所有点到其中心的距离并进行排序，只保留距离中心最近的点；
504)将上一步骤中得到的点与步骤s4标识为1的特征点数据进行拼接整合，得到精简后的最终点云数据。
s6：对精简后的点云数据进行精度评估。
为了评估精简点云的精度，需要计算原始点云与精简点云之间的几何误差，本发明求取精简点云与原始点云的最大误差和平均误差。
最大误差：
平均误差：
式中，d(g，s*)表示原始曲面s上采样点g到精简点云曲面s*上投影点g*的欧式距离；假设采样点g的法向矢量为np，则d(g，s*)＝np*(g*-g)。
本发明中邻域点个数k和曲率系数λh都为固定值，不需要根据模型不同进行修改。距离系数λd需要修改，但这个参数也只是对提取的特征点个数有影响。距离系数的改变不会导致算法不稳定，当λd比较大时，会得到特征区域中比较多的特征点，如果λd设置比较小，检测出的特征点就是比较尖锐的特征。本发明中优选曲率系数λh＝200，邻域点个数取k＝15。
为证明本发明方法的有效性，本实施例对两种激光点云数据重建实物进行了处理，一种为兔子图像，一种为椅子图像，并同时采用包围盒法、随机采样法与本发明方法的效果进行了对比，如图3～图7所示。三种方法在精简率相差不多的情况下，随着精简率的提高，包围盒法和随机采样法的精简后模型的细微特征变模糊，特别是在精简到68％左右时，有孔洞出现，且模型的细节特征丢失严重。本发明方法的精简效果图清晰完整，无细节特征丢失的情况，在精简率为69.23％和90.91％时也能较好地保持模型细节特征且避免孔洞区域的出现。
以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的工作人员在本发明揭露的技术范围内，可轻易想到各种等效的修改或替换，这些修改或替换都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以权利要求的保护范围为准。