基于物体表面正射投影的悬崖点云植被过滤方法

李彩林，陈文贺，王江妹，田鹏艳，姚吉利

(山东理工大学 建筑工程学院，山东 淄博 255049)



基于物体表面正射投影的悬崖点云植被过滤方法

李彩林，陈文贺，王江妹，田鹏艳，姚吉利

(山东理工大学 建筑工程学院，山东 淄博 255049)

悬崖陡坡是地形图重要的地貌要素，该要素在生态环境工程建设和地质灾害的防治等方面具有非常重要的作用，而悬崖植被的遮挡不利于悬崖观测和数据处理.提出一种悬崖植被过滤方法，该方法主要基于物体表面正射投影原理，首先将原始采集的悬崖三维点云转换到以悬崖表面为xy平面、垂直于悬崖面的方向作为z轴的空间直角坐标系中；然后在xy平面上划分格网，分格网拟合曲面建立局部数字地面模型(DTM)，通过设置合理的距离阈值剔除植被点；最后将过滤保留下来的悬崖岩石点经过投影反变换，映射到原始空间坐标系中.对实际悬崖点云数据进行实验分析，结果表明，基于物体表面正射投影的悬崖点云植被过滤方法是可行的、有效的.

三维激光扫描；悬崖点云；植被过滤；物体表面正射投影；DTM

0 引言

地形三维激光扫描仪具有获取数据速度快、精度高、无需接触即可获取大面积密集点云数据的一系列优点，广泛应用于工业设备形变测量、山体滑坡、实体变形监测等方面，其中悬崖、滑坡体的形变监测是其中相对活跃的一个研究领域[1-4].针对悬崖、滑坡体形变监测，目前常用的手段是利用三维激光扫描仪获取悬崖、滑坡体表面三维点云数据.在获取地势陡峭、表面复杂、植被覆盖广的滑坡体表面数据时，会出现非地形点比重大、点云密度分布极不均匀的情况[5]，这种情况加大了点云数据处理的工作量，同时增加了工作难度.过滤即保留地面点，去除非地面点.通过点云过滤，得到相应的数字产品.目前在点云滤波方面已经有了大量的研究，也取得不少滤波方面的研究成果.李等[5]研究了基于密集点云的双重滤波方法，主要运用到了局部区域增长法；Liu等[6]对居民区的Lidar点云数据进行树冠过滤提取，经过区域增长算法和多次反射方法取得较理想的效果；Sui等[7]提出基于改进的数学形态学的LiDAR点云数据滤波算法；文献[8]提出一种网格分块与移动最小二乘法的植被分类算法；Brodu等[9-10]提出利用不同尺度下点云的维度信息进行分类的算法；文献[11]提出一种点云数据全自动滤波算法；Kobler等[12]提出一种迭代过滤的算法.以上算法多适用于机载雷达数据，同时过滤选择的基准面是地面，所以以高程突变信息进行植被、建筑物滤波，而悬崖的坡面很陡，植被沿坡面生长，且与岩体裂隙交错，植被激光点云数据沿高程不具有突变信息，因此，对地面激光扫描仪获得的悬崖点云数据过滤并不适用.针对悬崖植被点云过滤存在的问题，笔者在已有研究基础上，提出了一种基于物体表面正射投影的悬崖点云植被过滤方法(Orthographic Projection based on Object Surface，简称OPOS)，OPOS的基本技术思想是将悬崖点投影到以悬崖表面作为投影面、悬崖表面法向作为z轴的三维空间直角坐标系中，然后在xy平面上分块进行曲面拟合建立局部DTM模型，采用距离阈值法进行悬崖上植被及附属物体的过滤，最后将悬崖岩石点经过投影反变换，映射到原空间坐标系.通过实验，并结合具体工程实例，表明OPOS方法是有效的、可行的.

1 基于物体表面正射投影原理

激光点云中每个数据点都具有空间三维坐标(x,y,z)，点云数据处理的实质就是依据数据点间的空间几何关系，在相应规则下实现空间变换与判断.对于较大区域悬崖点云数据，空间三维坐标可能会包含非常多的数据点，既影响数据处理效率，也会使数据处理结果的可靠性得不到保障，为计算方便，有必要将陡峭地形转化为缓和地形.因此，本文提出一种物体表面正射投影变换的方法.

基于物体表面正射投影原理是在原始点云数据中取悬崖上距离相隔较大的三个点建立三维空间直角坐标系，并依据六参数数学模型和罗德里格矩阵求出正射投影变换参数，悬崖上两个空间点的坐标转换的关系系数，并由此进行坐标转换，完成悬崖点云数据的正射投影.

1.1 用空间三点建立平面直角坐标系

A、B、C为悬崖原始数据上所取的三点，如图1所示，并按逆时针编号,其坐标系经常是国家测量坐标系，A、B、C投影后对应点为a、b、c.投影后的三维坐标系选择是a为原点，a到b的方向为x轴方向，建立右手三维空间直角坐标系o-xyz(物体表面正射坐标)，如图2所示.此时，3个点的z坐标全为0，其坐标表达为：a(0,0,0)，b(S1,0,0)，c(S3cosθ,S3sinθ,0).式中的Si(i=1,2,3)是由空间坐标计算的空间距离，θ为Δabc的顶角.物体表面正射投影关系就是两个三维坐标系之间的关系，也就是根据三个公共点，求解6个三维坐标转换参数.

图1 原始点云数据取点样图

1.2 物体表面正射投影关系的确定

根据参考文献[13-14]，物体表面正射投影关系为

(1)

图2 空间三点建立直角坐标系

式中：R=(I+S)(I-S)-1；S为反对称矩阵；I为3阶单位矩阵.S的表达形式为

(2)

式(1)中，ΔX、ΔY、ΔZ为平移参数；R为旋转矩阵，由三个旋转角的函数构成，也可由反对称矩阵S的三个元素a、b、c构成，根据文献[14]，其计算公式为

(3)

式中,

u2=XB-XA+S1,

v2=YB-YA,w2=ZB-ZA,

u3=ZC-ZA+S3cosθ,

v3=YC-YA+S3sinθ,

ΔH=u3v2w2-u2v3w2,w3=0.

(4)

在这里,ΔX、ΔY、ΔZ和a、b、c称为物体表面正射投影参数，简称投影参数.

1.3 物体表面正射投影反变换

由投影坐标系的坐标计算测量坐标系的坐标的过程称为物体表面正射投影反变换，公式(4)反变换为

(5)

2 植被过滤

点云过滤的目的，就是获取悬崖岩石表面点，并过滤掉植被等非悬崖岩石表面点.由于区域比较大，为了实现对点云数据的快速处理，有必要进行分块处理并建立对应的索引，同时为进一步精确处理，对每一分块区域进一步细分.因为悬崖表面在局部细节上表现为空间二次曲面的特性，由局部地形点拟合的曲面能较好地反映出该区域的地形特征，则选取每个细分区域的地形最低点拟合局部曲面建立DTM，进而通过计算每个悬崖点云数据点到建立的DTM曲面的距离，通过经验阈值，过滤植被点云.

(1)在物体表面正射投影坐标系内，将平面分成2.5米方格，采用二十五叉树的方法对各方格建立相应的索引，如图3所示.

图3 空间索引原理图

(2)对每一方格，进行0.5米分区处理，在每一个分区中选取小区的最低点，每一个分块由25个最低点拟合曲面，建立一个DTM模型.

(3)对所有的方格建立DTM模型，形成分块的DTM数据.

(4)求取每个网格内悬崖点云数据点到建立的DTM曲面的距离，根据经验阈值(植被点云的过滤效果与所选距离阈值密切相关，当选取过大时，会导致部分植被点云未被过滤干净；当选取过小时，会导致少量地面点云被误过滤掉.本文在确定距离阈值时，结合当前悬崖点云数据的点位精度和曲面拟合精度，通过多组实验验证，距离阈值一般采用3倍点位精度，可以得到比较好的过滤效果)，将允许范围内的点云数据保留，删除允许范围外的点云数据.

通过上述步骤即可以获取悬崖表面的岩石点云数据.

3 实验及分析

3.1 数据来源

为了验证本文方法的可行性，在IDL环境下实现了该算法.本实验所用数据来自泸溪县城东南沅江江岸的悬崖石刻工程，测绘悬崖立面图的目的是进行石刻工程论证、设计、调查悬崖岩石产状，以开发本地区的旅游资源，崖壁南北距离约1.5公里，如图4.在自主开发的软件EEXLT环境中进行点云数据预处理，主要包括离群点去除、点云去噪、数据精简，相关预处理方法也比较常规，离群点去除是通过计算每个点到其k近邻点的平均距离，然后根据平均距离值大小判断离群点；点云去噪是运用基于曲率特征的方法去除噪声点；数据精简的方法则运用了基于法向矢量距离累加的方法.预处理后，将悬崖点云数据分成十二块区域分别使用正射投影法进行植被过滤，取得较好效果.

图4 悬崖展示图

3.2 植被过滤

本次实验中使用正射投影坐标转换的原理将原始点云数据(见图5(a))坐标转换成正射投影面坐标系中的坐标[13-14]，如图5(b).为区分哪些是悬崖表面岩石点，哪些是非岩石点，实验时使用自主开发的EEXLT软件中专门浏览地形的圆柱浏览器CYLVIEWER，将悬崖上植被和探照灯等地物信息暴露，如图5(c)，红色为悬崖表面岩石点云数据，蓝色为悬崖表面上非岩石点云数据.

再对正射投影后的点云数据进行金字塔分区坐标排序，金字塔坐标分区即对点云数据以不同的边长分区，可进行多重滤波提升植被过滤效果.建立平面分辨率为2.5 m、0.5 m、0.1 m的金字塔索引数据库，选取每分区内的最低点建立DTM，过滤掉悬崖上的植被等非岩石点，见图6(a).将处理后的点云数据进行正射投影的坐标反变换，获取悬崖表面岩石点云数据，见图6(b).

3.3 EEXLT与LasaStudio软件的处理效果比较

比较中使用同一点云数据，见图5(a).打开LasaStudio.exe，新建工程，填写工程属性，导入航空点云数据，选择控高度制网格大小1 m，迭代次数2次，阈值18%，将悬崖原始点云数据直接进行过滤，结果见图7.图8为导入经正射投影变换后的点云数据效果图.经过自主开发的EEXLT软件处理的点云数据图为图6(b)所示.

图5 EEXLT点云处理

图6 植被过滤及点云反变换后结果图

由7、图8分析可以得出，有无运用正射投影方法进行点云植被过滤，其结果相差较大，不经过正射投影，运用LasaStudio软件进行点云植被过滤的效果不是很理想.同时,由图6(b)、图8结果可以看出，运用本文的点云植被过滤方法相较于LasaStudio基于双边滤波的过滤方法，效果相对较好.由此，针对悬崖地形的植被过滤，本文探究的基于物体表面正射投影的悬崖点云植被过滤方法，滤除非地形点、保留地形点的精度及效果更好.

3.4 算法分析比较

为进一步验证本文方法的有效性，将本文方法与目前较为经典的一种迭代过滤算法进行比较

分析[12]，采用两组悬崖点云数据进行植被过滤实验，原始数据分别如图9(a)和图10(a)所示.

利用迭代过滤方法和本文方法对宝剑山点云进行植被过滤后，地形点和植被点叠加显示效果图分别如图9(b)、(c)所示， 图中的白色部分表示悬崖表面地形点，灰色部分表示植被点.图11为区域放大图， 可以看出， 两种方法的植被分类结果差别较大，结合宝剑山悬崖植被实际分布情况，表明经过迭代过滤方法得到的植被大致呈条状分布，与宝剑山悬崖植被实际分布情况不符.通过图11(a)放大图可以进一步看出，迭代过滤方法将很多的植被点归类为地形点，同时将很多地形点错误地归类为植被点，达不到悬崖植被过滤的要求，同时也验证了迭代过滤算法一般适合于过滤基准面为地面的情况，不适合悬崖植被过滤的情况.而利用笔者提出的算法可以准确地对悬崖表面地形点和植被点进行分类，较好地实现悬崖植被的过滤，见图9(c)和图11(b)所示.

图7 LasaStudio未经正射投影点云处理效果图

图8 LasaStudio经正射投影点云数据处理效果图

图9 宝剑山点云植被过滤效果图比较

图10 白鸽山点云植被过滤效果图比较

利用迭代过滤方法和本文方法对第二组数据(白鸽山点云数据)进行植被过滤，过滤后的结果分别如图12所示.可以看出，分类结果差别依旧很大，图12(a)植被占据了大部分，不符合白鸽山植被的实际分布情况.因此，结合植被分布实际情况，得出本文算法的植被过滤结果比较理想.

图11 宝剑山植被过滤局部放大图比较

图12 白鸽山植被过滤局部放大图比较

综合两组悬崖植被点云数据分析可以得出，迭代过滤算法不适用于悬崖植被点云过滤，因为其依据地面为过滤基准面；而本文的算法由于坐标转换则较好地过滤了悬崖植被点云.

4 结论

笔者基于原始悬崖点云数据，首先选取相隔距离较大的三个点构成三维空间直角坐标系，并将悬崖坐标进行转换，完成悬崖点云数据的正射投影；然后，进行金字塔坐标分区，在每一区内选择高程最低点作为该区的DTM点，并进行DTM的建立，实现植被点云过滤的目的；最后通过反变换将植被过滤后的悬崖点云数据坐标还原，最终实现悬崖点云数据的过滤.

本文方法的一个重点是对悬崖点云进行正射投影和进行金字塔分区.对悬崖点云进行正射投影，使之能够像平地般进行点云拟合过滤，而金字塔分区则得到多重滤波的结果，使植被过滤的效果得到大大的提高.

本文方法虽可实现对陡崖的植被等非地面点的过滤，但需依据实际地形条件选择合适的分辨率(分区的大小)来建立DTM.若分辨率过小，则部分地面点可能会被过滤掉；若分辨率过大，则非地面点不能得到有效的过滤，达不到过滤的效果.

[1] PROKOP A,PANHOLZER I-L.Assessing the capability of terrestrial laser scanning for monitoring slow moving landslides[J].Natural hazards and earth system sciences,2009,9:1921-1928.

[2] SUI L, WANG X,ZHAO D,et al.Application of 3D laser scanner for monitoring of landslide hazards[C]//Isper Symposiumr.2008：277-282.

[3] OPPIKOFER T,JABOYEDOFF M,BLIKRA I, et al.Characterization and monitoring of the aknes rockslide using terrestrial laser scanning[J]．Natural hazards and earth system sciences,2009,9:1003-1019．

[4] 褚宏亮，殷跃平，曹峰,等.大型崩滑灾害变形三维激光扫描监测技术研究[J].水文地质工程地质，2015,42(3):128-134.

[5] 李永强,吴立新.陡坡密集点云双重滤波方法研究[J].地理与地理信息科学,2011,27(1):7-10.

[6] LIU J P,SHEN J,ZHAO R,et al.Extraction of individual tree crowns from airborne LiDAR data in human settlements[J]. Mathematical and computer modelling,2011,10:1-12.

[7] SUI Lichun, ZHANG Yibin, LIU Yan, et al. Filtering of airborn LIDAR point cloud data based on the adaptive mathematical morphology[J].Acta Geodaetica et Cartographica Sinica,2010,39(4):390-396.

[8] LIU Changjun, LI Jian, DING Liuqiang, et al. A point clouds filtering algorithm based on grid partition and moving least squares[J].Procedia engineering,2012(3):476-483.

[9] BRODU N, LAGUE D. 3D terrestrial LiDAR data classification of complex natural scenes using multi-scale dimensionality criterion:applications in geomorphology[J].ISPRS journal of photogrammetry and remote sensing,2012(68):121-134.

[10]LAGUE D,BRODU N,LEROUX J. Accurate 3D comparison of complex topography with terrestrial laser scanner:application to the rangitikei Canyon (NZ)[J].ISPRS journal of photogrammetry and remote sensing,2013(82):10-26.

[11]李健，方宏远，崔雅博,等. LIDAR点云数据全自动滤波算法研究[J].郑州大学学报(工学版),2016,37(1):92-96.

[12]KOBLER A, PFEIFER N,PETER O,et al.Repetitive interpolation: a robust algorithm for DTM generation from aerial laser scanner data in forested terrain[J].Remote sensing of environment,2007,108:9-23.

[13]姚吉利.三维坐标转换的静态滤波模型[J].武汉大学学报(信息科学版),2005,30(9):825-828.

[14]姚吉利,韩保民,杨元喜.罗德里格矩阵在三维坐标转换严密解算中的应用[J].武汉大学学报(信息科学版),2006,31(12):1094-1096.

The Filtration Method of Cliff Point Cloud Vegetation Based on the Surface Orthographic Projection

LI Cailin, CHEN Wenhe, WANG Jiangmei, TIAN Pengyan, YAO Jili

(School of Architecture Engineering, Shandong University of Technology, Zibo 255049, China)

Cliff and steep slope are important landscape elements of topographic map, and these elements play a very important role in the construction of the ecological environment and prevention of geological disasters, etc. However, it is unfavorable to observe and process data because of vegetation occlusion on cliff. In this paper, we present a cliff vegetation filtration method based on the principle of surface orthographic projection. Firstly, transform the original three dimensional point cloud of cliff to the spatial cartesian coordinate system, whosexyplane is the cliff face andz-axis is perpendicular to the direction of the cliff surface. Then the grid on thexyplane is divided to establish local grid Digital Terrain Model (DTM) by fitting surface, and the vegetation points can be extracted through setting a reasonable distance threshold. Finally, after inverse projection transformation, cliff rocky points preserved are mapped to the original spatial coordinate system. The experimental analysis using actual cliff point cloud data shows that the cliff point cloud vegetation filtering method based on the surface orthographic projection is feasible and effective.

3D laser scanning; cliff point cloud; vegetation filtration; surface orthographic projection; DTM

2016-07-20；

2016-09-19

国家自然科学基金资助项目(41601496)；山东省高等学校科技计划项目(J15LN32)；山东省自然科学基金资助项目(ZR2014DL001)；四川省地理国情监测工程技术研究中心资助项目(GC201512)；山东理工大学青年教师发展支持计划项目(4072-114016)

李彩林(1985—)，男，安徽安庆人，山东理工大学讲师，博士，主要研究方向为数字摄影测量与计算机视觉、数字图像处理，E-mail:licailin@whu.edu.cn．

1671-6833(2016)06-0077-06

P237

A

10.13705/j.issn.1671-6833.2016.06.005