CN114241217B - 一种基于圆柱特征的树干点云高效提取方法 - Google Patents

Abstract

本发明属于地面激光雷达点云数据处理技术领域，具体涉及一种基于圆柱特征的树干点云高效提取方法。本发明利用地面激光雷达获取的森林植被内部三维激光点云数据，基于树干的圆柱特征，经历地面点滤波、归一化分层、树干定位和树干提取处理，建立从地面激光雷达森林点云中提取树木树冠以下树干点云的树干点云提取方法。相比现有的树干点云提取方法，本发明能够准确定位树干，并精细提取树干点云，同时，算法效率高，鲁棒性强，适用于森林等大场景中树干点云的提取，能应用于林业资源调查、森林场景重建、森林结构参数提取等诸多领域。

Description

一种基于圆柱特征的树干点云高效提取方法

技术领域

本发明属于地面激光雷达点云数据处理技术领域，具体涉及一种基于圆柱特征的树干点云高效提取方法。

背景技术

激光雷达技术是一种主动遥感技术，它能够无接触的测量物体的三维表面形态，并以海量坐标点进行表达，即点云数据。地面激光雷达，包括固定式激光雷达、车载激光雷达，以及背包式激光雷达，能够获取高精度、高密度的林下点云数据，为林业资源调查、森林场景重建、森林结构参数提取等应用提供了高效便捷的技术手段。

可以将树大致分为两部分，包括树干和树冠。相比机载激光雷达，地面激光雷达能够获取详细的树干点云，这使得树干点云提取成为可能。从树干点云中能够直接获取树木位置、胸径、树干高度等参数，树干点云还可以间接辅助森林单木分割，进一步提取树高、冠幅、生物量，叶面积指数LAI等参数，因此，树干点云提取具有巨大的现实意义。

树干点云提取的关键步骤包括树干定位和树干提取。树干定位是指从森林场景点云中确定树干所在的位置，这方面的研究较多，例如，Tao等(2015)通过欧氏聚类进行树干定位，但聚类方法缺少语义信息，易将其他地物纳入树干，Safaie等(2021)将分层的树干段点云投影到平面，转换为二维图像后，通过Hough变换检测圆形特征，显然，将三维点云转换为二维图像，不但繁琐，还丢失了树干段的高度、倾角信息，降低了树干定位的正确性和完整性，Burt等(2019)通过随机采样一致性RANSAC算法进行圆柱拟合来定位树干，该方法充分基于树干的几何特征，并且有一定的抗噪能力，但相关参数较多，不好确定，此外RANSAC算法具有随机性，估算稳健的模型需要足够的迭代次数，因此，算法效率也较低。

树干提取是指在树干定位的基础上，提取定位点所在的树干点云。大部分研究在树干定位后先进行单木分割，之后再进行树干提取、枝叶分离、参数提取等操作，现有研究表明在树干定位后可以直接提取树干点云，Burt等(2019)在树干定位后，先抽取距离树干一定水平距离内的点，然后通过一系列方法去除地面、邻近低矮植被、细枝和叶，仅保留树干和主枝，最后同样用RANSAC算法逐段提取树干点云，可见该方法对树干的提取还不够直接，效率同样较低。Safaie等(2021)同样先抽取距离树干一定水平距离内的点，然后按一定高度进行水平切片，将切片转换为二值图像后，计算每一切片的点密度，在此基础上得到切片点密度随高度变化的高度密度直方图，基于树干与冠层交界处点密度发生突变这一现象，确定分割高度，将该高度以下且邻近树干位置的点划分为树干点，由上述算法思想可以看出，该方法存在上部树干提取不完整和下部过提取地面点的问题。

发明内容

针对上述存在问题或不足，为解决现有树干点云提取方法存在效率低，以及精度较低(如上部树干欠提取、下部树干过提取)的问题；本发明提供了一种基于圆柱特征的树干点云高效提取方法，以地面激光雷达场景点云为处理对象，先通过圆柱拟合及判断方法准确高效地定位树干，再通过树干点微圆柱搜索方法精细高效地提取场景树干点云。

一种基于圆柱特征的树干点云高效提取方法，包括以下步骤：

步骤1、场景点云地面点滤波：

将场景点云中的地面点与非地面点分离；优选的，采用CSF算法。

步骤2、归一化分层：

对步骤1获取的非地面点数据按照相对地面高度等水平间距分层，具体做法如下：

2-1、遍历每个非地面点，搜索离其最近的地面点，二者的高度差即为该非地面点的相对地面高度H。

2-2、按如下公式确定每个非地面点所属的层号，对非地面点数据对应分层为L_i：

式中，dx为分层间距,取20～30cm；floor表示向下取整，L_i∈(0,1,2,3,…)；

步骤3、树干定位：

3-1、在步骤2得到的分层点云中，从下至上抽取至少2层；再对每层点云进行欧氏距离聚类，聚类距离阈值大于两倍该层点云的平均最近邻距离f₁，小于30cm。

3-2、对步骤3-1得到的每一类簇，均对其依次进行以下判断：

条件1：类簇高度E满足，E＞dz-2f₂，f₂为当前类簇点云的平均最近邻距离。

条件2：对类簇进行圆柱模型拟合，拟合误差≤0.15(取决于树干表面的粗糙程度以及树干段与圆柱形的接近程度)。

3-2-1、构建初始圆柱模型：

首先，从当前类簇点云的最高点向下和最低点向上分别按垂直距离D，5cm≤D≤10cm，分别水平截取一段点云。

然后，将截取的两段点云分别投影到水平面，使用最小二乘法拟合得到两个圆，再将两个圆沿垂直方向分别移动到所属分段的中间高度位置。

最后，以两个圆圆心的连线作为圆柱模型轴线，两个圆的平均半径作为圆柱模型半径，从而得到初始圆柱模型C′＝(L′(x′₀,y′₀,z′₀,l′,m′,n′),r′)，C′包含圆柱轴线L′上的一点(x′₀,y′₀,z′₀)、圆柱轴线L′的方向向量(l′,m′,n′)和圆柱半径r′七个参数(两个圆心连线的中点位于类簇点云的中间高度位置，可以作为定位点来表征树干段的位置，原理如图2)。

3-2-2、圆柱模型优化：

采用Levenberg-Marquardt非线性最小二乘法对步骤3-2-1所得初始圆柱模型C′进行优化，得到优化圆柱模型C″＝(L″(x″₀,y″₀,z″₀,l″,m″,n″),r″)，以及相应的圆柱轴线L″上的一点x″₀,y″₀,z″₀、圆柱轴线L″的方向向量l″,#″,n″和圆柱半径r″七个参数。

Levenberg-Marquardt非线性最小二乘法可以在给定初始模型参数的前提下，进一步提高优化后模型与类簇点云的吻合度。

3-2-3、将各类簇分别采用相应的优化圆柱模型拟合计算拟合误差；

拟合误差用于表征圆柱模型与类簇点云的匹配度，计算公式如下：

式中，N为当前类簇点云的点数，d_i为当前类簇中点/到相应优化圆柱模型轴线的距离，r为当前类簇相应优化圆柱模型的半径；

如果依次满足以上两个条件则判定该类簇为树干段点云，则取位于L″上且在该类簇中间高度位置的一点作为定位点。以上过程对每个判定为树干段的类簇点云生成了一个定位点。

3-3、合并定位点：

步骤3-2对每个判定为树干段的类簇点云均生成了一个定位点；先将步骤3-2得到的所有定位点投影到水平面，再通过欧氏聚类将属于同一棵树的定位点合并，当某棵树的定位点至少有两个时，把这些点中的最低点作为该树唯一的定位点，并保存该定位点的坐标和相应的模型半径，用于该棵树后续树干点云提取；其中，聚类距离阈值为0.1～0.5m；最终得到场景点云中各棵树的定位点。

当某棵树的定位点只有一个时，认为该点是错误的判断，将其去除。每棵树通常会产生多个定位点，属于同一棵树的定位点的水平距离较为接近，因此先将步骤3-2得到的所有定位点投影到水平面，再通过欧氏聚类将属于同一棵树的定位点合并，当某棵树的定位点至少有两个时，把这些点中的最低点作为该树唯一的定位点，并保存该定位点的坐标和相应的模型半径，用于该棵树后续树干点云提取。

步骤4、树干点云提取：

算法思想类似于在树干定位点处生成一个半径略大于树干的水平圆环，圆环沿树干上下移动，圆环半径随树干粗细变化，环内的点被吸收进树干点，当圆环碰到地面或树干分叉时停止移动，树干提取结束(原理如图3)。具体过程如下：

4-1、将步骤3得到的定位点坐标和模型半径赋值给(p_x,p_y,p_z,p_r)。

4-2、以圆柱中心为(p_x,p_y,p_z)，圆柱半径为R(1.5p_r～2p_r)，圆柱高度为h，生成一个微圆柱，5cm≤h≤10cm。

4-3、先为场景点云构建K-D tree，再通过K-D tree按如下公式确定的搜索半径来快速搜索微圆柱外接球内的点；最后，在这些外接球内点中遍历搜索位于微圆柱内部的点；

由于搜索对象为千万级别的整个场景点云，因此，先为整个场景点云构建K-Dtree，再按确定的搜索半径来搜索微圆柱外接球内的点，以提高搜索效率；最后，在这些外接球内点中遍历搜索位于微圆柱内部的点。

4-4、确定微圆柱内部的点是否为树干点；

将微圆柱柱内点投影到二维平面，再进行最小二乘圆拟合，得到拟合圆心(p_x`,p_y`)和拟合半径p_r`，最终得到一个轴线铅垂，轴线上一点为(p_x`,p_y`,p_z)，半径为p_r`的圆柱模型，用步骤3-2-3中的方法计算拟合误差；

若拟合误差小于0.25，且|(p_r`-p_r)/p_r|小于0.3，则认为该柱内点为树干点，此时，提取柱内点并更新定位参数中的p_x＝p_x`，p_y＝p_y`，p_z＝p_z(operator)pace，p_r＝p_r`，跳到步骤4-2；其中pace为步长，0<pace≤h。

否则，柱内点不为树干点，提取结束。

operator分别取+和-，operator取+或-需分别以步骤3得到的定位点坐标和模型半径赋初值各自进行，不论先后但不能混用；operator为+是从定位点向上的树干点云的提取，operator为-是从定位点向上的树干点云的提取。

步骤3和步骤4所涉及的原理：

树干定位的原理：树干点云在被水平等距切割后，表现为多个上下平面平行的正(或斜)圆柱形态，可以通过圆柱拟合得到一个逼近该点云的模型(如图2所示)。通过模型参数和模型拟合误差，可以判断用于拟合模型的点云是否为树干点云。通常，树干分布在归一化分层的下部，同时又考虑到某一层树干处可能存在遮挡或点云缺失，为了提高树干定位的速度和完整度，由下至上抽取多层进行树干定位。由于树干被切分为多个树干段，因此会被多次定位，通过对投影到水平面的定位点进行聚类可以合并属于同一棵树的定位点，当某棵树的定位点只有一个时，认为该点是错误的判断，将其去除。定位点可能出现在冠层树干上，保留每棵树的最低定位点以确保定位点位于枝下树干。

树干提取的原理：

树干提取是在树干定位的基础上，从定位点开始通过微圆柱上下搜索树干点云的过程(如图3所示)，其本质仍然是圆柱拟合，只是树干段更短，可以将其投影到平面用圆拟合替代。在树干分叉处和接近地面处，圆拟合误差会变大，或拟合半径较上一次拟合发生较大变化，因此可作为树干提取的终止条件。

本发明首先对点云进行地面点滤波，将地面激光雷达获取的树林点云分为地面点和非地面点；再进行归一化分层处理，使得相对地面高度在同一范围内的点被划分到同一层；然后基于树干段的圆柱特征，由下至上选取若干层用于树干定位，通过聚类进一步合并属于同一棵树的定位点并去除错误的定位点；最后，以初始点云为处理对象，在树干定位点处上下搜索提取树干点云，直到碰到树干分叉或地面。其流程如图1所示。

本发明利用地面激光雷达获取的森林植被内部三维激光点云数据，基于树干的圆柱特征，经历地面点滤波、归一化分层、树干定位和树干提取处理，建立从地面激光雷达森林点云中提取树木树冠以下树干点云的树干点云提取方法。

综上所述，相比现有的树干点云提取方法，本发明能够准确定位树干，并精细提取树干点云，同时，算法效率高，鲁棒性强，适用于森林等大场景中树干点云的提取，能应用于林业资源调查、森林场景重建、森林结构参数提取等诸多领域。

附图说明

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

图2为圆柱拟合原理图；a.圆柱拟合3D原理图；b.沿轴线投影平面原理图；

图3为树干提取原理图；

图4为实施例树干提取结果图；a.plot1整体结果图；b.plot2局部结果图；

图5为实施例树干提取精度验证结果图；a.plot1树干提取精度验证结果图；b.plot2树干提取精度验证结果图。

具体实施方式

以下通过两个样方作为实施例并结合附图对本发明作进一步详细说明。本实施例开发环境Microsoft Visual Studio 2013、外部库PCL1.8.0、编程语言C++。

一种基于圆柱特征的树干点云高效提取方法，包括以下步骤：

步骤1、以成都电子科技大学清水河校区内的两块树林样方作为研究区域，研究区域地理位置坐标为东经103°55'11.56"～103°56'28.80"，北纬30°44'21.02"～30°45'35.91"，其中，plot1的优势树种为楠木，plot2的优势树种为玉兰。利用Leica ScanStationC10分别获取两块样方的地面激光雷达数据，数据采集时间为2020年6月，扫描方式为多站扫描，即在样方中心和四周进行5站式扫描，以避免树与树间的遮挡。

扫描系统具体参数参见表1。并使用Cyclone软件对获取的外业数据进行拼接、抽稀、去噪等预处理。以下为对两块样地采用同样方法和参数进行树干点云提取的详细步骤。

表1 Leica ScanStation C10技术参数

根据上述实验数据，使用CSF算法将预处理后的样方点云分为地面点和非地面点。

步骤2、对步骤1所得非地面点数据进行归一化分层，分层间距设置为0.3m。

步骤3、树干定位：

首先对每层点云进行欧氏距离聚类，聚类阈值由点云密度确定，本实施例中设置为0.05m，然后，通过判定条件进行树干定位；

条件1设置为类簇高度E大于0.2m；

条件2设置拟合误差小于0.1；

最后，将属于同一棵树的水平距离接近的定位点合并，合并的欧氏聚类阈值设置为0.5m。

步骤4、树干点云提取：

根据步骤3得到的定位点坐标和模型半径，在定位点处生成微圆柱，微圆柱随树干移动，并动态地调整大小，直到碰到树干分叉或地面，以此提取环内的树干点云。其中，微圆柱高度h设置为0.1m，设置圆拟合误差小于0.25，设置半径变化百分比小于0.3，移动步长pace设置为0.1m。

步骤5、精度检验。

经过以上步骤1-4，最终得到树干点云提取结果如图4所示。分别从树干定位和树干提取两个方面进行精度检验，树干定位结果见下表，其中polt1的树干全被识别到，但也有7棵错误识别的树，经过目视检验，主要是圆柱灯杆。而plot2只有一棵树干细而弯的小树未被识别。

树干提取精度通过比较提取的树干高度与实测树干高度，作线性回归分析得到。分别在两块样方正确定位的树干中随机抽取25％用于精度检验，plot1抽取了52棵，plot2抽取了18棵，检验结果如图5所示，其中，plot1树干提取高度与实测高度的RMSE为0.266m，线性回归直线的R²为0.956，直线斜率为1.016，plot2的RMSE为0.338m，R²为0.797，斜率为1.154。以上结果表明两块样方的树干点云提取效果都不错，且树干提取高度均表现为高估，原因是在分支较小或分支与树干连接处较平缓时，算法会向上多提取一段点云，plot2的高估现象更严重，是因为玉兰树有较多的细小分支。

该实施例根据本发明提出的方法，对两块树林样方的地面激光雷达数据进行处理，依照方法步骤，采用相同的参数实现了对两块样方树干点云的提取，并对提取结果从树干定位和树干提取两个角度进行了精度检验，都取得了较好的精度，说明本发明方法在不同植被组成的树林环境均能取得良好的提取效果。此外，电脑配置为Windows 7，Intel(R)Core(TM)i5-7500，8GB RAM，plot1有两千万个点，plot2有五百万个点，两块样地均在几分钟内提取完毕，表明该方法效率高，适用于森林等大场景的树干点云提取。

表2树干定位结果

综上可见，本发明基于树干的圆柱特征，经历地面点滤波、归一化分层、树干定位和树干提取处理，建立从地面激光雷达森林点云中提取树木树冠以下树干点云的树干点云提取方法。相比现有的树干点云提取方法，本发明能够准确定位树干，并精细提取树干点云，同时，算法效率高，鲁棒性强，适用于森林等大场景中树干点云的提取，能应用于林业资源调查、森林场景重建、森林结构参数提取等诸多领域。

Claims (2)

1.一种基于圆柱特征的树干点云高效提取方法，其特征在于，包括以下步骤：

步骤1、场景点云地面点滤波：

将场景点云中的地面点与非地面点分离；

步骤2、归一化分层：

对步骤1获取的非地面点数据按照相对地面高度等水平间距分层，具体做法如下：

2-1、遍历每个非地面点，搜索离其最近的地面点，二者的高度差即为该非地面点的相对地面高度H；

2-2、按如下公式确定每个非地面点所属的层号，对非地面点数据对应分层为L_i：

式中，dz为分层间距，取20～30cm；floor表示向下取整，L_i∈(0，1，2，3，...)；

步骤3、树干定位：

3-1、在步骤2得到的分层点云中，从下至上抽取至少2层；再对每层点云进行欧氏距离聚类，聚类距离阈值大于两倍该层点云的平均最近邻距离f₁，小于30cm；

3-2、对步骤3-1得到的每一类簇，均对其依次进行以下判断：

条件1：类簇高度E满足，E＞dz-2f₂，f₂为当前类簇点云的平均最近邻距离；

条件2：对类簇进行圆柱模型拟合，拟合误差≤0.15；

3-2-1、构建初始圆柱模型：

首先，从当前类簇点云的最高点向下和最低点向上分别按垂直距离D，5cm≤D≤10cm，分别水平截取一段点云；

然后，将截取的两段点云分别投影到水平面，使用最小二乘法拟合得到两个圆，再将两个圆沿垂直方向分别移动到所属分段的中间高度位置；

最后，以两个圆圆心的连线作为圆柱模型轴线，两个圆的平均半径作为圆柱模型半径，从而得到初始圆柱模型C′＝(L′(x′₀，y′₀，z′₀，l′，m′，n′)，r′)，C′包含圆柱轴线L′上的一点(x′₀，y′₀，z′₀)、圆柱轴线L′的方向向量(l′，m′，n′)和圆柱半径r′七个参数；

3-2-2、圆柱模型优化：

采用Levenberg-Marquardt非线性最小二乘法对步骤3-2-1所得初始圆柱模型C′进行优化，得到优化圆柱模型C″＝(L″(x″₀，y″₀，z″₀，l″，m″，n″)，r″)，以及相应的圆柱轴线L″上的一点x″₀，y″₀，z″₀、圆柱轴线L″的方向向量l″，m″，n″和圆柱半径r″七个参数；

3-2-3、将各类簇分别采用相应的优化圆柱模型拟合计算拟合误差；

拟合误差用于表征圆柱模型与类簇点云的匹配度，计算公式如下：

式中，N为当前类簇点云的点数，d_i为当前类簇中点i到相应优化圆柱模型轴线的距离，r为当前类簇相应圆柱模型的半径；

如果依次满足条件1和条件2则判定该类簇为树干段点云，则取位于L″上且在该类簇中间高度位置的一点作为定位点；

3-3、合并定位点：步骤3-2对每个判定为树干段的类簇点云均生成了一个定位点；先将步骤3-2得到的所有定位点投影到水平面，再通过欧氏聚类将属于同一棵树的定位点合并，当某棵树的定位点至少有两个时，把这些点中的最低点作为该树唯一的定位点，并保存该定位点的坐标和相应的模型半径，用于该棵树后续树干点云提取；其中，聚类距离阈值为0.1～0.5m；最终得到场景点云中各棵树的定位点；

当某棵树的定位点只有一个时，认为该点是错误的判断，将其去除；

步骤4、树干点云提取：

4-1、将步骤3得到的定位点坐标和模型半径赋值给(p_x，p_y，p_z，p_r)；

4-2、以圆柱中心为(p_x，p_y，p_z)，圆柱半径为R(1.5p_r～2p_r)，圆柱高度为h，生成一个微圆柱，5cm≤h≤10cm；

4-3、先为场景点云构建K-D tree，再通过K-D tree按如下公式确定的搜索半径来快速搜索微圆柱外接球内的点；最后，在这些外接球内点中遍历搜索位于微圆柱内部的点；

4-4、确定微圆柱内部的点是否为树干点；

将微圆柱柱内点投影到二维平面，再进行最小二乘圆拟合，得到拟合圆心(p_x`，p_y`)和拟合半径p_r`，最终得到一个轴线铅垂，轴线上一点为(p_x`，p_y`，p_z)，半径为p_r`的圆柱模型，用步骤3-2-3中的方法计算拟合误差；

若拟合误差小于0.25，且|(p_r`-p_r)/p_r|小于0.3，则认为该柱内点为树干点，此时，提取柱内点并更新定位参数中的p_x＝p_x`，p_y＝p_y`，p_z＝p_z(operator)pace，p_r＝p_r`，跳到步骤4-2；其中pace为步长，0＜pace≤h；

否则，柱内点不为树干点，提取结束；

operator分别取+和-，operator取+或-需分别以步骤3得到的定位点坐标和模型半径赋初值各自进行，不论先后但不能混用；operator为+是从定位点向上的树干点云的提取，operator为-是从定位点向上的树干点云的提取。

2.如权利要求1所述基于圆柱特征的树干点云高效提取方法，其特征在于：所述步骤1中采用CSF算法分离场景点云中的地面点与非地面点。

Images