CN106600617A - 基于曲率从Lidar点云数据提取建筑物轮廓线的方法 - Google Patents

Abstract

本发明公开了一种基于曲率从Lidar点云数据提取建筑物轮廓线的方法，涉及Lidar数据处理技术领域，应用于数字城市建设过程中的轮廓线提取及三维模型建立的生产过程中。该方法包括步骤：点云数据去噪处理；点云数据的动态网格划分；点云法矢和曲率的求解；点云数据的区域分割；点云数据的轮廓线提取。所以，该方法在地理信息行业数字城市轮廓线提取及三维建模过程中，为自动计算提取规则建筑轮廓线结构提供了一种可靠的解决方案，以改善现阶段纯人工提取建筑轮廓线的工作方式，提高数字城市三维数据生产的效率。

Description

基于曲率从Lidar点云数据提取建筑物轮廓线的方法

技术领域

本发明涉及Lidar数据处理技术领域，尤其涉及一种基于曲率从Lidar点云数据提取建筑物轮廓线的方法，应用于数字城市建设过程中的轮廓线提取及三维模型建立的生产过程中。

背景技术

建筑物作为城市区域重要的地表特征之一，在现代城市规划与建设，灾害预防，自然资源勘查与利用，数字城市建设等很多领域都有着重要的应用。

近几年来，随着城镇数字化进程的加快，基于快速且准确的数据源的城市建筑物研究逐渐成为一种新的发展趋势。

在数字城市三维模型的重建和大型工程施工中的质量监测等测绘领域，常常需要提取建筑物的轮廓线，目前，在建筑物轮廓线的提取工作中，大多仍然采用纯人工的方式，这种方式不仅费时费力，需要大量的人力消耗，而且工作效率低下，提取得到的轮廓线的精度差。

Lidar(LIght Detection And Ranging)技术使数据获取和处理向自动化的方向发展，研究如何利用Lidar数据自动提取建筑物的轮廓线的方法在数字城市三维建模过程中将具有十分重要的现实意义。

发明内容

本发明的目的在于提供一种基于曲率从Lidar点云数据提取建筑物轮廓线的方法，从而解决现有技术中存在的前述问题。

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

一种基于曲率从Lidar点云数据提取建筑物轮廓线的方法，包括如下步骤：

S1，对散乱分布的点云数据进行去噪处理，得到规则的点云数据；

S2，对所述规则的点云数据，通过设定容差e和阈值t，采用动态层法和局部重划法相结合的点云数据划分方法，对点云特征突变区部分网格进行动态调整，得到动态网格；

S3，将点云数据经过坐标转化法转换后，在每个网格内进行局部曲面拟合，求解点云法矢和曲率；

S4，以所述动态网格为基础，利用S3得到的法矢和曲率，根据曲率、法矢与曲面的关系，采用K-means方法，对点云数据进行区域分割，得到分割的点云数据块；

S5，采用八叉树的方法对所述点云数据块的边界点进行提取，利用DirectX绘制线性结构的函数，绘制边界线，完成点云数据的轮廓线提取。

优选地，S2包括如下步骤：

S21，定义点(Xmax+e，Ymax+e，Zmax+e)和点(Xmin-e，Ymin-e，Zmin-e)为对角点,表面平行于坐标平面的空间六面体为点云空间包围盒，其中，Xmax，Ymax,Zmax,Xmin,Ymin,Zmin分别为点云数据中最大点与最小点的x、y、z坐标值，e为容差值；

S22，沿坐标轴X,Y,Z方向，按等间隔t，将所述点云空间包围盒划分成多个空间六面体网格则在所述点云空间包围盒内，沿X,Y,Z方向，以t为步长分段，得到在三个方向上划分的段数m,n,l，用如下公式表示为：

，

则，根据上式，可知

对于空间内任意一个点P，都会属于某一个位置的基本正方体(i,j,k)内，坐标(x_p,y_p,z_p)满足：

S23，统计所述空间六面体网格内点云的个数并计算点云密度，并判断所述空间六面体网格是否为空，如果是，则删除。

S24，计算每一个网格的中心点O_i，并在网格的数据集中找到距离中心点O_i最近的点P_i；

S25，以所述P_i为基础，对网格内的数据分别在XYZ方向进行排序，根据数据的分布状况判断是否需要对网格进行调整，如果需要，则进行调整，否则，无需处理。

优选地，S25中，所述根据数据的分布状况判断是否需要对网格进行调整，具体为：

a.如果网格内的数据分布密度大于设定的密度阈值，点云个数超过了抛物面拟合的需要，点云数据之间在X,Y,Z方向的角度偏差小于设定的角度偏差阈值，则无需对网格进行调整；

b.如果网格内的数据分布密度大于设定的密度阈值，点云个数超过了抛物面拟合的需要，但是点云数据之间在X,Y，Z方向的角度偏差超出了角度偏差阈值，则根据点云偏差的突变，将原网格在角度偏差突变的方向上进行细分， 添加网格；

c.如果网格内数据分布均匀，点云个数满足抛物面拟合的需要，并且点云数据之间在X,Y,Z方向的角度偏差在角度偏差阈值范围内，则无需对网格进行调整；

d.如果网格内的数据分布密度小于设定的密度阈值，点云个数不能满足抛物面拟合的需要，但是点云数据之间在X,Y,Z方向的角度偏差在角度偏差阈值范围内，则以原网格为中心，边长2t扩大网格，或沿某一方向将网格与邻近网格合并，得到新的网格，并对所述新的网格进行重新判断，直到满足上述c中的情况。

优选地，S3具体为：每一个网格内，求出距离网格中心点最近的点集中的点Pi，利用Pi的k近邻拟合抛物面S，通过计算S的曲率特性来确定该数据点的曲率特性；

按照如下方法进行实施：每一个网格内点Pi的k近邻，记为Nd(Pi)拟合平面估算测量点法矢，目标函数为：可得到如下对称的半正定3阶矩阵:

其中:Cv是最小特征值对应的特征向量单位化，即为曲面在点Pi的法矢Ni，多次实验表明，当k取为24～30，曲率计算能获得较好效果；

在点Pi处建立局部坐标系(u,v,h)，且坐标原点为Pi，坐标轴h设为法矢Ni,坐标轴u和v可以在过点Pi垂直于Ni的平面内任意确定，将点Pi的k近邻由全局坐标(xi,yi,zi)转化为局部坐标系(ui,vi,hi)，并用转换后数据拟合如下面的抛物面：

S(u,v)＝a₀+a₁u+a₂v+a₃uv+a₄u²+a₅v²，

抛物面拟合通过求解目标函数实现，应用奇异值分解法可获得拟合抛物面的最小二乘解，上述方法包括如下步骤：

S31，将每一个划分好的网格邻域内点的坐标代入抛物面S方程中,设每个网格内有k个点，则得到如下方程组：

AX＝Z，

其中

X＝{a₀,a₁,a₂,a₃,a₄,a₅}^T，Z＝{h₀,h₁,.....h_k}^T；

S32，采用奇异值分解法求解S31中的线性方程组,由X＝A'Z求解曲面方程的系数a,b,c；

S33，求出系数a,b,c后，对拟合抛物面S(u,v,h)＝S(u,v,au²+buv+cv²)求一、二阶偏导数，得到Pi点的平均曲率H和高斯曲率K：

Su＝(1,0,2au+bv)，Suu＝(0,0,2a)，Sv＝(0,1,bu+2cv)，Svv＝(0,0,2c)，Suv＝(0,0,b)，

Pi点的法矢为：

由抛物面的第一、第二基本形式，可以得到：

E＝Su·Su，F＝Su·Sv，G＝Sv·Sv，L＝Suu·n，M＝Suv·n，N＝Svv·n，

则得到Pi点的平均曲率(H)和高斯曲率(K):

优选地，S5包括如下步骤：

S51，对于分割后的每块点云数据，在读取数据的过程中，得到待提取边界点云数据的最大与最小x、y、z坐标值，根据所述坐标值获得点云数据的最小外包围盒；

S52，判断S51中的每一个所述最小外包围盒内点云的个数，如果所述最小外包围盒内点云个数大于1，则对其进行八叉树网格化分，直到划分后得到的网格内点云个数不大于1，统计所述得到的网格和其相邻的空间中的26个网格是空网格的个数，得到统计结果；

S53，根据得到的统计结果判断所述划分后得到的网格是否为边界网格，如果该网格为边界网格，则返回S52，如果不是边界网格，则不用对其进行处理，并将点的内存释放；

S54，提取所述边界网格中的点云数据，得到边界点云，对所有的所述边界点云，根据网格的空间关系，进行索引排序，得到排序后的点云数据以及点云数据的个数，利用DirectX绘制线性结构的函数，绘制边界线，从而得到点云数据的轮廓线特征。

本发明的有益效果是：本发明实施例提供的基于曲率从Lidar点云数据提取建筑物轮廓线的方法，包括步骤：点云数据去噪处理；点云数据的动态网格划分；点云法矢和曲率的求解；点云数据的区域分割；点云数据的轮廓线提取。在地理信息行业数字城市轮廓线提取及三维建模过程中，为自动计算提取规则建筑轮廓线结构提供了一种可靠的解决方案，以改善现阶段纯人工提取建筑轮廓线的工作方式，提高数字城市三维数据生产的效率。

附图说明

图1是本发明利用曲率进行点云数据的分割和轮廓线提取的方法流程图；

图2是本发明采用八叉树的方法进行分割后点云数据的边界点提取方法示 意图；

图3是本发明采用八叉树的划分方法。

具体实施方式

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

如图1所示，本发明实施例提供了一种基于曲率从Lidar点云数据提取建筑物轮廓线的方法，包括如下步骤：

S1，对散乱分布的点云数据进行去噪处理，得到规则的点云数据；

S2，对所述规则的点云数据，通过设定容差e和阈值t，采用动态层法和局部重划法相结合的点云数据划分方法，对点云特征突变区部分网格进行动态调整，得到动态网格；

S3，将点云数据经过坐标转化法转换后，在每个网格内进行局部曲面拟合，求解点云法矢和曲率；

S4，以所述动态网格为基础，利用S3得到的法矢和曲率，根据曲率、法矢与曲面的关系，采用K-means方法，对点云数据进行区域分割，得到分割的点云数据块；

S5，采用八叉树的方法对所述点云数据块的边界点进行提取，利用DirectX绘制线性结构的函数，绘制边界线，完成点云数据的轮廓线提取。

其中，S2包括如下步骤：

S21，定义点(Xmax+e，Ymax+e，Zmax+e)和点(Xmin-e，Ymin-e，Zmin-e)为对角点,表面平行于坐标平面的空间六面体为点云空间包围盒，其中，Xmax，Ymax,Zmax,Xmin,Ymin,Zmin分别为点云数据中最大点与最小点的x、y、z坐标值，e为容差值；

S22，沿坐标轴X,Y,Z方向，按等间隔t，将所述点云空间包围盒划分成多个空间六面体网格则在所述点云空间包围盒内，沿X,Y,Z方向，以t为步长分段，得到在三个方向上划分的段数m,n,l，用如下公式表示为：

则，根据上式，可知

对于空间内任意一个点P，都会属于某一个位置的基本正方体(i,j,k)内，坐标(x_p,y_p,z_p)满足：

S23，统计所述空间六面体网格内点云的个数并计算点云密度，并判断所述空间六面体网格是否为空，如果是，则删除。

S24，计算每一个网格的中心点O_i，并在网格的数据集中找到距离中心点O_i最近的点P_i；

S25，以所述P_i为基础，对网格内的数据分别在XYZ方向进行排序，根据数据的分布状况判断是否需要对网格进行调整，如果需要，则进行调整，否则，无需处理。

S25中，所述根据数据的分布状况判断是否需要对网格进行调整，具体为：

a.如果网格内的数据分布密度大于设定的密度阈值，点云个数超过了抛物面拟合的需要，点云数据之间在X,Y,Z方向的角度偏差小于设定的角度偏差阈值，则无需对网格进行调整；

b.如果网格内的数据分布密度大于设定的密度阈值，点云个数超过了抛物面拟合的需要，但是点云数据之间在X,Y，Z方向的角度偏差超出了角度偏差阈值，则根据点云偏差的突变，将原网格在角度偏差突变的方向上进行细分，添加网格；

c.如果网格内数据分布均匀，点云个数满足抛物面拟合的需要，并且点云数据之间在X,Y,Z方向的角度偏差在角度偏差阈值范围内，则无需对网格进行调整；

d.如果网格内的数据分布密度小于设定的密度阈值，点云个数不能满足抛物面拟合的需要，但是点云数据之间在X,Y,Z方向的角度偏差在角度偏差阈值范围内，则以原网格为中心，边长2t扩大网格，或沿某一方向将网格与邻近网格合并，得到新的网格，并对所述新的网格进行重新判断，直到满足上述c中的情况。

S3具体为：每一个网格内，求出距离网格中心点最近的点集中的点Pi，利用Pi的k近邻拟合抛物面S，通过计算S的曲率特性来确定该数据点的曲率特性；

按照如下方法进行实施：每一个网格内点Pi的k近邻，记为Nd(Pi)拟合平面估算测量点法矢，目标函数为：可得到如下对称的半正定3阶矩阵:

其中:Cv是最小特征值对应的特征向量单位化，即为曲面在点Pi的法矢Ni，多次实验表明，当k取为24～30，曲率计算能获得较好效果；

在点Pi处建立局部坐标系(u,v,h)，且坐标原点为Pi，坐标轴h设为法矢Ni,坐标轴u和v可以在过点Pi垂直于Ni的平面内任意确定，将点Pi的k近邻由全局坐标(xi,yi,zi)转化为局部坐标系(ui,vi,hi)，并用转换后数据拟合如下面的抛物面：

S(u,v)＝a₀+a₁u+a₂v+a₃uv+a₄u²+a₅v²，

抛物面拟合通过求解目标函数实现，应用奇异值分解法可获得拟合抛物面的最小二乘解，上述方法包括如下步骤：

S31，将每一个划分好的网格邻域内点的坐标代入抛物面S方程中,设每个网格内有k个点，则得到如下方程组：

AX＝Z，

其中

X＝{a₀,a₁,a₂,a₃,a₄,a₅}^T，Z＝{h₀,h₁,.....h_k}^T；

S32，采用奇异值分解法求解S31中的线性方程组,由X＝A'Z求解曲面方程的系数a,b,c；

S33，求出系数a,b,c后，对拟合抛物面S(u,v,h)＝S(u,v,au²+buv+cv²)求一、二阶偏导数，得到Pi点的平均曲率H和高斯曲率K：

Su＝(1,0,2au+bv)，Suu＝(0,0,2a)，Sv＝(0,1,bu+2cv)，Svv＝(0,0,2c)，Suv＝(0,0,b)，

Pi点的法矢为：

由抛物面的第一、第二基本形式，可以得到：

E＝Su·Su，F＝Su·Sv，G＝Sv·Sv，L＝Suu·n，M＝Suv·n，N＝Svv·n，

则得到Pi点的平均曲率(H)和高斯曲率(K):

S5包括如下步骤：

S51，对于分割后的每块点云数据，在读取数据的过程中，得到待提取边界点云数据的最大与最小x、y、z坐标值，根据所述坐标值获得点云数据的最小外包围盒；

S52，判断S51中的每一个所述最小外包围盒内点云的个数，如果所述最小外包围盒内点云个数大于1，则对其进行八叉树网格化分，直到划分后得到的网格内点云个数不大于1，统计所述得到的网格和其相邻的空间中的26个网格是空网格的个数，得到统计结果；

S53，根据得到的统计结果判断所述划分后得到的网格是否为边界网格，如果该网格为边界网格，则返回S52，如果不是边界网格，则不用对其进行处理，并将点的内存释放；

S54，提取所述边界网格中的点云数据，得到边界点云，对所有的所述边界点云，根据网格的空间关系，进行索引排序，得到排序后的点云数据以及点云数据的个数，利用DirectX绘制线性结构的函数，绘制边界线，从而得到点云数据的轮廓线特征。

具体实施例：

本发明实施例提供了一种基于曲率从Lidar点云数据提取建筑物轮廓线的方法，可以采用如下步骤进行实施：

S1，点云数据去噪处理；

针对散乱点云没有明显的几何分布特征,呈散乱无序状态,没有规律性,因此为了较快、较好的对其进行处理，先对点云数据进行去噪等处理。

S2，点云数据的动态网格划分；

本发明针对散乱点云分布的不规则性，以及直接网格划分时只点的微分特性无法代表所有点的微分特性，本发明采用动态网格划分的方法，对点云数据进行网格化建立点云之间的拓扑关系。

设点云数据中最大与最小x、y、z坐标值分别为是{Xmax，Ymax,Zmax,Xmin,Ymin,Zmin}，且给定容差e，经多次实验，容差设定为点间距的3倍。并随机抽取3对邻近点云，估算点云数据的密度。

网格划分的具体步骤如下：

S21，首先定义以点(Xmax+e，Ymax+e，Zmax+e)和(Xmax-e，Ymax-e，Zmax-e)为对角点,且表面平行于坐标平面的空间六面体为点云的空间包围盒。

S22，根据点云密度设定阈值t,将点云包围盒沿坐标轴方向按等间隔t划分成空间六面体网格,则在包围盒内，沿X,Y,Z方向，以t为步长分段。可以得到，在三个方向上划分的段数分别为m,n,l：

那么对于空间内任意一个点P，它都会属于某一个位置的基本正方体(i,j,k)内，它的坐标(x_p,y_p,z_p)满足：

S23，经过直接网格划分后，统计每一个网格内点云的个数并计算点云密度，将空网格删掉。求出每一个网格的中心点O_i，在网格数据集中找到离中心点O_i最近的点P_i。对网格内数据在XYZ方向进行排序，判断数据的分布状况。

a.如果网格内的数据分布比较密，点云个数超过了抛物面拟合的需要，但是点云数据之间在X,Y,Z方向的角度偏差在一定的允许范围内，则不必对网格进行调整。

b.如果网格内的数据分布比较密，点云个数超过了抛物面拟合的需要，但是点云数据之间在X,Y,Z方向的角度偏差超出了允许范围，则根据点云偏差的突变，将原网格在一定方向上进行细分，添加网格。

c.如果网格内数据分布均匀，点云个数满足抛物面拟合的需要，并且点云数据之间在X,Y,Z方向的角度偏差在一定的允许范围内，则不必对网格进行调整。

d.如果网格内的数据分布比较稀疏，点云个数不满足抛物面拟合的需要，但是点云数据之间在X,Y,Z方向的角度偏差在一定的允许范围内，则判断邻近网格内数据的分布情况，以原网格为中心，边长2t扩大网格，或沿某一方向将网格与邻近网格合并。然后重新判断，直到满足需要。

S3，点云法矢和曲率的求解；

法矢和曲率是曲面的基本特性,也是曲面特征识别的重要依据之一。空间数据点曲面曲率计算方法一般有五种，分别是坐标转化法、正交步进法、交叉 曲面片法、曲面三角化法和TW法。

本发明将点云数据经过坐标转化法转换后，在每个网格内进行局部曲面拟合。其基本思想是：每一个网格内，求出距离网格中心点最近的点集中的点Pi，利用Pi的k近邻拟合抛物面S,通过计算S的曲率特性来确定该数据点的曲率特性。

每一个网格内点Pi的k近邻,记为Nd(Pi)拟合平面估算测量点法矢,目标函数为

经简单计算，可得到对称的半正定矩阵(3阶矩阵):

式中:Cv最小特征值对应的特征向量单位化,即为曲面在点Pi的法矢Ni。多次实验表明，当k取为24～30,曲率计算一般都能获得较好效果。

在点Pi处建立局部坐标系(u,v,h),且坐标原点为Pi,坐标轴h设为法矢Ni,坐标轴u和v可以在过点Pi垂直于Ni的平面内任意确定。将点Pi的k近邻由全局坐标(xi,yi,zi)转化为局部坐标系(ui,vi,hi)，并用转换后数据拟合如下面的抛物面。

S(u,v)＝a₀+a₁u+a₂v+a₃uv+a₄u²+a₅v²；

抛物面拟合通过求解目标函数实现。应用奇异值分解法可获得拟合抛物面的最小二乘解。

具体实现过程如下：

S31，将网格邻域内点的坐标代入式抛物面S方程中,设每个网格内有k个点，则得到方程组

AX＝Z

式中有

X＝{a₀,a₁,a₂,a₃,a₄,a₅}^T，Z＝{h₀,h₁,.....h_k}^T

S32，采用奇异值分解法求解线性方程组,由X＝A'Z求解曲面方程的系a,b,c。

S33，通过对拟合抛物面S(u,v,h)＝S(u,v,au²+buv+cv²)求一、二阶偏导数，得到Pi点(局部坐标系的原点)的平均曲率与高斯曲率：

Su＝(1,0,2au+bv) Suu＝(0,0,2a) Sv＝(0,1,bu+2cv) Svv＝(0,0,2c) Suv＝(0,0,b)

Pi点的法矢为：

由抛物面的第一、第二基本形式，可以得到：

E＝Su·Su，F＝Su·Sv，G＝Sv·Sv，L＝Suu·n，M＝Suv·n，N＝Svv·n

则得到Pi点的平均曲率(H)和高斯曲率(K):

S4，点云数据的区域分割；

通过计算得到每一个空间网格的平均曲率、高斯曲率和法矢等微分几何信息后，通过曲率、法矢与曲面的关系，采用K-means的方法，对点云数据进行区域分割。

S5，点云数据的轮廓线提取；

对于每一个网格根据它的里面是否有扫描点分为两种情况：网格里面没有 扫描点，网格里面有扫描点。有扫描点网格称为实网格，没有扫描点的网格称为虚网格。

八叉树剖分算法作为一种场景组织方法,广泛应用于计算机图形学系统中,本发明采用八叉树的方法提取分割后点云数据的边界点，可以显著提高提取速度和节约计算机资源。本发明对于分割后的点云数据，采用八叉数算法，提取点云数据的边界点。具体算法描述如下：

S51在导入数据的过程中，得到待提取边界点数据的最大与最小x、y、z坐标值，获得点云数据的最小外包围盒，并将其进行空间八划分，获得八个子包围盒，其中，八叉树划分方法如图2和图3所示。

S52判断每一个网格内点云的个数，如果网格内点个数大于1，则对其进行八叉树网格化分，并判断和其相邻的空间中的26个网格是空网格的个数，判断该网格是否为边界网格。

S53如果该网格为边界网格，则返回S52，如果不是边界网格，则不用对其进行处理，并将点的内存释放。

S54根据网格的空间关系，对其进行索引排序，得到点云数据的边界点。

由于上述步骤对每一个子包围盒点云的个数以及空间相邻空网格的个数进行判断，及时删除在八叉树空间剖分过程中产生的空叶节点和非边界网格，从而减小了存储空间,加快了网格模型的读取速度。

将提取出的边界点，利用其网格空间关系，进行索引排序，得到了排序后的点云数据以及点云数据的个数，利用DirectX绘制线性结构的函数，绘制边界线，从而得到点云数据的轮廓线特征。

通过采用本发明公开的上述技术方案，得到了如下有益的效果：本发明实施例提供的基于曲率从Lidar点云数据提取建筑物轮廓线的方法，包括步骤：点云 数据去噪处理；点云数据的动态网格划分；点云法矢和曲率的求解；点云数据的区域分割；点云数据的轮廓线提取。在地理信息行业数字城市轮廓线提取及三维建模过程中，为自动计算提取规则建筑轮廓线结构提供了一种可靠的解决方案，以改善现阶段纯人工提取建筑轮廓线的工作方式，提高数字城市三维数据生产的效率。

本说明书中的各个实施例均采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似的部分互相参见即可。

本领域人员应该理解的是，上述实施例提供的方法步骤的时序可根据实际情况进行适应性调整，也可根据实际情况并发进行。

上述实施例涉及的方法中的全部或部分步骤可以通过程序来指令相关的硬件来完成，所述的程序可以存储于计算机设备可读取的存储介质中，用于执行上述各实施例方法所述的全部或部分步骤。所述计算机设备，例如：个人计算机、服务器、网络设备、智能移动终端、智能家居设备、穿戴式智能设备、车载智能设备等；所述的存储介质，例如：RAM、ROM、磁碟、磁带、光盘、闪存、U盘、移动硬盘、存储卡、记忆棒、网络服务器存储、网络云存储等。

最后，还需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、商品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、商品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所 述要素的过程、方法、商品或者设备中还存在另外的相同要素。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视本发明的保护范围。

Claims (5)

1.一种基于曲率从Lidar点云数据提取建筑物轮廓线的方法，其特征在于，包括如下步骤：

S1，对散乱分布的点云数据进行去噪处理，得到规则的点云数据；

S2，对所述规则的点云数据，通过设定容差e和阈值t，采用动态层法和局部重划法相结合的点云数据划分方法，对点云特征突变区部分网格进行动态调整，得到动态网格；

S3，将点云数据经过坐标转化法转换后，在每个网格内进行局部曲面拟合，求解点云法矢和曲率；

S4，以所述动态网格为基础，利用S3得到的法矢和曲率，根据曲率、法矢与曲面的关系，采用K-means方法，对点云数据进行区域分割，得到分割的点云数据块；

S5，采用八叉树的方法对所述点云数据块的边界点进行提取，利用DirectX绘制线性结构的函数，绘制边界线，完成点云数据的轮廓线提取。

2.根据权利要求1所述的基于曲率从Lidar点云数据提取建筑物轮廓线的方法，其特征在于，S2包括如下步骤：

S21，定义点(Xmax+e，Ymax+e，Zmax+e)和点(Xmin-e，Ymin-e，Zmin-e)为对角点,表面平行于坐标平面的空间六面体为点云空间包围盒，其中，Xmax，Ymax,Zmax,Xmin,Ymin,Zmin分别为点云数据中最大点与最小点的x、y、z坐标值，e为容差值；

S22，沿坐标轴X,Y,Z方向，按等间隔t，将所述点云空间包围盒划分成多个空间六面体网格则在所述点云空间包围盒内，沿X,Y,Z方向，以t为步长分段，得到在三个方向上划分的段数m,n,l，用如下公式表示为：

$\left\{\begin{array}{c}m=INT\left(\frac{(X\max +e)-(X\min -e)}{t}\right)+1\\ n=INT\left(\frac{(Y\max +e)-(Y\min -e)}{t}\right)+1\\ l=INT\left(\frac{(Z\max +e)-(Z\min -e)}{t}\right)+1\end{array}\right.,$

则，根据上式，可知

对于空间内任意一个点P，都会属于某一个位置的基本正方体(i,j,k)内，坐标(x_p,y_p,z_p)满足：

$\left\{\begin{array}{c}xp\∈(\left(X\min -e+i\×t\right),\left(X\min -e+\left(i+1\right)\×t\right))\\ yp\∈(\left(Y\min -e+j\×t\right),\left(Y\min -e+\left(j+1\right)\×t\right))\\ zp\∈(\left(Z\min -e+k\×t\right),\left(Z\min -e+\left(k+1\right)\×t\right))\end{array}\right.;$

S23，统计所述空间六面体网格内点云的个数并计算点云密度，并判断所述空间六面体网格是否为空，如果是，则删除。

S24，计算每一个网格的中心点O_i，并在网格的数据集中找到距离中心点O_i最近的点P_i；

S25，以所述P_i为基础，对网格内的数据分别在XYZ方向进行排序，根据数据的分布状况判断是否需要对网格进行调整，如果需要，则进行调整，否则，无需处理。

3.根据权利要求2所述的基于曲率从Lidar点云数据提取建筑物轮廓线的方法，其特征在于，S25中，所述根据数据的分布状况判断是否需要对网格进行调整，具体为：

a.如果网格内的数据分布密度大于设定的密度阈值，点云个数超过了抛物面拟合的需要，点云数据之间在X,Y,Z方向的角度偏差小于设定的角度偏差阈值，则无需对网格进行调整；

b.如果网格内的数据分布密度大于设定的密度阈值，点云个数超过了抛物面拟合的需要，但是点云数据之间在X,Y，Z方向的角度偏差超出了角度偏差阈值，则根据点云偏差的突变，将原网格在角度偏差突变的方向上进行细分，添加网格；

c.如果网格内数据分布均匀，点云个数满足抛物面拟合的需要，并且点云数据之间在X,Y,Z方向的角度偏差在角度偏差阈值范围内，则无需对网格进行调整；

d.如果网格内的数据分布密度小于设定的密度阈值，点云个数不能满足抛物面拟合的需要，但是点云数据之间在X,Y,Z方向的角度偏差在角度偏差阈值范围内，则以原网格为中心，边长2t扩大网格，或沿某一方向将网格与邻近网格合并，得到新的网格，并对所述新的网格进行重新判断，直到满足上述c中的情况。

4.根据权利要求1所述的基于曲率从Lidar点云数据提取建筑物轮廓线的方法，其特征在于，

S3具体为：每一个网格内，求出距离网格中心点最近的点集中的点Pi，利用Pi的k近邻拟合抛物面S，通过计算S的曲率特性来确定该数据点的曲率特性；

按照如下方法进行实施：每一个网格内点Pi的k近邻，记为Nd(Pi)拟合平面估算测量点法矢，目标函数为：可得到如下对称的半正定3阶矩阵:

其中:Cv是最小特征值对应的特征向量单位化，即为曲面在点Pi的法矢Ni，多次实验表明，当k取为24～30，曲率计算能获得较好效果；

在点Pi处建立局部坐标系(u,v,h)，且坐标原点为Pi，坐标轴h设为法矢Ni,坐标轴u和v可以在过点Pi垂直于Ni的平面内任意确定，将点Pi的k近邻由全局坐标(xi,yi,zi)转化为局部坐标系(ui,vi,hi)，并用转换后数据拟合如下面的抛物面：

S(u,v)＝a₀+a₁u+a₂v+a₃uv+a₄u²+a₅v²，

抛物面拟合通过求解目标函数实现，应用奇异值分解法可获得拟合抛物面的最小二乘解，上述方法包括如下步骤：

S31，将每一个划分好的网格邻域内点的坐标代入抛物面S方程中,设每个网格内有k个点，则得到如下方程组：

AX＝Z，

其中

X＝{a₀,a₁,a₂,a₃,a₄,a₅}^T，Z＝{h₀,h₁,.....h_k}^T；

S32，采用奇异值分解法求解S31中的线性方程组,由X＝A'Z求解曲面方程的系数a,b,c；

S33，求出系数a,b,c后，对拟合抛物面S(u,v,h)＝S(u,v,au²+buv+cv²)求一、二阶偏导数，得到Pi点的平均曲率H和高斯曲率K：

Su＝(1,0,2au+bv)，Suu＝(0,0,2a)，Sv＝(0,1,bu+2cv)，

Svv＝(0,0,2c)，Suv＝(0,0,b)，

Pi点的法矢为：

由抛物面的第一、第二基本形式，可以得到：

E＝Su·Su，F＝Su·Sv，G＝Sv·Sv，L＝Suu·n，M＝Suv·n，N＝Svv·n，

则得到Pi点的平均曲率(H)和高斯曲率(K):

$K=\frac{LN-M^2}{EG-F^2},H=\frac{EN-2FM+GL}{2(EG-F^2)}.$

5.根据权利要求1所述的基于曲率从Lidar点云数据提取建筑物轮廓线的方法，其特征在于，S5包括如下步骤：

S51，对于分割后的每块点云数据，在读取数据的过程中，得到待提取边界点云数据的最大与最小x、y、z坐标值，根据所述坐标值获得点云数据的最小外包围盒；

S52，判断S51中的每一个所述最小外包围盒内点云的个数，如果所述最小外包围盒内点云个数大于1，则对其进行八叉树网格化分，直到划分后得到的网格内点云个数不大于1，统计所述得到的网格和其相邻的空间中的26个网格是空网格的个数，得到统计结果；

S53，根据得到的统计结果判断所述划分后得到的网格是否为边界网格，如果该网格为边界网格，则返回S52，如果不是边界网格，则不用对其进行处理，并将点的内存释放；

S54，提取所述边界网格中的点云数据，得到边界点云，对所有的所述边界点云，根据网格的空间关系，进行索引排序，得到排序后的点云数据以及点云数据的个数，利用DirectX绘制线性结构的函数，绘制边界线，从而得到点云数据的轮廓线特征。