CN110570428B - 一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法及*** - Google Patents

Abstract

本发明公开了一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法及***，着重解决从倾斜影像生成的彩色三维点云中如何准确获取建筑物屋顶三维面片的问题。一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，包括如下步骤：地面过滤和植被过滤；去除影像点云中的地面和植被；建筑物单体化；剔除非建筑物点，并实现建筑物聚类；去除建筑物立面；屋顶面片的分割、修补和细化。

Description

一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法 及***

技术领域

本发明涉及影像应用技术领域，尤其涉及一种从大规模影像密集匹配点云分割建筑物屋 顶面片的方法及***。

背景技术

城市作为时空信息的重要主体正发生着日新月异的变化。城市建筑物三维模型作为城市 地形地貌、地面上建筑物等等一系列城市空间当中信息的可视化表达，展示了目标的几何形 态、属性、空间位置及纹理等各种信息，是城市三维地理信息***的基础。快速提取建筑物和进行建筑物变化检测在GIS数据库更新、土地利用、数字化城市等诸多方面具有重要的作 用。

本发明方法以倾斜摄影测量***获取的三维彩色点云为数据来源，研究建筑物的屋顶面 片分割技术。下面分别对倾斜摄影测量***、基于机载LiDAR的建筑物屋顶提取技术进行详 细介绍。

倾斜摄影测量***是结合计算机视觉的倾斜摄影测量技术作为在传统摄影测量技术的基 础上所发展出的新技术，能够同时获取多个角度的地面影像，获取建筑物顶部影像及侧面的 丰富纹理。该技术的优势在于能够全方位全角度地提现地面物体真是的情况与形态，捕获地物尤其是建筑物高精度的侧面纹理信息，更能够对影像中的地物进行定位、建模，构建真实 的城市三维模型，其应用领域变得越来越广泛。

目前，倾斜成像***基本上由“1+N”相机构成，即一台下视相机、N台倾斜相机，从城 市三维建模的角度分析，“1+N”倾斜成像***能够获取覆盖完备的影像。目前最常见的是4 台倾斜相机，即前后、左右各2台。倾斜成像***基本都集成了高精度的POS***，进行影 像采集的同时记录影像曝光时刻相机的位置和姿态数据。

为了进一步改善倾斜成像***影像覆盖的完备性，采用的技术策略是增加倾斜相机的数 量，比如Track Air公式的OCTOBLIQUE MIDAS，增加至8台，由一个下视相机与八个倾斜角 度为45°的倾斜相机，以360°全景进行拍摄，额外添加的四个相机能够为影像全方位无死 角的覆盖所摄区域创造了有利条件。

综上，多视角倾斜相机一般搭载在航空飞行器上，多角度获取地物影像。从多视角倾斜 摄影技术特点上来看，主要有以下4个特点：

(1)同一时刻五镜头拍摄的影像均不存在重叠；

(2)在整个测区中，测视影像和垂直影像间的重叠关系是复杂未知的；

(3)同一航线具有重叠区的垂直影像和后视影像间存在约为180°的旋转；

(4)不存在交叉航线时，具有重叠区域的垂直影像和左右视影像间存在约为±90°的旋 转。

基于图像的稠密点云三维重建技术，又称多视图立体(multi-view stereo,MVS)，使 用的是从SFM算法中得到的定标图像，即已通过计算或其他方法得到图像所对应相机的参数。 目前，后处理商业软件主要有法国Acute3D公司所生产的ContextCapture软件，香港科技大学权龙教授建立的Altizure平台，等。基于图像的稠密点云三维重建技术主要流程如下：

(1)特征点检测与匹配；

(2)相机参数求解；通常使用SFM意义下的相机求解方法，在不精确的初始化情况下， 利用bundle adjustment约束，求得精确的相机参数，并得到稀疏的点云模型；

(3)稠密点云的三维重建；现阶段主要流行的是PMVS算法，即利用SFM得到的稀疏点 云和相机参数，通过PMVS基于面片生成的三维重建算法，生成稠密匹配下的三维点云模型。

(4)网格空洞修补和添加纹理映射。

结合计算机视觉的倾斜摄影测量技术突破了遥感技术、传统摄影测量与测绘行业的限制， 发展前景广阔。该技术在我国的摄影测量与遥感领域的市场应用的潜力巨大，对于“智慧城 市”和“数字城市”的构建具有极大的推动作用。在本发明方法中，倾斜摄影技术是研究的基础，通过该手段获取的稠密彩色三维点云是建筑物屋顶轮廓线提取的对象。

通常来说，基于机载LiDAR的建筑物屋顶提取技术主要为2步：建筑物单体化和屋顶3D 面片生成。建筑物单体化方面，常用两种思路。1)直接从原始点云中识别建筑物；2)分别 过滤掉地面、植被之后得到建筑物点云。而屋顶面片提取技术按照建筑立面去除、面片分割 进行。

发明内容

为了解决上述的技术问题，本发明的目的是提供一种从大规模影像密集匹配点云分割建 筑物屋顶面片的方法及***，着重解决从倾斜影像生成的彩色三维点云中如何准确获取建筑 物屋顶三维面片的问题。

为了达到上述的目的，本发明采用了以下的技术方案：

一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，包括如下步骤：

1)地面过滤和植被过滤；去除影像点云中的地面和植被；

2)建筑物单体化；剔除非建筑物点，并实现建筑物聚类；

3)去除建筑物立面；

4)屋顶面片的分割、修补和细化。

作为优选，所述地面过滤，利用布料滤波算法将点云数据中的地面点与非地面点分离， 得到非地面点。

作为优选，所述地面过滤的具体实现过程如下：

1)将去除离群点的原始激光点云做180度翻转；

2)初始化布料格网；用定义的网格分辨率确定布料粒子总数；布料的初始位置被放置在 点云最高点的上方；

3)将激光点云和网格粒子投影到一个水平面，发现每个粒子的最近邻激光点，并将激光 点的高程记为交点高程值；

4)对每个粒子，计算它由于重力作用后的位置；如果它的高程值大于交点高程值，它可 以继续移动；如果小于等于交点高程值，则这个粒子的高程值被设置为交点高程值，并被标 记为不可移动；

5)对于每个格网粒子，计算由于粒子间的相互作用力造成的位置偏移；

6)迭代执行上述步骤4)至5)，直到所有粒子的最大高程方差足够小或者超过了定义 的最大迭代次数；

7)计算粒子和点云之间的距离；

8)分离地面和非地面点；对于任何一个激光点，如果与对应网格粒子点的距离小于阈值， 则被识别为地面点；否则，为非地面点。

作为优选，所述地面过滤的具体实现过程还包括一后处理步骤：搜索每个可移动粒子的 四近邻粒子，如果发现有不可移动粒子，则比较两者对应激光点的高程；如果最小高程差小 于阈值，则可移动粒子继续移动，直到被设置为不可移动；迭代执行这个步骤，直到遍历所 有可移动粒子。

作为优选，所述植被过滤利用绿信比G_r和植被系数指标VI识别绿色或者近绿色点，再 利用基于法向量的变化检测算法对识别点集进行优化，进而剔除植被点。

作为优选，所述绿信比即绿色波段的DN值占RGB的DN值之和的比值；绿信比和植被系 数的定义如下式所示：

式中：VI和G_r分别为植被系数和绿信比；R、G、B为点云的3个颜色通道值，R′、G′、B′分别为R、G、B所占比例，R′、G′、B′定义如下式所示：

作为优选，所述建筑物单体化，利用欧式聚类算法进行连通域分析得到初步建筑物单体 化结果，再利用基于高程阈值的约束项进一步优化。

作为优选，建筑物立面的去除，利用基于法向量的点特征f去除建筑物立面，得到建筑 物屋顶点云；点特征f的定义如下公式所示：

f＝1-|n_p·e_z|

式中，n_p为点p处的法向量，e_z为单位列向量(0,0,1)。

作为优选，所述屋顶面片的分割是指，对单栋建筑物屋顶点云，利用基于点间距的区域 增长算法实现屋顶面片分割；所述屋顶面片的修补是指，对屋顶面片分割后的输出，利用基 于PCA的面片修补算法修补孔洞，并平整面片；所述屋顶面片的细化是指，对屋顶面片分割后的输出，利用二维格网，实现屋顶面片的去重处理。

一种从大规模影像密集匹配点云分割建筑物屋顶面片的***，采用如上所述的方法，包 括：

第一模块，用于使用布料滤波算法分离地面点，得到非地面点；

第二模块，用于利用绿信比和植被系数得到初始植被点，基于法向量的变化检测进一步 约束植被点；

第三模块，用于执行欧式聚类算法和高程阈值约束，从而实现建筑物单体化；

第四模块，用于实现建筑物立面消除算法，和实现屋顶面片分割，获得初始的屋顶面片；

第五模块，用于实现基于PCA的面片修补算法，和基于二维格网的面片去重，得到平整 的屋顶面片。

借用LiDAR点云处理的一些思想，本发明提出了一种从影像点云中全自动提取建筑物三 维轮廓线的算法框架。由于影像点云噪声很多，相比LiDAR点云，算法更加复杂。按照流程， 本框架主要由三部分组成：建筑物单体化、屋顶面片分割、边界线生成。其中，每一部分都由2-3个基本算法组成。下面详细介绍每一部分的内容。

1、建筑物单体化

为了实现从大规模影像点云中准确提取出建筑物三维轮廓线，建筑物单体化是必不可少 的一步。因为屋顶轮廓线的提取是以单体建筑物为单元，在此基础上，再进行轮廓线的提取。 借用LiDAR点云“自下而上”分离建筑物的思想，本方法依次通过分离地面、植被来实现建筑物单体化。

1.1地面过滤

从数据性质上分析，影像点云类似于机载激光点云，因此可以利用机载激光点云中地面 识别算法进行地面去除。这里采用Wuming Zhang于2016提出的布料滤波算法，论文详见Wuming Z,Jianbo Q,Peng W,et al.An Easy-to-Use Airborne LiDAR DataFiltering Method Based on Cloth Simulation[J].Remote Sensing,2016,8(6):501-。该算法参数设置较少、算法泛化能力强、能得到稳定的地面过滤效果。

布料滤波算法模拟了一个简单的物理过程。假设一块布料在空中由于重力作用掉落下来。 如果布料足够柔软，能够依附于物体表面，那么布料的最后形状就是场景的数字表面模型。 现在将点云场景翻转180度，使点云中的地面部分位于上方。不同的是，这里布料具有一定的刚度，那么布料最后形成的形状就是数字地形模型，即地面。分析布料的粒子和相应激光 点云的相互作用就可以确定布料的最后形状，从而用它将地面点和非地面点分离。基于这种 思想，该算法用到了三维计算机图形学中的布料模拟技术，即通过计算机编程来模拟布料。 这里的布料是一种由大量有固定质量但无空间大小的粒子相连形成的格网。格网的顶点是这些物理粒子，格网的边是一种服从霍克定律的“弹簧”。

该算法的具体实现过程如下：

1)将去除离群点的原始激光点云做180度翻转；

2)初始化布料格网。用户自定义的网格分辨率确定布料粒子总数；布料的初始位置被放 置在点云最高点的上方；

3)将激光点云和网格粒子投影到一个水平面。发现每个粒子的最近邻激光点，并将激光 点的高程记为交点高程值；

4)对每个粒子，计算它由于重力作用后的位置。如果它的高程值大于交点高程值，它可 以继续移动；如果小于等于交点高程值，则这个粒子的高程值被设置为交点高程值，并被标 记为不可移动；

5)对于每个格网粒子，计算由于粒子间的相互作用力造成的位置偏移；

6)迭代执行4)-5)，直到所有粒子的最大高程方差足够小或者超过了用户定义的最大 迭代次数；

7)计算粒子和点云之间的距离；

8)分离地面和非地面点。对于任何一个激光点，如果与对应网格粒子点的距离小于阈值， 则被识别为地面点；否则，为非地面点；

为了使算法适用于有较大梯度变化的地形，加入了一个后处理算法。原理如下：搜索每 个可移动粒子的四近邻粒子，如果发现有不可移动粒子，则比较两者对应激光点的高程。如 果最小高程差小于阈值，则可移动粒子继续移动，直到被设置为不可移动。迭代执行这个步骤，直到遍历所有可移动粒子。

1.2植被过滤

本算法没有采用传统激光点云的思想。从数据内容分析，数据源是彩色三维点云，具有 正确的RGB信息；由于植被通常是绿色或者近绿色，因此植被的绿光波段DN(DigitalNumber)值一般比红光和蓝光波段DN值高。结合这两点，本发明提出使用绿信比G_r来识别植被。再 利用植被系数指标VI进行优化，从而实现了剔除植被点。绿信比即绿色波段的DN值占RGB 的DN值之和的比值；植被系数VI的定义公式是经验公式的结果。绿信比和植被系数的定义 如下式所示。

式中：VI和G_r分别为植被系数和绿信比；R、G、B为点云的3个颜色通道值，R′、G′、B′分别为R、G、B所占比例，G′、G′、B′定义如下式所示：

值得注意的是，现实世界中存在一些特殊的建筑。比如，在棚户区中存在一些建筑物屋 顶为蓝色；某些居民楼或者学校的建筑物墙面上有爬山虎，幼儿园墙面为彩色等情况。为了 防止植被点过滤时，将非植被点误剔除，本发明采用基于法线变化的方法来进行优化。

首先，估计点的法向量，这里我们采用常见的基于协方差矩阵的计算方法。如果一个区 域内的点足够多，那么它就能构造出一个表面，那么我们就可以用它来估计法线。如公式(1) 所示，我们以点云中的点p为中心，建立一个邻域N_p。其中，P为原始点集，p∈P,d(p，q) 表示两点之间的距离，r是以p为圆心的搜索半径。

N_p＝{q|q∈P，d(p，q)＜r} (1)

然后，我们构造邻域内的协方差矩阵C_p，定义如公式(2)所示。p是邻域p中所有 点的中心。计算协方差矩阵的特征值，并排序λ₁<λ₂<λ₃，最小特征值λ₁对应的特征向量被 视为点p的法向量。

M.Pauly,Point Primitives for Interactive Modeling and Processing of3D Geometry. Konstanz,Germany:Hartung-Gorre,2003.中提到利用法向量的协方差矩阵的特征值可以计 算法向量的最大变化，这里我们采用这种方法。

我们构造邻域

内所有点法向量的协方差矩阵，如公式(3)所示。同C_p，计算该 协方差矩阵的特征值，并排序/>

可以定量表示该点p在高斯球面上法线的最大变化。

因此，我们可以设定阈值Tⁿ，如果点的

值小于Tⁿ，则把它视为非植被点。

1.3建筑物单体化

对于激光点云，过滤掉地面、植被之后，通常剩余点包含建筑物(含围墙)、车辆和电 力设施，比如电线、电杆等。但是，影像点云中还包含一些块状的噪声点。利用欧式聚类算 法，可以剔除大部分非建筑物点，并实现建筑物聚类。

欧式聚类算法是一种基于点间距的聚类算法，算法实现步骤如下：

1)设置点间距的阈值T，所有点都被设置为未访问点；

2)设置点集P为空。以点云中任意一个未访问点为起点，标记为访问点，寻找与该点的 空间距离小于或等于T的所有点，将添加到P；

3)以点集P中任意一个未访问点为起点，重复2)，直到点集P不能增加新点。点集P为一个聚类；

4)重复执行2)-3)，直到所有点都被标记为访问点。

完成聚类工作。

但是，还是会有部分非建筑物点被识别为建筑，比如残留地面、类平面的绿色植被、连 接建筑物的墙体等，这些误分类的存在会影响屋顶面片提取的准确性。为了去除这些误分类， 首先找到每类点集的高程最小值Zmin。结合建筑屋顶常识，设置高程阈值Tz＝Zmin+4。经过 这些步骤，能够实现较高精度的建筑物单体化。

2、屋顶面片分割

经过建筑物单体化之后，我们可以以单栋建筑为操作单元进行屋顶面片提取。相比机载 激光点云，影像点云中建筑物信息更加丰富，不仅包括建筑物立面信息，还有可能包括建筑 物部分室内信息。

2.1建筑物立面去除

影像点云数据与机载激光点云数据不同，大部分建筑都有完整的立面信息。常识告诉我 们，建筑物立面通常都是垂直于水平面的。基于这一点，我们计算点云每个点的特征f，定 义如公式(4)所示。如果点位于立面，那么其法向量中的z值极有可能等于0，即特征值f 极有可能等于1.

f＝1-|n_p·e_z| (4)

式中，n_p为点p处的法向量，e_z为单位列向量(0,0,1)。

由于影像点云中建筑物外表结构比较清晰，通过特征值f只能剔除建筑物墙面点，对于 非墙面点，比如雨棚等构筑物不能有效去除。建筑物的先验知识告诉我们，这些构筑物呈小 块状零星分布，不同块之间有较大距离。因此，我们可以利用2.1.3中提到的欧式聚类算法 有效去除之。

2.2屋顶面片的分割

根据先验知识，建筑物屋顶通常都是由各种平面或者类平面组成，这为我们提供了思路。 可以使用区域增长算法分割屋顶面片。

经典区域增长算法的工作原理是基于点法线之间的角度比较，即通过设置曲率变化和法 线角度变化的阈值，合并在平滑度上接近的点，从而分割面片。因此，该算法每一次输出都 被认为是共面点。但是，这种方法在此处有一些缺点。该算法会将两平行的平面识别为同一个面片，这样会导致屋顶面片识别出错。同时，由于本方法已经去除建筑物立面点，在计算 点法向量时，如果考虑点数较多，会造成法向量计算出错。

结合这些信息，我们提出一种结合点间距的区域增长算法来分割屋顶面片。首先，设置 点间距的阈值，通过迭代的计算最近邻点的点间距，我们可以对点集聚类，并统计每类点的高 程均值。为了防止孔洞的影响，我们将高程近似的聚类合并。然后，我们以每一类点作为输入，按照经典区域增长算法的原理，识别出不同面片。反复迭代，直到遍历所有点。

2.3屋顶面片的修补

由于计算误差和照相机摄影角度的原因，影像点云往往会存在部分点与现实世界不一致 甚至部分点缺失等情况，比如屋顶平面与立面之间的过渡区域为不规则曲面、屋顶面片中有 孔洞等。同时，建筑物屋顶通常具有一些附属物，比如水箱、通风设施、太阳能热水器等。经过立面去除后，仍然会有点残留。这些情况将造成分割出的屋顶面片边界不规则甚至不完 整。通过观察大量建筑屋顶发现，建筑物屋顶通常都是由矩形、三角形、扇形组成的规则面 片，面片内的点分布无规律，但是边界点有严格的形状约束。受到这一点启发，本发明提出 了一种新的基于PCA的面片修补方法，既能使面片规则化，又能保持面片之间的拓扑关系。

该算法主要分为两步：屋顶面片及附属物的合并，和面片的修补。算法原理如下：

1)将所有点标识位设为0。如公式(9)所示，找到点数最多的面片P_max，将其标志位改为1。

P_max＝{P_i|max{P_i)，i∈n} (9)

式中：P_i表示一个面片内的点数，n表示屋顶的面片数。

2)计算所有面片与P_max的最小距离d。为了准确找出附属物点集，选取标识位为0，且 0.5<d<2.0的面片。这些面片是与P_max相关的附属物，统一称为P′。将P′的标识位改为1；

3)利用主成分分析算法将P′和P_max转换到主平面。将P′内所有点的z坐标都设为P_max的高程均值。完成投影。

4)在主平面上，基于平均点间距

进行补点。

i.找到点集中坐标的最值X_min、X_max、X_min和X_max；

过点

做垂直于x轴的直线L_i,将与L_i距离在/>

内的点都投影在该直线 上。计算相邻点间距dis，如果/>

则在两点之间***/>

个点。其中， i为1,2,3，…，直到X_min+i*d>X_max。

ii.以(0，Y_min)为起始点，重复执行ii步，直到Y_min+i*d>Y_max。

iii.从x和y坐标的最大值出发，重复执行ii～iii步。

注：为了不破坏面片边界的拓扑关系，设置阈值T₁，比如T₁＝4.0。

5)将P′从主平面还原回世界坐标系，得到屋顶面片P。

6)迭代执行步骤1)～5)，直到点的标识位都为1。

经过本算法的处理，理论上，我们可以获得规则且较完整的建筑物屋顶面片。

2.4屋顶面片的细化

经过实验发现，由于照相机摄影角度的原因，影像点云通常含有建筑物内部点。为了避 免将建筑物内部点识别为屋顶面片，我们提出一种基于二维格网的滤波方法。

原理如下：首先，找到x、y轴上的最值X_min、X_max、Y_min、Y_max；。按照公式(10) 计算每个点所属网格的行列索引号。统计每个网格内的点数，并按照高程值对点进行排序。 如果网格内存在相邻点的高程之差大于阈值Tz，比如Tz＝2.5，则舍弃高程较低的那部分点。

式中：row_i表示点i的行索引，column_i表示其列索引，width表示网格的宽度，x_i表示 点i的x坐标，y_i表示其y坐标。

本发明由于采用了以上的技术方案，与基于机载LiDAR的建筑物屋顶提取技术相比，本 发明技术方案具有以下明显优点：

1)影像点云获取技术需要的设备便宜、简单，并能够获得精确逼真的模型。倾斜成像系 统主要由高精度的POS***和多个相机组成，造价不高。与之相反，机载激光雷达设备一个 激光扫描仪就要几十万，再加上惯导***等设备，成本高昂。同时，基于图像的三维重建技术能够生成精确逼真的稠密点云模型，满足用于建筑物屋顶轮廓线的提取。

2)影像点云数据具有准确的rgb信息，能够很好的减少算法时间复杂度。与LiDAR复杂 的植被识别算法相比，我们能够通过简单的绿信比和植被指数的计算很好的识别植被，从而 减少了算法时间复杂度。

3)提出了一种全自动建筑物屋顶轮廓线提取算法，通过简单的参数设置能够自动获取屋 顶轮廓线，从而解决了传统房屋地籍测量中耗时长、人力物力投入大的问题。

附图说明

图1是本发明方法的流程图；

图2是实施例1中原始影像点云的示意图；

图3是实施例1中建筑物单体化的示意图；

图4是实施例1中屋顶面片点云的示意图；

图5是RANSAC分割算法中屋顶面片点云的示意图；

图6是RANSAC分割算法和本发明方法的对比示意图。

具体实施方式

下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中自始至终相同或 类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的 实施例是示例性的，旨在用于解释本发明，而不能理解为对本发明的限制。

在本发明的描述中，需要理解的是，术语“中心”、“纵向”、“横向”、“长度”、 “宽度”、“厚度”、“上”、“下”、“前”、“后”、“左”、“右”、“竖直”、“水 平”、“顶”、“底”“内”、“顺时针”、“逆时针”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的 装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限 制。

此外，术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、“第二”的特征可以明示 或者隐含地包括一个或者更多个该特征。在本发明的描述中，除非另有说明，“多个”的含义是两个或两个以上，除非另有明确的限定。

在本发明中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”、“固 定”等术语应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接； 可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以根据具体情况理解上述术语 在本发明中的具体含义。

在本发明中，除非另有明确的规定和限定，第一特征在第二特征之“上”或之“下”可以包括第一和第二特征直接接触，也可以包括第一和第二特征不是直接接触而是通过它们之 间的另外的特征接触。而且，第一特征在第二特征“之上”、“上方”和“上面”包括第一特征在第二特征正上方和斜上方，或仅仅表示第一特征水平高度高于第二特征。第一特征在 第二特征“之下”、“下方”和“下面”包括第一特征在第二特征正下方和斜下方，或仅仅表示第一特征水平高度小于第二特征。

以下结合附图和实施例对本发明技术方案进行具体说明。

为了测试该技术方案的正确性，选取一组由无人机采集的倾斜影像经过运动恢复结构算 法后生成的彩色点云数据进行实验。原始数据收集于中国的某个小区，一共21个txt文件，每个文件大小不一，共包含大约38593万个点，大小约为17.6G。场景内建筑物类型复杂， 包括一个幼儿园、12栋居民楼、2栋复杂结构的商业建筑、大量绿化设施和车辆，以及道路和路灯等交通设施。

由于原始影像点云点数较多，点密度较大，因此，需要对数据进行抽稀。同时，通过观 察发现，场景中地面以下存在大量分离的块状噪声点。本发明按文件迭代的执行抽稀和聚类 算法，聚类算法去除块状噪声之后，抽稀比例按1:50执行，最终合并为一个文件，大小为 359M，点数约为772万。如图2所示。

如图1所示，实施流程包括以下步骤：

步骤a，利用布料滤波算法将点云数据中的地面点与非地面点分离，得到非地面点。

布料滤波算法模拟了一个简单的物理过程。假设一块布料在空中由于重力作用掉落下来。 如果布料足够柔软，能够依附于物体表面，那么布料的最后形状就是场景的数字表面模型。 现在将点云场景翻转180度，使点云中的地面部分位于上方。不同的是，这里布料具有一定的刚度，那么布料最后形成的形状就是数字地形模型，即地面。分析布料的粒子和相应激光 点云的相互作用就可以确定布料的最后形状，从而用它将地面点和非地面点分离。基于这种 思想，该算法用到了三维计算机图形学中的布料模拟技术，即通过计算机编程来模拟布料。 这里的布料是一种由大量有固定质量但无空间大小的粒子相连形成的格网。格网的顶点是这些物理粒子，格网的边是一种服从霍克定律的“弹簧”。

该算法的具体实现过程如下：

1)将去除离群点的原始激光点云做180度翻转；

2)初始化布料格网。用户自定义的网格分辨率确定布料粒子总数；布料的初始位置被放 置在点云最高点的上方；

3)将激光点云和网格粒子投影到一个水平面。发现每个粒子的最近邻激光点，并将激光 点的高程记为交点高程值；

4)对每个粒子，计算它由于重力作用后的位置。如果它的高程值大于交点高程值，它可 以继续移动；如果小于等于交点高程值，则这个粒子的高程值被设置为交点高程值，并被标 记为不可移动；

5)对于每个格网粒子，计算由于粒子间的相互作用力造成的位置偏移；

6)迭代执行4)-5)，直到所有粒子的最大高程方差足够小或者超过了用户定义的最大 迭代次数；

7)计算粒子和点云之间的距离；

8)分离地面和非地面点。对于任何一个激光点，如果与对应网格粒子点的距离小于阈值， 则被识别为地面点；否则，为非地面点；

为了使算法适用于有较大梯度变化的地形，加入了一个后处理算法。原理如下：搜索每 个可移动粒子的四近邻粒子，如果发现有不可移动粒子，则比较两者对应激光点的高程。如 果最小高程差小于阈值，则可移动粒子继续移动，直到被设置为不可移动。迭代执行这个步骤，直到遍历所有可移动粒子。

具体实施时，可参见相关文献Wuming Z,Jianbo Q,Peng W,et al.An Easy-to-Use Airborne LiDAR Data Filtering Method Based on Cloth Simulation[J].RemoteSensing, 2016,8(6):501-。

步骤b，利用绿信比G_r和植被系数指标VI识别绿色或者近绿色点，再利用基于法向量的 变化检测算法对识别点集进行优化，进而剔除植被点。

实施例步骤b中，由于植被通常是绿色或者近绿色，因此植被的绿光波段DN(Digital Number)值一般比红光和蓝光波段DN值高。结合这两点，本发明提出使用绿信比G_r和植被 系数指标VI识别植被。绿信比和植被系数的定义如下式所示：

式中：VI和G_r分别为植被系数和绿信比；R、G、B为点云的3个颜色通道值，R′、 G′、B′分别为R、G、B所占比例，R′、G′、B′定义如下式所示：

通过设置绿信比G_r和植被系数指标VI的阈值，比如G_r>0.35，VI>0，可以有效识别出绿 色或者近绿色点集。

但是，现实世界中往往存在一些特殊颜色的建筑，比如实施例中有近绿色墙面的幼儿园。 为了防止植被点过滤时，将非植被点误剔除，需要采用基于法向量的变化检测算法进一步筛 选识别点集，确定植被点。

基于法向量的变化检测算法原理如下：

首先，估计点的法向量，这里我们采用常见的基于协方差矩阵的计算方法。如果一个区 域内的点足够多，那么它就能构造出一个表面，那么我们就可以用它来估计法线。如公式(1) 所示，我们以点云中的点p为中心，建立一个邻域N_p。其中，P为原始点集，p∈P,d(p，q)表 示两点之间的距离，r是以p为圆心的搜索半径。

N_p＝{q|q∈P，d(p，q)＜r} (1)

然后，我们构造邻域内的协方差矩阵C_p，定义如公式(2)所示。p是邻域p中所有点的中心。计算协方差矩阵的特征值，并排序λ₁<λ₂<λ₃，最小特征值λ₁对应的特征向量被视 为点p的法向量。

利用M.Pauly,Point Primitives for Interactive Modeling and Processingof 3D Geometry.Konstanz,Germany:Hartung-Gorre,2003.中提到的基于法向量的协方差矩阵 的特征值可以量化法向量的最大变化。

具体的，构造邻域

内所有点法向量的协方差矩阵，如公式(3)所示。同C_p，计算该协方差矩阵的特征值，对特征值进行排序，/>

可以定量表示该点p在高斯球面上法向量的最大变化。

通过设定阈值Tⁿ，可以有效识别出近绿色点集中的非植被点。

步骤c,利用欧式聚类算法进行连通域分析得到初步建筑物单体化结果，再利用基于高程 阈值的约束项进一步优化。

实施例步骤c中，通过步骤b获得点云，通常包含建筑物(含围墙)、车辆和电力实施， 比如电线、电杆等。通过欧式聚类算法进行连通域分析，设置点间距阈值和最小点数阈值可 以提出大部分非建筑点，实现初步建筑物聚类。

具体的，欧式聚类算法是一种基于点间距的聚类算法，算法实现步骤如下：

1)设置点间距的阈值T，所有点都被设置为未访问点；

2)设置点集P为空。以点云中任意一个未访问点为起点，标记为访问点，寻找与该点的 空间距离小于或等于T的所有点，将添加到P；

3)以点集P中任意一个未访问点为起点，重复2)，直到点集P不能增加新点。点集P为一个聚类；

4)重复执行2)-3)，直到所有点都被标记为访问点。

执行完欧式聚类算法后，实施例中还是会有部分非建筑物点被识别为建筑，比如残留地 面、类平面的绿色植被、连接建筑物的墙体等，这些误分类的存在会影响屋顶面片提取的准 确性。为了去除这些误分类，这里需要设置基于高程阈值的约束项。

具体的，找到每类点集的高程最小值Zmin，结合建筑屋顶常识，设置高程阈值，比如 z>Zmin+3，可以实现过滤大部分非建筑物点集。之后，再次执行欧式聚类算法。经过这些步 骤，能够实现高准确率的建筑物单体化。如图3建筑物单体化所示。

以单栋建筑物点云为处理单元，执行后续处理步骤。

步骤d，利用基于法向量的点特征f去除建筑物立面，得到建筑物屋顶点。

实施例步骤d中，为了去除单栋建筑物的建筑物立面，得到建筑物屋顶点云，计算点云 中每个点的点特征f。

具体的，由于建筑物立面通常都是垂直于水平面的，基于这一点，设计点特征f的定义 如公式(4)所示。如果点位于立面，其特征值f极有可能等于1.

f＝1-|n_p·e_z| (4)

式中，n_p为点p处的法向量，e_z为单位列向量(0,0,1)。

实施例中，点特征f只能剔除建筑物墙面点，对于非墙面点，比如雨棚等构筑物不能有 效去除。采用欧式聚类算法去除之。

步骤e，对单栋建筑物屋顶点云，利用基于点间距的区域增长算法实现屋顶面片分割。

实施例步骤e中，建筑物屋顶通常都是由各种平面或类平面组成，这为我们提供了思路。 可以使用区域增长算法分割屋顶面片。

具体的，经典区域增长算法的工作原理是基于点法线之间的角度比较，即通过设置曲率 变化和法线角度差异的阈值，合并在平滑度上接近的点，从而分割面片。因此，该算法每一 次输出都被认为是共面点。但是，这种方法在此处有一些缺点。该算法会将两平行的平面识别为同一个面片，这样会导致屋顶面片识别出错。同时，由于本方法已经去除建筑物立面点， 在计算点法向量时，如果考虑点数较多，会造成法向量计算出错。

结合这些信息，我们提出一种结合点间距的区域增长算法来分割屋顶面片。首先，设置点 间距的阈值，通过迭代的计算最近邻点的点间距，我们可以对点集聚类，并统计每类点的高 程均值。为了防止孔洞的影响，我们将高程近似的聚类合并。然后，我们以每一类点作为输入，按照经典区域增长算法的原理，识别出不同面片。反复迭代，直到遍历所有点。

步骤f，对步骤e的输出，利用基于PCA的面片修补算法修补孔洞，并平整面片。

实施例步骤f中，步骤e得到的屋顶面片通常边界不规则甚至不完整，甚至会有孔洞存 在。通过执行基于PCA的面片修补方法，既能修补孔洞，平整面片，又能保持面片之间的拓 扑关系。

具体的，该算法主要分为两步：屋顶面片及附属物的合并，和面片的修补。算法原理如 下：

1)将所有点标识位设为0。如公式(9)所示，找到点数最多的面片P_max，将其标志位改为1。

P_max＝{P_i|max{P_i)，_i∈n) (9)

式中：P_i表示一个面片内的点数，n表示屋顶的面片数。

2)计算所有面片与P_max的最小距离d。为了准确找出附属物点集，选取标识位为0，且 0.5<d<2.0的面片。这些面片是与P_max相关的附属物，统一称为P。将P的标识位改为1；

3)利用主成分分析算法将P和P_max转换到主平面。将P内所有点的z坐标都设为P_max的 高程均值。完成投影。

4)在主平面上，基于平均点间距

进行补点。

i.找到点集中坐标的最值X_min、X_max、Y_min和Y_max；

ii.过点

做垂直于x轴的直线L_i,将与L_i距离在/>

内的点都投影在该 直线上。计算相邻点间距dis，如果/>

则在两点之间***/>

个点。其中， i为1,2,3，…，直到X_min+i*d>X_max。

iii.以(0，Y_min)为起始点，重复执行ii步，直到Y_min+i*d＞Y_max。

iv.从x和y坐标的最大值出发，重复执行ii～iii步。

注：为了不破坏面片边界的拓扑关系，设置阈值T₁，比如T₁＝4.0。

5)将P从主平面还原回世界坐标系，得到屋顶面片P。

迭代执行步骤1)～5)，直到点的标识位都为1。

步骤g，对步骤e的输出，利用二维格网，实现屋顶面片的去重处理。

实施例步骤g中，由于照相机摄影角度和建筑物窗户的原因，步骤f得到的输出面片有 部分是建筑物内部面片。为此，这里采用基于二维格网的滤波方法剔除这部分点集。

原理如下：首先，找到x、y轴上的最值X_min、X_max、Y_min、Y_max。按照公式(10)计 算每个点所属网格的行列索引号。统计每个网格内的点数，并按照高程值对点进行排序。如 果网格内存在相邻点的高程之差大于阈值Tz，比如Tz＝2.5，则舍弃高程较低的那部分点。

式中：row_i表示点i的行索引，column_i表示其列索引，width表示网格的宽度，x_i表示点i的x坐标，y_i表示其y坐标。

本实施例中建筑物屋顶分割面片结果如图4所示。实验表明本发明的分割算法能够有效 的分割出数据集中由矩形组合形成的屋顶面结构。以下通过实验数据对分割结果做定量的统 计。

表1屋顶分割实验结构统计

其中，BN，表示建筑物的编号；SP，表示建筑物屋顶分割后得到的屋顶面片数据；RP， 表示数据集中建筑物固有的屋顶面片数目；Completeness，表示对建筑物屋顶面片分割的完 整度，即分割得到的屋顶面片数目与固有的屋顶面片数目之比；

实验结果显示，本发明所提出的屋顶面片分割方法在绝大多数建筑物中都能够得到较好 的分割结果，算法稳定性好，平均分割结果准确率为97.37％。对比实验结果和原始数据发现， f和g组实验结果出现欠分割的原因是两个屋顶面片处于同一水平面上，且相邻，因此造成无法分割；m和n组实验数据结构相似，出现欠分割的原因是两个屋顶面片之间的间隔距离 过小，导致面片合并，以及由于面片面积较小，点数过少，因此造成的分割不足；h组出现 过分割的原因是该面片中有一个较大的孔洞，并且其附近附属物过多，因此造成建筑立面去除后，孔洞面积加大，从而导致了过分割。这些实验结果足以说明本发明的算法对由凸多边 形组合形成的建筑物屋顶结构具有很好的适应性。

具体实施时，采用模块化方式实现相应***。具体包括以下模块：

模块1，用于使用布料滤波算法分离地面点，得到非地面点；

模块2，用于利用绿信比和植被系数得到初始植被点，基于法向量的变化检测进一步约 束植被点；

模块3，用于执行欧式聚类算法和高程阈值约束，从而实现建筑物单体化；

模块4，用于实现建筑物立面消除算法，和实现屋顶面片分割，获得初始的屋顶面片；

模块5，用于实现基于PCA的面片修补算法，和基于二维格网的面片去重，得到平整的 屋顶面片。

各模块具体实现可参加相应步骤，这里不再详述。

本方法提出了一个框架流程，涉及多个算法的线性组合。如表2和表3所示，整个框架 表面上统计多达18个参数，实际上，需要调节的参数只有4个，分别为点特征f的阈值Tf、区域增长算法中的平坦度阈值SmoothnessThreshold和曲率阈值CurvatureThreshold、邻域点数或半径。下面详细分析各个参数。

表2建筑物单体化技术中参数

表3屋顶面片分割技术中参数

表2中，植被过滤部分，由于算法是基于点的颜色信息来过滤绿色植物，因此Gr、Vi为常量；对于点法线的变化阈值Tⁿ，大量实验发现，对其基于平均值进行归一化处理之后，Tⁿ＝0总能取得理想的效果，因此可以认为其固定不变。布料滤波部分，使用默认参数就可以满足需求。因为通过设置单体化技术中的高程阈值T_z可以优化地面过滤结果。单体化中，结合常识建筑物高度总大于3m，因此高程阈值T_z＝2.5为常量；最小点数minP和点间距阈值dis通常都是常量值，除非建筑物规格或点密度发生大的变动，分别为minP＝10000，dis＝0.5。

表3中，立面去除部分，点特征f的阈值T1通常约等于1，这里需要通过对样例建筑物 的测试设置合适值。对于欧式聚类算法，单体化和面片分割中都有参与，点间距的阈值都采 用dis＝0.5，而最小点数minP在这里有变化，通常为500左右。面片分割中的区域增长算法， 需根据影像点云的质量对SmoothnessThreshold、CurvatureThreshold两参数做适当调整。 面片修补部分，Td主要用于防止破坏建筑的拓扑结构，根据房屋建筑常识，它一般有一个最小值，即Td＝3；点间距阈值T₁通常为常量值，T₁＝0.05。面片优化部分，格网宽通常设置 为点间距的3倍。

值得注意的是，算法中多次涉及到了法线的估计，其方法都采用了基于邻域的主成分分 析法。其中，由于影像点云点的坐标与真实情况相比会发生偏移，因此邻域设置必须较大， 通常最近邻点数设置为500或者邻域半径为0.5米左右，这里需要根据影像点云的质量作出 调整。

本方法是一系列算法的组合应用，里面涉及多个子算法。现有技术的对比文献中，都只 是涉及到部分原理重叠,比如区域增长算法、欧式聚类算法、二维格网化、主成分分析法等。

下面详细分析各个对比文件，找出具有真实对比意义的现有技术。

CN107220987A记录了一种基于主成分分析的建筑物屋顶边缘快速检测方法，其目的是利 用主成分分析算法辅助检测屋顶面片边缘，而本发明利用主成分分析算法进行屋顶面片的修 补。与本发明核心要点有明显区别。

CN108010092B记录了一种基于低空摄影测量数据建模后的城市高密度区太阳能利用潜 力评估方法。该文献中，在建筑物单体化部分，其试图通过在参数化建模工具grasshopper 里通过volvox插件根据影像点云中的彩色信息分离出建筑物点云。但是它有缺陷，其没有考 虑建筑物颜色与地面或植被相似的情况。与本发明核心要点有明显区别。

CN107944384A(协方差矩阵-区域增长)公开了一种基于三维Voronoi图的建筑物屋顶面 片分割方法，处理对象为机载LiDAR点云。该方法利用三维Voronoi图来代替kdtree建立点 邻域，避免了邻域设置时的参数设置，有效解决了构建点云数据空间拓扑关系时邻域不易控 制的问题。但是，由于处理对象为机载LiDAR点云，它采用的植被过滤算法完全区别于本本 发明，与本发明核心要点有明显区别。

CN107230251A公开了一种从倾斜成像数据和激光雷达数据创建3D城市模型的技术。该 文献涉及一种将倾斜成像数据和激光雷达数据匹配融合来创建3D城市模型的方法和混合3D 成像装置的设计。而本发明描述了一种处理彩色三维点云的算法，两者没有对比性。与本发 明核心要点有明显区别。

CN108090957A公布了一种基于BIM的测绘地形的方法。该文献介绍了一种以人机交互为 主的建立三维点云模型的方法，而本发明涉及一种全自动的建筑物屋顶面片识别算法，与本 发明核心要点有明显区别。

CN106846494A(格网法_GPU)公布了一种倾斜摄影三维建筑物模型自动单体化算法。该 算法通过将建筑物轮廓线的顶点信息传入GPU着色器，通过将轮廓线内的建筑物标记高亮， 使其实现建筑物单体化操作。该文献采用了GPU技术，其技术路线与本发明存在很大不同。 与本发明核心要点有明显区别。

CN108074232A(格网-体元)和CN108109139A相同，公开了一种基于体元分割的机载LiDAR 建筑物检测方法。该方法通过将原始机载LiDAR点云数据规则化为灰度3D体元数据集，基于 连通性和辐射特性相似性准则，将灰度3D体元数据分割并标记为若干个3D连通区域；然后， 基于建筑物屋顶和立面的特性，完成基于体元分割的机载LiDAR建筑物检测。该方法采用了 体元技术，与本发明核心要点有明显区别。

CN104809689B公开了一种基于轮廓的建筑物点云模型地图配准方法。该文献的目的是实 现点云模型与卫星底图之间的配准。与本发明核心要点有明显区别。

CN102411778B公开了一种机载激光点云与航空影像的自动匹配方法。该方法直接从机载 LiDAR点云中提取建筑物轮廓线，并以此获得配准基元的建筑物角特征，然后在航空影像近 似外方位元素的辅助下进行点云与影像间同名角特征的自动匹配；之后，采用光束法区域网 平差和循环迭代策略，实现航空影像和点云数据的整体最优配准。与本发明核心要点有明显 区别。

CN101726255B公开了一种从三维激光点云数据中提取感兴趣建筑物的方法。该方法企图 直接利用欧式聚类算法从非地面点中分离建筑物，从理论上来说这有一定局限性。因为现实 场景中往往地面上还会有植被、大型汽车等非建筑物点，因此简单的执行欧式聚类算法稳健性不高。同时，该文献的目的是通过设置建筑物边特征来检索建筑物。因此，与本发明核心 要点有明显区别。

CN107545602A(空间拓扑关系-建筑物建模)公布了一种基于LiDAR点云的空间拓扑关系 约束下的建筑物建模方法。该文献重点在建筑物屋顶几何基元的空间拓扑关系处理，在建筑 物单体化和屋顶分割方面几乎没有详细介绍。与本发明核心要点有明显区别。

CN104036544B公开了一种基于机载LiDAR数据的建筑物屋顶重建方法。该文献的重点在 如何准确获取屋顶面片矢量边界上，而建筑物屋顶面片分割部分比较简单，直接采用区域生 长算法就获得了屋顶面片，然后利用最小二乘法拟合得到屋顶面片所在平面的平面方程。与 本发明核心要点有明显区别。

CN105572687B公布了一种基于车载激光雷达点云制作建筑物数字线划图的方法。该方法 是对车载激光雷达数据进行建筑物立面检测与提取，与本发明没有对比性。

文献城区LiDAR点云的树木提取提出了一种基于区域增长与梯度分割的树木提取算法， 与本发明没有对比性。

CN106970375A公布了一种机载激光雷达中自动提取建筑物信息的方法。该方法本质上利 用平面拟合方法检测出所有建筑物点云平面，然后利用三维形态学腐蚀运算去除部分误分类 点云，从而得到准确的建筑物点云。该方法在机载LiDAR点云中适用，然而影像点云不同于 激光点云，它不仅点云密度大，而且在包含屋顶点云的同时，还包含完整的建筑物立面点。因此，如果不采取其它措施，直接进行平面拟合方法，无法通过少量的迭代检测出建筑物点 云。因此，与本发明核心要点有明显区别。

基于局部约束的建筑物屋顶点云平面分割方法介绍了一种加入点云点法向量约束的 RANSAC算法，它在一定程度上解决了传统RANSAC算法在同一平面上不同屋顶面的分离问题。 本发明的面片分割算法采用了区域增长算法，与本发明核心要点有明显区别。

期刊论文《顾及建筑物屋顶结构的改进RANSAC点云分割算法》利用坡度和高差的三角形 区域生长方法，对复杂建筑物的不同结构层次进行分解，然后提出一种浮动一致集阈值的 RANSAC算法检测建筑物屋顶平面，有一定适用性。本发明的面片分割算法采用了区域增长算法，与本发明核心要点有明显区别。

研究生毕业论文《机载LiDAR点云数据滤波及建筑物点群分割研究》改进了随机霍夫变 换对建筑物屋顶面面片所含点云的提取方法。整个算法的假设是基于建筑物模型是由多个平 面构成，这种假设有一定局限性，因为真实场景中，往往有部分建筑物是由曲面构成的。与本发明核心要点有明显区别。

CN109242855A公开了基于多分辨率三维统计信息的屋顶分割方法、***及设备。该方法 通过对影像点云提取不同分辨率的三维点云特征统计信息，然后直接利用全局能量优化对三 维点云场景进行语义分类，得到建筑物的点云。其整个技术路线不同于本发明采用的方法。

CN108171720A公开了一种基于几何统计信息的倾斜摄影模型对象边界探测方法。该方法 通过计算建筑物的二维最小外包矩形，使它与模型求交，选择属于当前对象所在区域的点云 数据来得到单体化的建筑物。该方法过程复杂，计算量大。与本发明采用技术路线有很大差别。

CN108898144A公开了一种建筑物损毁状态检测方法。该方法直接以建筑物点云数据为输 入对象，不涉及建筑物分割及屋顶面检测等，没有对比性。

CN109461207A公开了一种点云数据建筑物单体化方法及装置。该方法在对体素进行分类 处理前，先对体素中的三维点进行了一个预处理，将三维点重新做了划分，使得每个体素所 包含的三维点更为相似。该方法用到了全局能量优化函数，计算量偏大。与本发明采用技术路线有很大差别。

CN109754020A公开了一种融合多层级渐进策略和非监督学习得地面点云提取方法。该方 法是一种地面点云提取算法，与本发明核心要点有明显区别。

CN109859315A公开了一种三维影像的地表植被分离方法。该方法是一种软件操作方法， 与本发明没有对比性。

CN109870106A公开了一种基于无人机图片的建筑物体积测量方法。该方法通过无人机图 片的密集匹配获取建筑物影像点云，然后对Delaunay三角化后的模型利用积分法计算建筑物 的体积。与本发明核心要点有明显区别。

CN105844629B(RGA-DP-RANSAC)公开了一种大场景城市建筑物立面点云自动分割方法。 该方法在提取机载LiDAR建筑物屋顶点云数据的方法为：采用渐进不规则三角网加密实现地 面点和非地面点的分离，然后以地面点高程为基准，滤除高程差小于2.0m地物点，然后采用 RANSAC面片检测算法分割建筑物屋顶点云。该方法相对简单，对拥有大量立面信息的影像点云不适用。与本发明核心要点有明显区别。

CN106097311A(渐进形态学滤波-区域增长-最小回环检测)公开了一种机载激光雷达数 据的建筑物三维重建方法。建筑物点云提取中，由于点云回波次数在影像点云中没有，因此 这里不适用。其次，它采用了传统区域增长算法来识别建筑物点云，这有一定局限性。建筑物屋顶分割部分，在采取基于点云空间分布特性的聚类增长分割算法和RANSAC算法检测出屋 顶平面之后，再利用面片法向量夹角约束和面片距离约束进行优化，合并相似面片。与本发 明核心要点有明显区别。

CN106600680A公开了一种建筑物体框模型批量精细三维建模方法。该方法是一套成熟的 生产与质检作业流程。与本发明核心要点有明显区别。

CN107644452A(机载LiDAR-屋顶面片分割)公开了一种机载LiDAR点云屋顶面片分割方 法及***。该方法与本发明的区别点在于：

1)以3D Voronoi图为机载LiDAR点云数据构建点云邻域***，包括为机载LiDAR点云 数据中任意一点构建基于3D Voronoi邻域的邻接关系。本发明采用kdtree。

2)该文献利用LiDAR点云的多次回波信息来区分建筑物点与植被点。本发明处理数据为 影像点云，没有回波信息，采用了其它方法。

CN104484668B公开了一种无人机多重叠遥感影像的建筑物轮廓线提取方法。该方法与本 发明的区别点在于：

1)采用颜色不变量对地面点中的植被滤除，没有考虑建筑物有颜色相近的情况。而本发 明在采用颜色过滤之后，加入了基于法向量的变化检测算法，保证了颜色过滤不会造成建筑 物点的损失。

2)利用基于平面拟合区域增长法得到建筑物点云，其本质上是基于假设建筑物是由平面 组成的。而本发明在过滤植被和地面之后，直接采用欧式聚类算法通过设置点间距阈值和最 小聚类点数来得到建筑物。其不基于这个假设，保证了可以识别任意结构的建筑。

3)立面分割部分，该文献采用了面片法向量的方向来分离建筑物立面。本发明虽然采用 点的法向量。

CN106023312B(立面去除-平面拟合-区域增长)公开了一种基于航空LiDAR数据的三维 建筑物模型自动重建方法。该方法与本发明的区别点在于：

1)抓住机载LiDAR点云中建筑物立面点密度小于屋顶点的方法来提出墙面点；

2)屋顶层重采样部分使用了支持向量机算法，增大了计算量。

3)“种子区域选取-屋顶面片生长-面片平整优化”部分，采用了利用点曲率阈值初步得 到初始平面参数，然后利用距离阈值、距离标准差来实现屋顶面片生长。

综上，本发明相对于现有技术的区别点在于：

在建筑物识别方面，采用了排除法。即将地物分为3类：建筑、地面、植被，剔除地面点、植被点，剩余就是建筑点。而参考文献中的很多方法是基于“建筑物是由平面面片组成”的假设，该假设有一定局限性，不能提取出非平面组成的建筑物。而本发明不存在这个问题。

在屋顶面片分割方面，这里在提取出原始屋顶面片点云之后，采用了基于PCA的面片修 补方法，不仅解决了原始影像点云中存在孔洞的问题，还对面片进行了平整处理。而参考文 献中，没有涉及对影像点云中孔洞的处理。

具体的，

1)地面过滤算法中，采用了最新的布料滤波算法。它具有算法稳健、参数设置简单，原 理容易理解等优点，同时，它场景适用性强，不仅能够处理平坦的地面，对有一定坡度的地 面同样适用。而参考文献中绝大部分用的是渐进形态滤波算法。

2)植被过滤算法中，发明了一种新的方法。首先，基于点的颜色信息将点云分为绿色点 和非绿色点。然后，计算通过协方差矩阵中的特征值度量绿色点中法向量的变化，以此来判 断哪些点属于建筑物点。在这前后两步的双重作用下，能够保证在剔除绿色植被的同时，保证不会误删除建筑物点。参考文献中，也有类似采用颜色信息来过滤植被的方法，但是它们 没有考虑建筑物颜色类似于植被的情况，因此算法有一定的局限性。

3)屋顶面片修补算法中，提出了一种基于PCA投影转换的面片修补方法，既能处理屋顶 的附属物，使面片平整，又能修复孔洞，保持面片之间的拓扑关系。而参考文献中并没有考 虑到对面片内孔洞的修复问题。

4)屋顶面片的细化部分中，采用了一种基于二维格网的识别重叠面片的方法。在机载 LiDAR点云或者无人机倾斜影像匹配点云中，由于激光雷达或者照相机的倾斜角度问题，点 云中往往包含一些建筑物内部点甚至建筑物内部面片点集，这里通过对格网点集进行高程聚 类，剔除点集中高程较低者的方式，完成建筑物内部点的删除，从而提高了建筑物屋顶面片 识别的准确性。

下面验证本发明方法与现有技术实例的比较。

1.植被过滤算法方面：

这里选取CN104484668B中提到的基于颜色不变量的植物滤波方法和本发明方法进行比 较，测试数据集用具体实施例中的数据集。

具体的，设点云中每个点的坐标为(x,y,z)，颜色三个通道为(R,G,B)，颜色不变量的对 于植被的阈值为Tg，利用绿色和蓝色通道定义的颜色不变量公式为：

其中，Ig(x,y,z)、I_b(x,y,z)表示点云在(x,y,z)点的绿色和蓝色分量值。ψ_g(x,y,z)表示 在(x,y,z)点的颜色不变量。当ψ_g＜T_g时，表示该点为植被点；否则，为非植被点。

表3植被过滤算法比较

方法	Pts
		本发明方法	216335
颜色不变理论	225893

其中，Pts表示过滤算法获得的植被点。由于正确植被点数无法量化，这里我们只能依据 肉眼视觉结果进行结果比较。对比实验表明，在肉眼观察识别出的正确植被点数相同的情况 下，由于本发明的方法加入了基于法向量的变化检测算法，很好的保护了近绿色建筑点集，因此使本方法的结果点数小于颜色不变理论。综上，本发明的算法稳定性较好并且由于现有 的基于颜色信息的植被过滤算法。

2、屋顶面片分割算法方面：

这里我们选取RANSAC点云分割算法与本实施例提出的算法进行对比。分割效果如图5所 示。结果如下：

表4屋顶分割实验结构统计

为了更好地对两种算法进行比较，这里选取了g、m、n三栋建筑进行单独分析。如图6 所示，为g、m、n三栋建筑物的屋顶面片提取结果。其中，(a-c)为RANSAC方法，(d-f)为 本发明方法。结果显示，RANSAC点云分割算法对面片点数敏感，不能有效提取出小面片。因此，本发明的算法稳定性较好且能够较好的保留屋顶面的细节结构。

在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体 示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或 者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不一定指的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任 何的一个或多个实施例或示例中以合适的方式结合。

尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的， 不能理解为对本发明的限制，本领域的普通技术人员在不脱离本发明的原理和宗旨的情况下 在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (9)

1.一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，包括如下步骤：

1)地面过滤和植被过滤；去除影像点云中的地面和植被；

2)建筑物单体化；剔除非建筑物点，并实现建筑物聚类；

3)去除建筑物立面；

4)屋顶面片的分割、修补和细化；

其中，所述屋顶面片的分割是指，对单栋建筑物屋顶点云，利用基于点间距的区域增长算法实现屋顶面片分割；所述屋顶面片的修补是指，对屋顶面片分割后的输出，利用基于PCA的面片修补算法修补孔洞，并平整面片；所述屋顶面片的细化是指，对屋顶面片分割后的输出，利用二维格网，实现屋顶面片的去重处理；

其中，所述基于PCA的面片修补算法的具体步骤如下：

1)将所有点标识位设为0；如下述公式(9)所示，找到点数最多的面片P_max，将其标识位改为1；

P_max＝{P_i|max{P_i}，i∈n} (9)

式中：P_i表示一个面片内的点数，n表示屋顶的面片数；

2)计算所有面片与P_max的最小距离d；为了准确找出附属物点集，选取标识位为0，且0.5<d<2.0的面片；这些面片是与P_max相关的附属物，统一称为P′；将P′的标识位改为1；

3)利用主成分分析算法将P′和P_max转换到主平面；将P′内所有点的z坐标都设为P_max的高程均值；完成投影；

4)在主平面上，基于平均点间距

进行补点；

5)将P′从主平面还原回世界坐标系，得到屋顶面片P；

6)迭代执行所述基于PCA的面片修补算法的步骤1)～5)，直到点的标识位都为1。

2.根据权利要求1所述的一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，所述地面过滤，利用布料滤波算法将点云数据中的地面点与非地面点分离，得到非地面点。

3.根据权利要求1所述的一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，所述地面过滤的具体实现过程如下：

1)将去除离群点的原始激光点云做180度翻转；

2)初始化布料格网；用定义的网格分辨率确定布料粒子总数；布料的初始位置被放置在点云最高点的上方；

3)将激光点云和网格粒子投影到一个水平面，发现每个粒子的最近邻激光点，并将激光点的高程记为交点高程值；

4)对每个粒子，计算它由于重力作用后的位置；如果它的高程值大于交点高程值，它可以继续移动；如果小于等于交点高程值，则这个粒子的高程值被设置为交点高程值，并被标记为不可移动；

5)对于每个格网粒子，计算由于粒子间的相互作用力造成的位置偏移；

6)迭代执行上述步骤4)至5)，直到所有粒子的最大高程方差足够小或者超过了定义的最大迭代次数；

7)计算粒子和点云之间的距离；

8)分离地面和非地面点；对于任何一个激光点，如果与对应网格粒子点的距离小于阈值，则被识别为地面点；否则，为非地面点。

4.根据权利要求3所述的一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，所述地面过滤的具体实现过程还包括一后处理步骤：搜索每个可移动粒子的四近邻粒子，如果发现有不可移动粒子，则比较两者对应激光点的高程；如果最小高程差小于阈值，则可移动粒子继续移动，直到被设置为不可移动；迭代执行这个步骤，直到遍历所有可移动粒子。

5.根据权利要求1所述的一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，所述植被过滤利用绿信比G_r和植被系数指标VI识别绿色或者近绿色点，再利用基于法向量的变化检测算法对识别点集进行优化，进而剔除植被点_。

6.根据权利要求5所述的一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，所述绿信比即绿色波段的DN值占RGB的DN值之和的比值；绿信比和植被系数的定义如下式所示：

式中：VI和G_r分别为植被系数和绿信比；R、G、B为点云的3个颜色通道值，R′、G′、B′分别为R、G、B所占比例，R′、G′、B′定义如下式所示：

7.根据权利要求1所述的一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，所述建筑物单体化，利用欧式聚类算法进行连通域分析得到初步建筑物单体化结果，再利用基于高程阈值的约束项进一步优化。

8.根据权利要求1所述的一种从大规模影像密集匹配点云分割建筑物屋顶面片的方法，其特征在于，建筑物立面的去除，利用基于法向量的点特征f去除建筑物立面，得到建筑物屋顶点云；点特征f的定义如下公式所示：

f＝1-|n_p·e_z|

式中，n_p为点p处的法向量，e_z为单位列向量(0,0,1)。

9.一种从大规模影像密集匹配点云分割建筑物屋顶面片的***，其特征在于，采用如权利要求1至8中任一项所述的方法，包括：

第一模块，用于使用布料滤波算法分离地面点，得到非地面点；

第二模块，用于利用绿信比和植被系数得到初始植被点，基于法向量的变化检测进一步约束植被点；

第三模块，用于执行欧式聚类算法和高程阈值约束，从而实现建筑物单体化；

第四模块，用于实现建筑物立面消除算法，和实现屋顶面片分割，获得初始的屋顶面片；

第五模块，用于实现基于PCA的面片修补算法，和基于二维格网的面片去重，得到平整的屋顶面片。

Images