CN108074232B - 一种基于体元分割的机载lidar建筑物检测方法 - Google Patents

Abstract

本发明提出一种基于体元分割的机载LIDAR建筑物检测方法，该方法为：读取原始机载LIDAR点云数据，形成原始机载LIDAR点云数据集；将原始机载LIDAR点云数据规则化为灰度3D体元数据集；基于连通性和辐射特性相似性准则，将灰度3D体元数据分割并标记为若干个3D连通区域；基于建筑物屋顶和立面的特性，依次检测建筑物屋顶和立面形成的3D连通区域，完成基于体元分割的机载LIDAR建筑物检测；该方法很好地利用了3D体元数据中各体元间隐含的邻域关系及建筑物的特性，有助于基于体元建模理论的机载LIDAR点云数据处理及应用的发展。

Description

一种基于体元分割的机载LIDAR建筑物检测方法

技术领域

本发明属于遥感数据处理技术领域，具体涉及一种基于体元分割的机载LIDAR建筑物检测方法。

背景技术

建筑物是3D地理信息产品中不可或缺的组成部分，因而自动、高精度及快速的建筑物目标检测成为了研究热点。机载激光雷达(Light Detection And Ranging，LIDAR)技术能够提供密集的、精确的、具有地理参考的真三维(3Dimension，3D)点云数据，并包含有回波信号的强度信息。因而机载LIDAR数据特别适合用于3D的目标检测。经典的建筑物检测方法依据其采用的数据结构可分为：基于离散点云、栅格格网及不规则三角网的建筑物检测方法。点云结构完全保留了机载LIDAR数据的原始特征，为真3D数据。但该结构中各激光点相对独立，没有明确地记录关于每个激光点的邻接信息，不能直接得到数据处理需要的邻域信息，由此导致数据处理算法设计困难，运行效率较低。另外，基于点云的建筑物检测结果不易实现矢量化。栅格格网及不规则三角网的同一平面(X，Y)坐标只能对应一个高程(Z)值，该类数据结构表达对3D LIDAR点云数据而言必然存在信息损失，进而影响基于该结构的目标检测结果的完整性。另外，基于栅格格网的建筑物检测方法的检测结果为2D形式。可见，经典建筑物检测方法所采用的数据结构均不利于发挥机载LIDAR真3D的技术优势。体元数据结构是一种真3D数据结构，用其表达LIDAR点云数据不会造成信息损失。同时该结构内部的体元间隐含有几何拓扑关系，因而基于该数据结构的数据处理算法设计相对容易。基于体元结构的机载LIDAR数据的分析多见于林业或地面目标的检测，本发明则创新性地将体元结构与建筑物目标检测相结合，提出了基于体元的3D建筑物目标检测方法。

发明内容

针对现有技术的不足，本发明提出一种基于体元分割的机载LIDAR建筑物检测方法。

一种基于体元分割的机载LIDAR建筑物检测方法，包括以下步骤：

步骤1：读取原始机载LIDAR点云数据，形成原始机载LIDAR点云数据集；

步骤2：将原始机载LIDAR点云数据规则化为灰度3D体元数据集；

步骤2.1：从原始机载LIDAR点云数据中剔除异常数据，得到去除异常数据集；

步骤2.2：将去除异常数据集规则化为灰度3D体元数据集；

步骤3：基于连通性和辐射特性相似性准则，将灰度3D体元数据分割并标记为若干个 3D连通区域；

步骤4：基于建筑物屋顶和立面的特性，依次检测建筑物屋顶和立面形成的3D连通区域，完成基于体元分割的机载LIDAR建筑物检测；

步骤4.1：基于面积特性、高程跳变特性及密度特性，检测建筑物屋顶的3D连通区域，完成建筑物屋顶检测；

步骤4.2：基于缓冲区分析，检测建筑物立面的3D连通区域，完成建筑物立面检测。

2、根据权利要求1所述的基于体元分割的机载LIDAR建筑物检测方法，其特征在于，所述步骤2.1具体包含如下步骤：

步骤2.1.1：统计原始机载LIDAR点云数据中各个激光点高程值的频次，并以直方图的形式可视化显示统计结果；

步骤2.1.2：确定与真实地形对应的最高高程阈值T_h和最低高程阈值T_l；

步骤2.1.3：针对原始机载LIDAR点云数据中各个激光点，若其高程值高于最高高程阈值T_h或低于最低高程阈值T_l，则该激光点为异常数据，进行剔除，否则保留该激光点，最终获得去除异常数据集。

所述步骤2.2具体包含如下步骤：

步骤2.2.1：用去除异常数据集的有向包围盒表示三维空间范围；

步骤2.2.2：根据去除异常数据集中激光点的平均点间距确定体元在x、y、z方向上的分辨率(Δx，Δy，Δz)，即体元大小；

步骤2.2.3：依据体元分辨率(Δx，Δy，Δz)对有向包围盒进行划分，得到3D体元格网，每一个3D体元格网单元即为体元；

步骤2.2.4：将去除异常数据集中各个激光点映射到3D体元格网，进而根据3D体元格网中包含的激光点的强度属性为各体元赋值，得到灰度3D体元数据集。

所述步骤3包括以下步骤：

步骤3.1：依次扫描灰度3D体元数据集中的非0值体元，直到第k个未被标记的非0值体元，其中，k＝1，2，…；

步骤3.2：采用深度优先策略遍历与第k个非0值体元3D连通且体元值相近的所有体元，并标记为L_l，其中，l为标记标签的索引，其中，l＝1，2，…；

步骤3.3：继续扫描灰度3D体元数据集中的未被标记的非0值体元，直到所有的体元都被标记，得到若干个3D连通区域；

所述体元值相近的体元即体元值之差小于灰度差阈值的两个体元。

所述步骤4.1具体包含如下步骤：

步骤4.1.1：基于面积特性剔除非建筑物屋顶的3D连通区域；

步骤4.1.2：基于高程跳变特性剔除非建筑物屋顶的3D连通区域；

步骤4.1.3：基于密度特性剔除非建筑物屋顶的3D连通区域。

所述步骤4.2具体包含如下步骤：

步骤4.2.1：检测各个建筑物屋顶的轮廓；

步骤4.2.2：在水平面上，以任一建筑物屋顶轮廓为中心，以内侧和外侧各一个体元为宽度建立缓冲区；

步骤4.2.3：对任一3D连通区域，若其位于缓冲区内部且其灰度值与该区域对应的建筑物屋顶的灰度值之差小于灰度差阈值，则判定为建筑物立面。

所述根据3D体元格网中包含的激光点的强度属性为各体元赋值的具体过程如下所示：

将含有激光点的体元赋值为激光点强度均值、将不含有激光点的体元赋值为0，进一步将各体元赋值离散化到{0，…，255}，得到各体元值。

本发明的有益效果：

本发明提出一种基于体元分割的机载LIDAR建筑物检测方法，该方法首先将机载LIDAR 点云数据规则化为体元数据，然后依据连通性和辐射特性相似性准则将体元数据分割并标记为若干个3D连通区域，最后利用建筑物屋顶及立面的特性依次检测由建筑物屋顶、立面形成的3D连通区域。体元相对于传统的点云、栅格格网及不规则三角网等建筑物检测中常用的离散LIDAR点云数据表达方式，其是一种简单明晰的空间结构表达且基于该表达的数据处理算法设计更容易。本发明提出的基于体元分割的建筑物检测方法，即是将点云数据中的点云分割及目标检测转换成基于体元分割的目标检测，很好地利用了3D体元数据中各体元间隐含的邻域关系及建筑物的特性，有助于基于体元建模理论的机载LIDAR点云数据处理及应用的发展。

附图说明

图1为本发明具体实施方式中基于体元分割的机载LIDAR建筑物检测方法的流程图；

图2为本发明具体实施方式中原始机载LIDAR点云数据；

其中，(a)为Area2点云数据，(b)为Area3点云数据，(c)为Area2点云数据对应的图像，(d)为Area3点云数据对应的图像；

图3为本发明具体实施方式中将原始机载LIDAR点云数据规则化为灰度3D体元数据集的流程图；

图4为本发明具体实施方式中LIDAR点云数据规则化所得灰度3D体元数据顶视图；

其中，(a)为Area2对应的体元数据，(b)为Area3对应的体元数据；

图5为本发明具体实施方式中将灰度3D体元数据分割并标记为若干个3D连通区域的流程图；

图6为本发明具体实施方式中LABEL(k，u，L_l)程序过程；

图7为本发明具体实施方式中灰度3D体元数据分割结果顶视图；

其中，(a)为Area2的分割结果，(b)为Area3的分割结果；

图8为本发明具体实施方式中Area3的灰度频率直方图；

图9为本发明具体实施方式中完成基于体元分割的机载LIDAR建筑物检测的流程图；

图10为本发明具体实施方式中缓冲区设置的示意图；

其中，黑色体元代表提取的建筑物屋顶轮廓，深灰色体元代表内侧缓冲区，浅灰色体元代表外侧缓冲区；

图11为本发明具体实施方式中完成基于体元分割的机载LIDAR建筑物检测结果；

其中，(a)为Area2的检测结果，(b)为Area3的检测结果。

具体实施方式

下面结合附图对本发明具体实施方式加以详细的说明。

本实施方式中，在CPU dual-core 2.4GHz、内存4GB、Windows 7旗舰版系统上使用MATLAB 7.11.0平台编程实现该方法，并进一步通过对该方法的精度评定验证方法的有效性。

一种基于体元分割的机载LIDAR建筑物检测方法，如图1所示，包括以下步骤：

步骤1：读取原始机载LIDAR点云数据，形成原始机载LIDAR点云数据集。

本实施方式中，采用国际摄影测量与遥感协会(International Society forPhotogrammetry and Remote Sensing，ISPRS)第III/4工作组提供的两组(Area2和Area3，如图2所示)专门用于目标分类算法测试的城区样本数据作为实验数据，以检验方法的有效性和可行性。实验数据由LeicaALS50机载激光扫描仪获取(航高500米，视场角45)。这两组数据中包含被树木围绕的高层城市住宅建筑物和带有小附属建筑物的居民区。点云数据密度为4点/m²。

本实施方式中，定义原始机载LIDAR点云数据集P＝{p_i(x_i，y_i，z_i)，i＝1，...，n}，其中，i是原始机载LIDAR点云数据的索引，n是原始机载LIDAR点云数据的个数，p_i是第i个原始机载LIDAR点云数据，其坐标为(x_i，y_i，z_i)。

步骤2：将原始机载LIDAR点云数据规则化为灰度3D体元数据集，具体流程如图3所示。

步骤2.1：从原始机载LIDAR点云数据中剔除异常数据，得到去除异常数据集。

步骤2.1.1：统计原始机载LIDAR点云数据中各个激光点高程值的频次，并以直方图的形式可视化显示统计结果。

步骤2.1.2：确定与真实地形对应的最高高程阈值T_h和最低高程阈值T_l。

步骤2.1.3：针对原始机载LIDAR点云数据中各个激光点，若其高程值高于最高高程阈值T_h或低于最低高程阈值T_l，则该激光点为异常数据，进行剔除，否则保留该激光点，最终获得去除异常数据集。

本实施方式中，去除异常数据集记做Q＝{q_i′(x_i′，y_i′，z_i′)，i′＝1，...，t}，其中，i′是去除异常数据集中数据的索引，t是去除异常数据集中数据的个数，q_i′是去除异常数据集中第i′个数据，其坐标为(x_i′，y_i′，z_i′)。

本实施方式中，最高高程阈值T_h和最低高程阈值T_l为常数，其取值需根据原始机载LIDAR点云数据的空间分布情况确定。

步骤2.2：将去除异常数据集规则化为灰度3D体元数据集。

步骤2.2.1：用去除异常数据集的有向包围盒表示三维空间范围。

本实施方式中，去除异常数据集Q的有向包围盒为长方体，其底面大小可通过求取去除异常数据集Q在XY平面上投影的最小外接矩形确定，其高可由(z_max-z_min)确定。

其中，z_max是去除异常数据集Q中激光点的z坐标的最大值，z_min是去除异常数据集Q中激光点的z坐标的最小值，z_max＝max{z_i′，i′＝1，...，t}，z_min＝min{z_i′，i′＝1，...，t}。

步骤2.2.2：根据去除异常数据集Q中激光点的平均点间距确定体元在x、y、z方向上的分辨率(Δx，Δy，Δz)，即体元大小。

本实施方式中，体元在x、y、z方向上的分辨率Δx、Δy、Δz的计算公式如式(1)所示：

其中，S_xy＝{(x_i′，y_i′)，i′＝1，...，t}为去除异常数据集Q在XOY平面上的投影所得二维点集， C(S_xy)是点集S_xy的凸壳，A(C(S_xy))是凸壳C(S_xy)面积。

步骤2.2.3：依据体元分辨率(Δxx，Δy，Δz)对有向包围盒进行划分，得到3D体元格网，每一个3D体元格网单元即为体元。

本实施方式中，基于体元分辨率(Δxx，Δy，Δz)就可以将有向包围盒划分为3D体元格网，用 3D体元阵列表示。设V是3D体元阵列中的体元集合，如式(2)所示：

V＝{v_j(r_j，c_j，l_j)，j＝1，…，m}， (2)

其中，j是体元索引；m是体元数；v_j是第j个体元的体元值；(r_j，c_j，l_j)是第j个体元在体元阵列中的坐标(行、列和层号)。X方向上的体元个数为R，Y方向上的体元个数为C，Z方向上的体元个数为L。其中，R、C、L由公式(3)所示：

其中，

为向上取整操作符，x_max＝max{x_i′，i′＝1，...，t}，x_min＝min{x_i′，i′＝1，...，t}， y_max＝max{y_i′，i′＝1，...，t}.y_min＝min{y_i′，i′＝1，...，t}。

由此可以得出，体元数m如式(4)所示：

m＝R*C*L (4)

步骤2.2.4：将去除异常数据集Q中各个激光点映射到3D体元格网，进而根据3D体元格网中包含的激光点的强度属性为各体元赋值，得到灰度3D体元数据集。

本实施方式中，将去除异常数据集Q中各个激光点映射到3D体元格网，进而根据3D体元格网中包含的激光点的强度属性为各体元赋值。其中，将含有激光点的体元赋值为激光点强度均值、将不含有激光点的体元赋值为0，如式(5)所示：

其中，

为向下取整操作符。进一步将各体元赋值离散化到{0，…，255}，得到各体元值。由此，得到灰度3D体元数据集，完成对去除异常数据集的规则化。

本实施方式中，由点云规则化所得体元数据顶视图如图4所示，灰度值代表不同的体元强度值，可以看出，建筑物屋顶的强度值相对一致且和其它目标的强度值存在差异。

步骤3：基于连通性和辐射特性相似性准则，将灰度3D体元数据分割并标记为若干个 3D连通区域，具体流程如图5所示。

步骤3.1：依次扫描灰度3D体元数据集中的非0值体元，直到第k个未被标记的非0值体元，其中，k＝1，2，…。

步骤3.2：采用深度优先策略遍历与第k个非0值体元3D连通且体元值相近的所有体元，并标记为L_l，其中，l为标记标签的索引，l＝1，2，…。

本实施方式中，针对体元值为v_k＝u的第k个体元，采用深度优先策略遍历与第k个非0 值体元3D连通且体元值相近的所有体元，如图6所示，调用程序LABEL(k，u，L_l)标记与第k 个非0值体元3D连通且强度值相近的体元的标签为L_l。体元值相近的体元即体元值之差小于强度差阈值Tg的两个体元。

步骤3.3：继续扫描灰度3D体元数据集中的未被标记的非0值体元，直到所有的体元都被标记，得到若干个3D连通区域。

本实施方式中，分割的目标就是将连通且辐射特性相似的体元合并到一个3D连通区域。假设3D体元阵列V共包含l个3D连通区域，分割的任务就是将l个标签分配给V中的各个体元，以便属于同一个3D连通区域的体元具有相同的标签，而属于不同的3D连通区域的体元则具有不同的标签。具体实施过程中，基于连通且强度值相近的体元为同一类别的思想，利用3D连通分量标记算法将V分割并标记为l个3D连通区域，得到的分割结果如图7所示。可以看出绝大多数属于单体建筑物的体元都被分割到一个3D连通区域中。

本实施方式中，在上述标记过程中应用不同的空间邻域尺度及不同的强度差阈值会得到不同的分割结果并进而影响后继的建筑物提取结果。最佳空间邻域尺度将在实验中确定。最佳强度差阈值由下述方案确定(以实验数据中Area3为例)：

统计3D体元阵列V中的非0值体元的灰度值的频率，并以直方图形式显示，如图8所示，图中包含4个正态分布，其峰值分布对应1、28、81和133。为了保证属于单体建筑的体元均被分割到一个3D连通区域，对应于建筑物的第三个正态分布被用来计算最优的强度差阈值T_g。第三个正态分布灰度值的范围为[61，103]，其均值μ和标准差σ分别是81.9和11.5。1.84σ被用作最优的强度差阈值T_g。选择乘数1.84的原因在于所有的建筑物的灰度等级都在 1.84σ范围内。

本实施方式中，参考数据采用了ISPRS第III/4工作组提供的建筑物标准数据(被准确分类为建筑物点集和非建筑物点集的实验数据)，以定量评价本发明方法的计算精度。

本发明提出的建筑物检测方法的成果是以体元形式表示的，而参考数据中的建筑物则是以离散的LIDAR激光点云数据表达的，为和参考数据做对比以评价本发明提出的方法精度，首先统计本方法所检测的建筑物体元内包含的原始机载LIDAR点云数据的个数，然后和参考数据进行比对进而用I类误差(将建筑物激光点错分为非建筑物激光点比例)、II类误差(将非建筑物激光点错分为建筑物激光点比例)、总误差(错分的建筑物激光点的比例)、正确率(正确检测的建筑物激光点数占检测结果中建筑物激光点总数的比例)、完整度(正确检测的建筑物激光点数占标准数据中建筑物激光点总数的比例)、质量及Kappa系数来定量评价本发明所提出的建筑物检测方法的有效性。

表1为本实施例中空间邻域尺度为6、8、26、56以及80时，基于2个实验数据的体元数据进行分割进而进行检测其中的建筑物，对应的建筑物检测结果的精度指标。该表中的数据旨在考查不同领域尺度对建筑物检测结果的影响，并由此确定最优邻域尺度。

表1不同邻域尺度的建筑物检测结果的精度指标

由表1可知，6、18、26、56和80的平均Kappa系数分别为34.5％、66.5％、71.3％、89.1％和87.3％，这说明：(1)56邻域对应最大的Kappa系数，因此，从Kappa系数指标来看，56 邻域是最佳邻域尺度。(2)邻域尺度的增加并不意味着检测精度的必然提高。该算法的思想是，建筑物信息可以通过基于体元阵列中定义的连通性和强度相似性来传播。以6邻域为例，建筑物信息只能向其上、下、左、右、前、后6个方向传播，由此导致只有平顶建筑物(例如， Area 2)上的体元才能被合并到一个3D连通区域并被正确检测，而位于尖顶建筑物(例如，Area 3)上的体元可能被划分到多个3D连通区域并因此被后继的面积、高程跳变等特性而误判。这可以解释为什么使用6邻域产生了最高的I类误差(见表1中第4列不同邻域尺度对应的I类误差)，也可以解释为什么6邻域的Area 2的Kappa系数远高于Area 3。随着邻域尺度的增加，传播的方向增加，更多的体元被分类为建筑物，使用18、26、56邻域比6邻域产生更好的结果。但如果邻域尺度太大，有些非建筑物体元可能被误判为建筑物，因此导致II类误差的增加(见使用80-邻接的II型误差面积2和3)，这可以解释为何80邻域对比56邻域的精度反而有所下降。(3)6、18、26、56和80的平均总误差分别为25.8％、14.9％、12.8％、5.1％和5.9％，这说明56邻域对应最小的平均总误差。因此，从总误差指标来看，56邻域也是最佳邻域尺度。

步骤4：基于建筑物屋顶和立面的特性，依次检测建筑物屋顶和立面形成的3D连通区域，完成基于体元分割的机载LIDAR建筑物检测，具体流程如图9所示。

步骤4.1：基于面积特性、高程跳变特性及密度特性，检测建筑物屋顶的3D连通区域，完成建筑物屋顶检测。

本实施方式中，建筑物屋顶的特性为：具有一定的面积；与其周围的地形存在高程差；和其它目标(如植被等)在空间分布上存在密度差异。根据以上面积、高程跳变及密度特性建立从l个3D连通区域中搜寻建筑物屋顶对应的3D连通区域的方案：

步骤4.1.1：基于面积特性剔除非建筑物屋顶的3D连通区域。

本实施方式中，令原始机载LIDAR点云数据集P中的最小建筑物面积为A_min，令原始机载LIDAR点云数据集P中的最大的建筑物面积为A_max。

原始机载LIDAR点云数据集P中的建筑物面积是指建筑物在水平平面的投影面积。

对任一3D连通区域，若其水平投影面积大于等于A_min并且小于等于A_max，则该3D连通区域判定为建筑物屋顶，予以保留，否则剔除该3D连通区域。

在本实施方式中，A_min和A_max为常数，由用户根据给定的原始机载LIDAR点云数据的情况定义。

步骤4.1.2：基于高程跳变特性剔除非建筑物屋顶的3D连通区域。

对任一3D连通区域，若其轮廓线的平均高程和其周围地形的平均高程大于给定的高程阈值T_e，则该3D连通区域判定为建筑物，予以保留，否则剔除该3D连通区域。其中，某个3D连通区域的周围地形的平均高程可用下述方案获取：利用结构元素[1 1 1；1 1 1；1 11]对该 3D连通区域做3D形态学膨胀操作，记膨胀处理结果的外轮廓线为C_s＝{v_a(r_a，c_a，l_a)，a＝ 1，...，q}，其中，s代表第s个3D连通区域，a是第s个3D连通区域的外轮廓线中各个体元的索引，q是第s个3D连通区域的外轮廓线包含的体元个数。

搜寻与其平面坐标相同且高程低于l_a的非0值体元，上述体元的平均高程即是第s个3D连通区域的周围地形的平均高程；

本实施方式中，高程阈值T_e为常数，其取值2米。

步骤4.1.3：基于密度特性剔除非建筑物屋顶的3D连通区域。

对任一3D连通区域，若其密度大于给定的密度阈值T_d，则该3D连通区域判定为建筑物屋顶，予以保留，否则剔除该3D连通区域；

本实施方式中，密度阈值T_d为常数，其取值可根据各个3D连通区域的密度分布情况确定。

步骤4.2：基于缓冲区分析，检测建筑物立面的3D连通区域，完成建筑物立面检测。

本实施方式中，建筑物立面的特性为：垂直于建筑物屋顶轮廓，位于建筑物屋顶轮廓的周围一定范围内。根据上述特点，确定搜寻建筑物立面对应的3D连通区域的方案：

步骤4.2.1：检测各个建筑物屋顶的轮廓。

步骤4.2.2：在水平面上，以任一建筑物屋顶轮廓为中心，以内侧和外侧各一个体元为宽度建立缓冲区，如图10所示。

步骤4.2.3：对任一3D连通区域，若其位于缓冲区内部且其灰度值与该区域对应的建筑物屋顶的灰度值之差小于灰度差阈值T_g，则判定为建筑物立面。

应用本发明方法所得的建筑物检测结果如图11所示，其中，(a)为Area2的检测结果， (b)为Area3的检测结果。其中，建筑物检测结果是建筑物体元形式，单个体元的大小为0.4m×0.4m×0.4m，如图10中的黑色立方体所示。Area2和Area3检测到的建筑物体元数分别为25933和37589个。上述建筑物体元直接可用作建筑物模型，其是一种新形式的3D建筑物体元模型。

表2为本实施方式中，以参考数据为标准对2个实验数据的56邻域尺度下的建筑物检测精度进行的定量评价。

表2建筑物检测结果的精度

由表2可知：建筑物检测的平均质量、平均完整度及平均正确率分别为89.2％、89.7％及 96.8％。从而验证了本发明提出的方法的有效性。

最后应说明的是：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明权利要求所限定的范围。

Claims (6)

1.一种基于体元分割的机载LIDAR建筑物检测方法，其特征在于，包括以下步骤：

步骤1：读取原始机载LIDAR点云数据，形成原始机载LIDAR点云数据集；

步骤2：将原始机载LIDAR点云数据规则化为灰度3D体元数据集；

步骤2.1：从原始机载LIDAR点云数据中剔除异常数据，得到去除异常数据集；

步骤2.2：将去除异常数据集规则化为灰度3D体元数据集；

步骤3：基于连通性和辐射特性相似性准则，将灰度3D体元数据分割并标记为若干个3D连通区域；

步骤4：基于建筑物屋顶和立面的特性，依次检测建筑物屋顶和立面形成的3D连通区域，完成基于体元分割的机载LIDAR建筑物检测；

步骤4.1：基于面积特性、高程跳变特性及密度特性，检测建筑物屋顶的3D连通区域，完成建筑物屋顶检测；

步骤4.2：基于缓冲区分析，检测建筑物立面的3D连通区域，完成建筑物立面检测；

步骤4.2.1：检测各个建筑物屋顶的轮廓；

步骤4.2.2：在水平面上，以任一建筑物屋顶轮廓为中心，以内侧和外侧各一个体元为宽度建立缓冲区；

步骤4.2.3：对任一3D连通区域，若其位于缓冲区内部且其灰度值与该区域对应的建筑物屋顶的灰度值之差小于灰度差阈值，则判定为建筑物立面。

2.根据权利要求1所述的基于体元分割的机载LIDAR建筑物检测方法，其特征在于，所述步骤2.1具体包含如下步骤：

步骤2.1.1：统计原始机载LIDAR点云数据中各个激光点高程值的频次，并以直方图的形式可视化显示统计结果；

步骤2.1.2：确定与真实地形对应的最高高程阈值T_h和最低高程阈值T_l；

步骤2.1.3：针对原始机载LIDAR点云数据中各个激光点，若其高程值高于最高高程阈值T_h或低于最低高程阈值T_l，则该激光点为异常数据，进行剔除，否则保留该激光点，最终获得去除异常数据集。

3.根据权利要求1所述的基于体元分割的机载LIDAR建筑物检测方法，其特征在于，所述步骤2.2具体包含如下步骤：

步骤2.2.1：用去除异常数据集的有向包围盒表示三维空间范围；

步骤2.2.2：根据去除异常数据集中激光点的平均点间距确定体元在x、y、z方向上的分辨率(Δx,Δy,Δz)，即体元大小；

步骤2.2.3：依据体元分辨率(Δx,Δy,Δz)对有向包围盒进行划分，得到3D体元格网，每一个3D体元格网单元即为体元；

步骤2.2.4：将去除异常数据集中各个激光点映射到3D体元格网，进而根据3D体元格网中包含的激光点的强度属性为各体元赋值，得到灰度3D体元数据集。

4.根据权利要求1所述的基于体元分割的机载LIDAR建筑物检测方法，其特征在于，所述步骤3包括以下步骤：

步骤3.1：依次扫描灰度3D体元数据集中的非0值体元，直到第k个未被标记的非0值体元，其中，k＝1,2,…；

步骤3.2：采用深度优先策略遍历与第k个非0值体元3D连通且体元值相近的所有体元，并标记为L_l，其中，l为标记标签的索引，其中，l＝1,2,…；

步骤3.3：继续扫描灰度3D体元数据集中的未被标记的非0值体元，直到所有的体元都被标记，得到若干个3D连通区域；

所述体元值相近的体元即体元值之差小于灰度差阈值的两个体元。

5.根据权利要求1所述的基于体元分割的机载LIDAR建筑物检测方法，其特征在于，所述步骤4.1具体包含如下步骤：

步骤4.1.1：基于面积特性剔除非建筑物屋顶的3D连通区域；

步骤4.1.2：基于高程跳变特性剔除非建筑物屋顶的3D连通区域；

步骤4.1.3：基于密度特性剔除非建筑物屋顶的3D连通区域。

6.根据权利要求3所述的基于体元分割的机载LIDAR建筑物检测方法，其特征在于，所述根据3D体元格网中包含的激光点的强度属性为各体元赋值的具体过程如下所示：

将含有激光点的体元赋值为激光点强度均值、将不含有激光点的体元赋值为0，进一步将各体元赋值离散化到{0,…,255}，得到各体元值。

Images