CN107274417B - 一种基于机载激光点云聚集关系的单木分割方法 - Google Patents

Abstract

本发明属于机载激光雷达森林遥感应用领域，具体涉及一种基于机载激光点云聚集关系的单木分割方法。本发明利用机载激光雷达获取森林植被冠层三维激光点云数据，经历林木判别、局部最大值法初始分割、删除错误的分割单树、循环聚类和最近邻点法处理，建立基于原始三维激光点云数据的从森林中分割出单棵树的单木分割方法。相比CHM单木分割方法，能避免信息丢失、减少错误分割和不适当的形状分割，分割效果优异，并且计算时间也得到了大幅度的缩减，能快速、准确的对森林激光点云数据分割出属于特定单棵树的回波激光点来，有助于单树的参数反演、单树的三维结构特征重建以及对森林植被参数反演的验证。

Description

一种基于机载激光点云聚集关系的单木分割方法

技术领域

本发明属于机载激光雷达森林遥感应用领域，具体涉及一种基于机载激光点云聚集关系的单木分割方法。

背景技术

机载激光雷达(airborne light detection and ranging，A-LiDAR)是通过搭载在航空飞机上的传感器所发出的激光来测定传感器与目标物之间距离的一种有效且主动的遥感技术，不仅可以精确、快速的获取地球表面几何信息(如三维地理坐标)、回波地物的辐射信息，而且可以应用于大尺度上，再加上其具有很强的穿透森林植被的能力,为森林结构参数的提取以及森林生态系统功能参数的反演提供了全新的技术手段。

在林业应用方面，利用机载激光雷达数据进行研究大致包含森林参数反演(树高、冠幅、叶面积指数等)、单木分割、林木判别以及森林3D模型重建。林木判别是从原始激光点云中识别出植被回波点的方法；而单木分割算法主要用来解决的问题是在林木判别后得到的植被回波点中分割出属于特定单棵树木的点云。在实际野外测量中，比如LAI值、树高、冠幅等实测参数往往都是针对单棵树木进行人工测量的，而在利用A-LiDAR或者其他遥感数据反演森林参数过程中，基本上都是以这些单木实测数据源作为验证数据，而且单木分割的位置精度、点云数据的密度等也影响着森林参数反演的精度,因而单木分割这个步骤往往是不能忽视的。

机载激光雷达数据不仅提供了高精度的地理位置数据，并且采用“自上而下”的方式扫描植被冠层，得到的回波数据中包含描述森林冠层的垂直分布特性的信息，使得在森林样方下的单木分割成为可能。

目前国际上存在的单木分割方法中，最常用的方法是根据植被冠层的冠层高度模型(Canopy Height Model，CHM)进行单木分割，大致过程是它首先通过建立的数字表面模型(Digital Surface Model，DSM)与数字高程模型(Digital Elevation Model，DEM)相减得到CHM，然后从CHM中定位出局部最大值，并自定义一个范围来扩展至被识别单棵树的冠幅，以此冠幅为基础提取出该棵树。使用CHM方法进行单木分割时，由于建立的CHM常常会使用插值点代替真实的三维空间点，会造成不可逆的信息丢失，再加上自定义的窗口大小对CHM进行处理，窗口大小与树高和冠幅大小相联系，因而在利用CHM进行分割时往往造成错误分割以及不适当的形状分割。

而在三维空间上直接对原始点云进行分割发展的越来越好，比如H.Lee等人(2010)，根据局部最大值的点与点之间的距离提出可适应性聚类方法直接对植被点云数据进行聚类，使其单木分割的总体精度达到了95.1％。Lu等人(2014)基于强度和激光点云数据的在空间分布的3D结构，提出了先识别树干的自低向上的方式分割阔叶林，结果表明该方法能从原始点云中检测出84％的树木来，而在检测出的树木中进行分割的正确率达到了97％。Li等人(2012)使用小光斑离散回波激光雷达点云数据，通过考虑两棵树之间的相对间隔，以及形状指数来判断激光点，结果表明该方法能从原始点云中检测出86％的树木来，正确分割率达到了94％。

发明内容

针对上述存在问题或不足，本发明提供了一种基于机载激光点云聚集关系的单木分割方法，基于原始点云数据，根据单棵树内的机载激光点云在三维空间上的聚集关系以及邻近的树与树之间的点云在三维空间上的聚集关系进行单木分割，提供了单木分割的新路径。

本发明具体技术方案如下：

步骤1、林木判别：

由于单木分割算法只是针对植被点云数据的，把回波点分类为植被回波点、地面回波点和其余回波点。

步骤2、局部最大值法初始分割：

在未编码(默认编码为0)的原始点云中寻找回波点类别属于植被的且坐标Z值最大的点，并记录该点的坐标序号，同时把该点编码设为k(k>＝1)，用半径为R(实测目标树木数据优势植被的最小冠幅)的圆标识该区域，并且把该区域内的所有植被回波点都设为编码k；然后在原始点云中去除已经编码了的点，树编码号变为k+1，以此循环下去，直到原始点云中所有植被点的编码都不为0。

步骤3、删除错误编码树：

当同一编码树的点云数量小于激光点云密度或两颗树之间的间隔小于3R，修改该树的编码为0。

步骤4、聚类：

A、使用步骤3得到编码数据，依据树编码号选取任意两颗树并确定两树之间的间隔I。

B、if判断，如果I小于maxR，maxR为实测目标树木数据优势植被的最大冠幅半径,则依次执行①、②、③、④。

①、由于在树冠边缘的不规则，高树相对矮树树高的边缘植被点与其顶点形成的垂直夹角θ是不同的；首先比较两树的高矮，以低树的中心点为原点，不考虑高程坐标Z，建立笛卡尔直角坐标系，在高树顶点与低树顶点连线上的点云作为高树的边缘植被点(如图4所示)。②、求最小高树边缘植被激光点与高树顶点形成的阈值正切值mintanθ(原理如图3),即以mintanθ为阈值tanθ。

③、记录与mintanθ比较的激光点：

为了减少聚类计算时间以及矮树边缘异值点对聚类结果的影响：记录两树顶点连线为长边长度及方向，矮树树冠幅为宽边的平行四边形所属激光点point2Peek，且该点的编码不为高树编码的点(如图4所示)。

④、根据多木判别原理进行聚类：首先判断正切值，如果步骤③中记录的激光点对应的正切值大于mintanθ，判断该点编码对应的矮树为高树的一部分，并把该矮树对应的点编码改为高树的编码，将平行四边形内所属激光点遍历(原理如图3)。

⑤、根据单木归属原理进行聚类：如果point2Peek对应的正切值小于mintanθ，那么判断该点为高树的冠幅点，修改该点的编码为高树的编码，遍历直到修改完该棵树的所有回波点(原理如图2)。

D、重复执行A、B步骤，直到遍历完所有组树。

步骤5：最近邻点法处理，把未被编码的激光点分到距离最近的已编码聚类的树中。

步骤④和步骤⑤的聚类原理

单木归属原理：在同一棵树中，假设在高度为h时,边缘植被点与顶点形成的垂直夹角为θ，边缘植被点以内的点与顶点形成的夹角为θ`；在高度低于h时，冠幅内的点与顶点形成的夹角为a、b(如图2所示)。在同一高度h，θ`大于θ角，在高度小于h时，a、b小于θ角。若已知高度大于h的点属于该棵树的冠层点，那么可以求出在高度为h时边缘点与顶点形成的最小角度θ；根据高度小于h的其余点和顶点形成的夹角a、b与θ的关系，可以判断其余点与该树的所属关系。

多木判别原理：在一组树中，与矮树顶点同一高度上，高树的边缘植被点与其顶点形成的垂直夹角为θ；高树顶点与矮树顶点形成的夹角为a；矮树其余点与高树顶点形成的夹角为b(如图3所示)。如果在矮树中的角a小于θ，且存在部分角b小于θ，那么矮树不划归为高树冠层部分；否则把矮树划归为高树。

本发明首先使用TerraScan软件进行林木判别，把原始机载激光雷达数据分类为地面回波点、植被回波点和其余回波点；然后只针对植被回波点利用局部最值法进行初始的单木判别，同时在此过程中删除错误判别的树；再利用单棵树内的点云以及邻近的树与树之间点云在三维空间上的聚集关系进行聚类处理，从而得到较为精确的分割结果；最后跟据最近邻点法把未分类的回波点划分到离其最近的单木中，从而提取出整片区域的树木来，其流程如图1所示。

本发明是通过利用机载激光雷达获取森林植被冠层三维激光点云数据，并根据这些三维点云数据在单棵树内的聚集关系以及邻近的树与树之间的聚集关系，经历林木判别、局部最大值法初始分割、删除错误的分割单树、循环聚类和最近邻点法处理，建立基于原始三维激光点云数据的从森林中分割出单棵树的单木分割方法。

综上所述，本发明相比CHM单木分割方法，能避免信息丢失、减少错误分割和不适当的形状分割，分割效果优异，并且计算时间也得到了大幅度的缩减，能快速、准确的对森林激光点云数据分割出属于特定单棵树的回波激光点来，有助于单树的参数反演(如叶面积指数、树高、生物量等)、单树的三维结构特征重建以及对森林植被参数反演的验证。

附图说明

图1为本发明的流程示意图；

图2为单木归属原理图；

图3为多木判别原理图；

图4为聚类过程中的确定两棵树的相对方向以及取值范围；

图5为森林样方林木判别结果图；

图6为本发明的操作界面图；

图7为初始识别和删除错误编码树后的结果图；

图8为单木分割最终结果图；

图9为使用实测样方内的林木棵数来验证分割精度的结果图；

图10为使用实测树高来验证分割精度的结果图。

具体实施方式

以下通过一个样方的实施并结合附图对本发明作进一步详细说明：

开发环境Microsoft Visual Studio 2010、开发平台框架.Net Framework3.5、三维可视化引擎DirectX 11、编程语言C#。

步骤1，以内蒙古呼伦贝尔市西北部，额尔古纳市东南部的上库力农场为研究区域，研究区域地理位置坐标为东经120°36′50.48″～120°52′56.53″，北纬50°21′11.08″～50°24′32.00″，并选取实地测量的30个样方作为单木分割验证区，其优势植被为白桦树。利用Leica LAS60扫描系统获取该区域的机载激光雷达数据，航飞时间在2012年8月16日至2012年9月25日，系统具体参数参见表1。并基于TerraScan软件把原始激光回波点分类为植被、地面、其余回波点，分类结果如图5。以下步骤为对样方编号为YF_1，中心坐标为(120.770998°，50.357450°)的单木分割详细步骤。

表1机载激光雷达系统参数

步骤2，在原始激光点云中寻找局部最大高程值点，对该点进行编码并以优势植被最小冠幅的近似圆标识该区域(样方实测最小冠幅为1m)，以此循环下去，直到所有的原始点云都被编码，最后得到初始识别的树木棵树为197(图6所示)。采用上述局部最值法进行单木分割，由于只是使用优势植被的最小冠幅R来标识各个区域，而R并不能代表所有树的冠幅，那么得到的分割结果将会存在多分现象、错误识别现象以及形状不规范性现象等。为了提高分割精度，需要其他的方法对初始分割结果进行后处理。

步骤3，在初始识别得到的结果中，如果所识别的树包含的激光回波点数小于5个点(这个点数阈值根据点云密度和最小冠幅面积而设定)，就把属于该棵树的所有编码设为未编码状态，编码值设为0。当两棵树之间的距离小于1.5米时，取消矮树的所有编码，最终得到修改后的树木棵树为45棵，如图6、图7所示。

步骤4，在实测中该样方单棵树的最大冠幅为4.2m，这里把最大半径阈值设为2.1m，根据激光点云在单树或者相邻树的聚集情况进行聚类处理，最后得到的分割棵树为38棵，而实测棵树为33棵。聚类结果如图7所示。

步骤5，根据距离把未被编码的回波点分到最近的已编码已聚类的树中，如图7所示。

步骤6，精度检验。分别从单木分割后得到的样方内树木棵树和分割后得到的点云数据进行单木树高的反演两个方面对分割结果的精度进行了检验。如图8、图9所示。在图8中，对实地测量的30个样方进行分割，分割后得到的树木棵树的精度分别是，绝对系数R²为0.7398，均方根误差RMSE为6.845棵。在图9中，对分割后得到数据进行了单木提取以及单木树高反演，树高反演精度分别是，绝对系数0.7912，均方根误差1.894米。

该实例根据本发明提出的方法，对白桦树样方的机载激光雷达点云数据进行分析，依照技术方案所述，实现了对30个样方的单木分割，并对分割后的结果进行了两种方法的精度检验，精度检验结果分别为(R²＝0.7398，RMSE＝6.845棵)、(R²＝0.7912，RMSE＝1.894米)，都取得了较好的精度，充分说明该方法的可适用性。由检验结果可知，在一般情况下，分割后的单木棵树比实测的稍多，是因为在分割时，如果单棵树存在比较大的分支，那么单棵树内的回波点云聚集情况将会改变为相邻树的点云聚集情况，造成多分，尤其是多枝丫的树出现多分的频率会更高。树高反演的均方根误差稍微显大，是因为树高反演的数据来源是经过单木分割后得到的点云数据，单木分割中本身就存在一定的误差(如错误分割、少分、多分等)，这样就导致单木分割中的误差传递到单木树高反演中，进一步带来偏差。该分割算法更适用于不存在大的分支的树种(针叶林)或者树木冠层形状比较规范的阔叶林等。

Claims (1)

1.一种基于机载激光点云聚集关系的单木分割方法，其特征在于：

步骤1、林木判别：

把回波点分类为植被回波点、地面回波点和其余回波点；

步骤2、局部最大值法初始分割：

将原始点云的编码值均设为0,0代表点云处在未编码状态；在未编码的原始点云中寻找回波点类别属于植被的且坐标Z值最大的点，并记录该点的坐标序号，同时把该点编码设为k，k≥1，用半径为R的圆标识目标树木区域，并且把该区域内的所有植被回波点都设为编码k，R为实测目标树木数据优势植被的最小冠幅；

然后在原始点云中去除已经编码了的点，树编码号变为k+1，以此循环下去，直到原始点云中所有植被点的编码都不为0；

步骤3、删除错误编码树：

当同一编码树的点云数量小于D或两颗树之间的间隔小于3R时，修改两颗树中高度低的树的编码为0,D等于激光点云密度乘以半径为R的圆面积；

步骤4、聚类：

A、使用步骤3得到编码数据，依据树编码号选取任意两颗树并确定两树之间的间隔I；

B、if判断，如果I小于maxR,maxR为实测目标树木数据优势植被的最大冠幅半径，则依次执行①、②、③、④；如果I≥maxR，则执行①、②、③、⑤；

①、比较两树的高矮，以低树的中心点为原点，不考虑高程坐标Z，建立笛卡尔直角坐标系，高树相对矮树树高的边缘植被点与其顶点形成的垂直夹角为θ,高树顶点与低树顶点连线上的点云作为高树的边缘植被点；

②、以mintanθ为阈值tanθ，即最小tanθ；

③、记录与mintanθ比较的激光点：

记录两树顶点连线为长边长度及方向，矮树冠幅为宽边的矩形所属激光点point2Peek，且该点的编码不为高树编码的点；

④、根据多木判别原理进行聚类：首先判断正切值，如果步骤③中记录的激光点对应的正切值大于mintanθ，判断该点编码对应的矮树为高树的一部分，并把该矮树对应的点编码改为高树的编码，将矩形形内所属激光点遍历；

⑤、根据单木归属原理进行聚类：如果point2Peek对应的正切值小于mintanθ，那么判断该点为高树的冠幅点，修改该点的编码为高树的编码，遍历直到修改完该棵树的所有回波点；

D、重复执行A、B步骤，直到遍历完所有组树；

步骤5：最近邻点法处理，把未被编码的激光点分到距离最近的已编码聚类的树中。

Images