CN114972743A - 一种基于半径扩展的层次级单棵树木提取方法 - Google Patents

Abstract

本发明公开了一种基于半径扩展的层次级单棵树木提取方法，包括对场景进行分类，提取出其中的树点并利用滤波算法进行去噪；通过建立局部坐标参考系来确定投影方向，从投影点中找到局部极大值平滑差值的正负交叉点，通过设置高度阈值找到潜在树顶点，通过欧式距离对潜在树顶点进行合并优化；利用点云数据与圆的关系进行自顶向下逐层分割实现单棵树木的提取。本发明解决了提取结果受噪声干扰、场景复杂性以及数据不完整性等因素影响的问题。

Description

一种基于半径扩展的层次级单棵树木提取方法

技术领域

本发明属于计算机视觉及人工智能技术领域，涉及一种基于半径扩展的 层次级单棵树木提取方法。

背景技术

树木作为室外场景中的重要元素之一，对于当地的气候条件、空气质量 和生态系统具有十分重要的作用，间接影响着居民的健康和生活质量。单棵 树木提取已被广泛用于各种应用中，例如监视树木的生长与健康、车辆和行 人的导航和定位、大规模室外场景的3D建模，以及城市绿量估计。因此， 对于树木的检测和提取是十分必要的。而MLS技术的快速发展也为我们带 来了更加丰富的树木结构信息，为树木提取提供了更多方法。

目前，国内外学者对于树木的处理方式主要分为两种：基于栅格化数据 的方法和基于点云数据的方法。基于栅格化数据的方法通常是把点云转换为 树冠高度模型(canopyheight model,CHM)，然后使用模板匹配、区域生长和 分水岭等2D图像处理技术进行单棵树木的提取。其中，分水岭算法是一种 传统的基于拓扑理论的图像分割算法，通过在相邻集水盆的汇合处建立分水 岭能够很好地实现单棵树冠分割。Koch等通过固定窗口检测出局部极大值， 然后利用分水岭进行树木检测，但是检测结果往往取决于窗口大小。考虑到树木提取结果很大程度上依赖于局部极大值的检测，Chen等提出了一种基 于标记点控制的分水岭算法，根据树冠和树高设置窗口大小，通过将动态窗 口提取的种子点作为标记进行分水岭法来实现更加精细的单棵树木分割。 Mongus等通过使用局部拟合曲面检测CHM中的树顶点，然后将树顶点作 为分水岭算法的标记，并基于高度、面积和形状进行区域合并来减少过分割 问题。Li等提出一种多尺度的局部极大值算法，以提高不同空间分辨率下获得的局部极大值的精度，从而将激光雷达点云分割成单棵树。因此，虽然分 水岭算法的分割结果比较稳定，但是容易产生多个局部极值引起过分割现 象。而且在生成CHM的过程中，会对原始点云数据进行平滑移除等操作， 只保留表面高程信息，导致无法充分利用点云的三维信息。

基于点云数据的方法主要依赖于点云数据特有的信息，包括属性特征、 局部几何特征和全局特征，并结合树木结构特点来完成检测。国内外学者主 要采用两种策略提取树木：1)基于聚类分类的策略；2)基于形状拟合的策 略。Morsdorf等使用聚类方法完成树木的分割，将通过数字表面模型(DSM) 获得的局部极大值作为k均值聚类的种子点进行单树提取，但该方法适用于 均匀分布的点云数据。Li等通过在树干上选取合适的种子点进行双重生长分 别提取出树干与树冠部分。Zhong等分析八叉树节点的水平直方图来检测树 干，利用改进Ncut方法对重叠树冠进行分割。Fan等首先从场景中提取杆状 物体，接着计算全局和局部特征，利用随机森林进行分类实现单棵树木的提 取。上述三种方法在提取树木之前首先均需要根据地面检测算法剔除地面， 其次在提取单棵树木时需要利用树木与其他物体的特征差异去除干扰物体。 但是由于这些特征的表达能力往往依赖于点云数据的质量且难以表示场景 整体的语义特征，因此在复杂场景中的应用效果有限。近年来，深度学习在 点云领域取得了突破性进展，其中PointNet第一个应用在点云数据上的端到 端网络，它通过引入空间变换网络和最大池化分别解决了点云数据的置换不 变性和无序性问题。Luo等人利用语义分割网络从原始点云中分割树点，并 提出了PDE-Net网络预测每个树簇指向树中心的方向来实现单棵树木的提 取。因此，基于深度学习的语义分割方法有效克服了点云数据的非结构化、 不均匀性等固有特性，为从原始点云中提取树点提供了有效的解决方案。

发明内容

本发明的目的是提供一种基于半径扩展的层次级单棵树木提取方法，解 决了树点提取不准确的问题。

本发明所采用的技术方案是，一种基于半径扩展的层次级单棵树木提取 方法，具体按照如下步骤实施：

步骤1，使用主成分分析算法PCA获取原始场景中每个点的维度特征， 通过将这些特征输入到PointNet网络来对场景进行分类，提取出其中的树点 并利用滤波算法进行去噪；

步骤2，对于去噪后的树点，通过建立局部坐标参考系来确定投影方向， 从投影点中找到局部极大值平滑差值的正负交叉点，通过设置高度阈值找到 潜在树顶点，通过欧式距离对潜在树顶点进行合并优化；

步骤3，确定树顶点的位置后，基于距离进行初步聚类来获得每个簇的 矩形和圆所在的区域，利用点云数据与圆的关系进行自顶向下逐层分割实现 单棵树木的提取。

本发明的特点还在于，

步骤1具体为：

步骤1.1，维度特征提取，将维度特征引入到语义分割中；

步骤1.2，基于语义分割的树木检测，在获得维度特征后，使用PointNet 网络对原始场景的点云数据进行分割；

步骤1.3，分割结果优化，依语义分割结果提取出所有的树点，使用PCL 中的直通滤波器和统计滤波器进行噪声点和离群点的剔除。

步骤1.1具体为：

步骤1.1.1，给定点云数据集P＝{p₁,p₂,...p_N}，N为点云数据的个数，给 定邻域点的范围[k_min,k_max]，然后将邻域k_i初始化为k_min；

步骤1.1.2，从三维点云中选取点p_i局部邻域点并构建协方差矩阵，依据 PCA算法计算得到其特征值λ₁、λ₂、λ₃，λ₁≥λ₂≥λ₃；

步骤1.1.3，利用特征值计算三个不同维度特征的值，

式(1)中，pro_1D、pro_2D、以及pro_3D分别表示该点属于线状点、平面点 和散乱点的概率；

步骤1.1.4，计算点p_i的香农熵值Entropy：

Entropy＝-pro_1Dln(pro_2D)-pro_2Dln(pro_2D)-pro_3Dln(pro_3D), (2)

记录此时的熵值，熵值越小，代表该点包含的信息越小，维度特征越单 一；

步骤1.1.5，如果k_i＜k_max，转到步骤1.1.2，令k_i＝k_i+Δk，Δk为邻域k_i的 增量，值为1，计算对应p_i维度特征，否则，转到步骤1.1.6；

步骤1.1.6，统计不同邻域半径下的维度特征对应的熵值，将熵值最小时 的维度特征值作为点p_i的最佳特征值。

步骤1.2具体为：

步骤1.2.1，将点云数据的空间坐标和维度特征组合成六维数据 (x,y,z,pro_1D,pro_2D,pro_3D)数据输入到网络中，经过一个共享多层感知机MLP 和两个T-Net将N*6维的数据转变为利于分割N*64维的数据，得到点云数 据的局部特征；

步骤1.2.2，通过一个共享MLP得到N*1024的高维数据，并利用maxpool 得到1*1024的全局特征；

步骤1.2.3，将64维的局部特征和1024维的全局特征进行拼接，得到 N*1088的数据，再经过一个共享的MLP，得到表示每个点属于树点和非树 点的概率的N*2维矩阵，完成场景语义分割，将点云数据分为树点和非树点。

步骤1.3中直通滤波器和统计滤波器进行噪声点和离群点的剔除具体 为：直通滤波器通过确定点云数据X轴,Y轴,Z轴的范围，快速去除离群点； 统计滤波器通过计算每个点到临近点的平均距离，并与给定均值和方差对比 剔除掉噪声点。

步骤2具体为：

步骤2.1，根据局部坐标参考系确定投影方向，利用PCA算法计算得到 投影点簇的最大的三个特征值和特征向量，选取X轴和Z轴所在平面作为投 影方向；

步骤2.2，对噪声点进行处理，剔除投影后的重复点，基于去噪的数据 提取局部极大值，根据局部极大值定位潜在树顶点，计算所有潜在峰值点之 间的欧式距离，将峰值点之间的距离按升序排序，重新计算峰值点对之间的 距离并进行距离评估，直到所有峰值点之间的距离均不小于距离阈值时合并 结束。

步骤2.1具体为：

步骤2.1.1，根据局部坐标参考系确定投影方向：假设检测到的树点集合 T＝{t_i|i＝1,2,...,N_t}，N_t为树点个数；定义平面模型ax+by+cz+d＝0，其中a＝b＝d＝0，c＝1；将T中的点云数据进行水平面投影；

步骤2.1.2，利用PCA算法计算得到投影点簇的最大的三个特征值和特 征向量，将最大特征值对应的特征向量v₁所在方向作为X轴，中间特征值对 应的特征向量v₂所在方向作为Y轴，最小特征值对应的特征向量v₃所在方向作 为Z轴来建立局部坐标参考系；

步骤2.1.3，选取X轴和Z轴所在平面作为投影方向，将场景数据沿着Y 轴进行投影，对所有的树点计算质心，

式(3)中，n是所有投影点的质心点；

将质心n所在的XOZ平面作为投影平面，则投影平面的法向量为v₂，坐 标为(a,b,c)，设投影后的数据集为T′＝{t′_i|i＝1,2,...N_t}，t′_i表示第i个点t_i投影到XOZ 平面上的点，对于任意点t_i(t_i ^x,t_i ^y,t_i ^z)∈T，t′_i(t′_i ^x,t′_i ^y,t′_i ^z)∈T′，

因为

与v₂同向，令l为向量

的长度，则

因此

因为

且向量

在v₂方向的投影长度为向量

的长度 l，所以

即

结合式(4)和式(5)，得到点t_i在XOZ平面上的投影点t′_i的坐标为：

步骤2.2具体为：

步骤2.2.1，对噪声点进行处理，剔除投影后的重复点：对于投影轮廓上 的点t′_i(t_i ^x′,t_i ^y,t′_i ^z)和t′_j(t′_j ^x,t′_j ^y,t′_j ^z)，点t′_j(t_j ^x′,t_j ^y,t′_j ^z)表示与点 t′_i(t_i ^x′,t_i ^y,t′_i ^z)不重复的点，如果t′_i ^x＝t′_j ^x且t′_i ^z＝t′_j ^z，则只保留其中一个点； 对于t′_i ^x＝t′_j ^x的所有点，只保留z值最大的点，基于去噪的数据提取局部极 大值，从去噪后的场景获得的局部极大值集中于树冠的外部轮廓位置；

步骤2.2.2，根据局部极大值定位潜在树顶点：从X坐标的最小值到X坐 标的最大值为所有局部极大值指定顺序，表示为M＝{m_i|i＝1,2,...,N_m}，N_m为局 部极大值的点数，计算点m_i(m_i ^x,m_i ^y,m_i ^z)在Z轴的差值DM_i，

步骤2.2.3，寻找平滑差值的正负变化点来定位潜在树顶点，先利用符号 函数来判断点m_i平滑差值的正负,

式(8)中，sign(m_i)为1时代表点m_i差值为正数，sign(m_i)为0时代表点 m_i差值为零，sign(m_i)为-1时代表点m_i差值为负数，因此，当sign(m_i)＞sign(m_i+1) 时为交叉点，通过检测其中Z轴坐标超过高度阈值的点来找到潜在峰值点；

设定树顶点位于整个场景高度的1/2以上，以避免其他明显不属于树顶 点数据的干扰，设置高度阈值为，

h_th＝(z_max-z_min)×1/2+z_min (9)

式(9)中，z_max和z_min为所有数据在Z轴上的最高点和最低点，高度小 于h_th的所有交叉点都被剔除；

对于符合上述条件的点，找到其邻域内的两个点，将其中Z坐标最大的 点即为潜在树顶点；

步骤2.2.4，计算所有潜在树顶点之间的欧式距离，将潜在树顶点之间的 距离按升序排序，对于距离最近的一对潜在树顶点，判断它们之间的距离是 否小于距离阈值d_th，如果小于该阈值，则把这两个潜在树顶点用其中心点 进行替换，重新计算峰值点对之间的距离并进行距离评估，直到所有峰值点 之间的距离均不小于距离阈值时合并结束。

步骤3具体为：

步骤3.1，对于给定的树顶点G＝{G_i|i＝1,2,...,N_g}，根据其指定半径的邻域 信息进行初始聚类得到聚类簇Clu＝{Clu_i|i＝1,2,...,N_g}，N_g表示树顶点的个数；

步骤3.2，计算簇所在圆的初始半径R＝{R_i|i＝1,2,...N_g}和圆心 O＝{O_i|i＝1,2,N._g.}.；

步骤3.3，根据高度将树点分为N_l层，每一层点组成集合 L＝{L_i|i＝1,2,...N_l}，N_l表示将树点分为的层数；

步骤3.4，令f＝1，基于初始半径和圆心对L_f层数据分割，转到步骤 3.2更新当前的圆心的半径；

步骤3.5，令f＝f+1，如果f＜N_l，则对L_f层的数据进行分割，转到 步骤3.2更新圆心和半径，继续执行步骤3.5；否则，直接对该层树点进行分 割，直至所有点处理完毕，实现单棵树木的提取。

本发明的有益效果是：

本发明一种基于半径扩展的自顶向下单颗树木提取方法，通过检测树顶 点，并将它作为分水岭算法中的标记点进行自顶向下的扩展，既可以避免直 接运用分水岭算法时将三维点云数据转化二维图像时的信息丢失问题，又可 以直接利用点云的三维结构信息。本发明能够有效地将单棵树木提取出来， 并且对于相连树木也能够进行正确提取；其次，该方法能够解决由于点云数 据的不完整性和密度不均匀，导致树木没有连接起来引起的过分割和欠分割 的问题。

附图说明

图1为本发明一种基于半径扩展的自顶向下单颗树木提取方法的过程 图；

图2为本发明PointNet语义分割网络架构图；

图3为本发明树木检测结果优化的示意图，其中，图3(a)为树簇的正 视图，图3(b)为树簇的俯视图，图3(c)为优化后树簇的正视图，图3 (d)为优化后树簇的俯视图；

图4为本发明确定树点投影方向示意图，其中，图4(a)为树簇的水平 面投影，图4(b)为建立局部坐标参考系，图4(c)为XOZ平面投影；

图5为本发明树点去噪前后的对比示意图，其中，图5(a)为去噪前的 点云数据，图5(b)为去噪前点云数据的局部极大值，图5(c)为去噪后 的点云数据，图5(d)为去噪后点云数据的局部极大值；

图6为本发明树木提取方法中潜在树顶点提取结果，其中，图6(a)为 点云数据提取潜在树顶点的正视图，图6(b)为点云数据提取潜在树顶点的 俯视图；

图7为本发明树木提取方法中最优树顶点的定位，其中，图7(a)为点 云数据提取最优树顶点的正视图，图7(b)为点云数据提取最优树顶点的俯 视图；

图8为本发明树木提取方法中基于半径扩展的层次级单颗树木提取示意 图，图8(a)为基于树顶点的初始聚类，图8(b)为基于树顶点的初始聚 类的结果，图8(c)为树点分层，图8(d)为半径扩展的方法，图8(e) 为自顶向下分层分割示意图，图8(f)为单颗树木提取结果。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

本发明一种基于半径扩展的层次级单棵树木提取方法，如图1所示，具 体按照如下步骤实施：

步骤1，使用主成分分析算法PCA获取原始场景中每个点的维度特征， 通过将这些特征输入到PointNet网络来对场景进行分类，提取出其中的树点 并利用滤波算法进行去噪；

步骤1.1，维度特征提取，将维度特征引入到语义分割中；

步骤1.1.1，给定点云数据集P＝{p₁,p₂,...p_N}，N为点云数据的个数，给 定邻域点的范围[k_min,k_max]，然后将邻域k_i初始化为k_min；

步骤1.1.2，从三维点云中选取点p_i局部邻域点并构建协方差矩阵，依据 PCA算法计算得到其特征值λ₁、λ₂、λ₃，λ₁≥λ₂≥λ₃；

步骤1.1.3，利用特征值计算三个不同维度特征的值，

式(1)中，pro_1D、pro_2D、以及pro_3D分别表示该点属于线状点、平面点 和散乱点的概率；

步骤1.1.4，计算点p_i的香农熵值Entropy：

Entropy＝-pro_1Dln(pro_2D)-pro_2Dln(pro_2D)-pro_3Dln(pro_3D), (11)

记录此时的熵值，熵值越小，代表该点包含的信息越小，维度特征越单 一；

步骤1.1.5，如果k_i＜k_max，转到步骤1.1.2，令k_i＝k_i+Δk，Δk为邻域k_i的 增量，值为1，计算对应p_i维度特征，否则，转到步骤1.1.6；

步骤1.1.6，统计不同邻域半径下的维度特征对应的熵值，将熵值最小时 的维度特征值作为点p_i的最佳特征值。

步骤1.2，基于语义分割的树木检测，在获得维度特征后，如图2所示， 使用PointNet网络对原始场景的点云数据进行分割；

步骤1.2.1，将点云数据的空间坐标和维度特征组合成六维数据 (x,y,z,pro_1D,pro_2D,pro_3D)数据输入到网络中，经过一个共享多层感知机MLP (multi-layerperceptron)和两个T-Net将N*6维的数据转变为利于分割N*64 维的数据，得到点云数据的局部特征；

步骤1.2.2，通过一个共享MLP得到N*1024的高维数据，并利用maxpool 得到1*1024的全局特征；

步骤1.2.3，将64维的局部特征和1024维的全局特征进行拼接，得到 N*1088的数据，再经过一个共享的MLP，得到表示每个点属于树点和非树 点的概率的N*2维矩阵，完成场景语义分割，将点云数据分为树点和非树点。

步骤1.3，如图3所示，分割结果优化，依语义分割结果提取出所有的 树点，使用PCL中的直通滤波器和统计滤波器进行噪声点和离群点的剔除， 直通滤波器通过确定点云数据X轴,Y轴,Z轴的范围，快速去除离群点；统计 滤波器通过计算每个点到临近点的平均距离，并与给定均值和方差对比剔除 掉噪声点。

步骤2，对于去噪后的树点，通过建立局部坐标参考系来确定投影方向， 从投影点中找到局部极大值平滑差值的正负交叉点，通过设置高度阈值找到 潜在树顶点，通过欧式距离对潜在树顶点进行合并优化；

步骤2.1，根据局部坐标参考系确定投影方向，利用PCA算法计算得到 投影点簇的最大的三个特征值和特征向量，选取X轴和Z轴所在平面作为投 影方向；

步骤2.1.1，根据局部坐标参考系确定投影方向：假设检测到的树点集合 T＝{t_i|i＝1,2,...,N_t}，N_t为树点个数；定义平面模型ax+by+cz+d＝0，其中a＝b＝d＝0， c＝1；如图4所示，将T中的点云数据进行水平面投影；

步骤2.1.2，利用PCA算法计算得到投影点簇的最大的三个特征值和特 征向量，将最大特征值对应的特征向量v₁所在方向作为X轴，中间特征值对 应的特征向量v₂所在方向作为Y轴，最小特征值对应的特征向量v₃所在方向作 为Z轴来建立局部坐标参考系；

步骤2.1.3，选取X轴和Z轴所在平面作为投影方向，将场景数据沿着Y 轴进行投影，对所有的树点计算质心，

式(3)中，n是所有投影点的质心点；

将质心n所在的XOZ平面作为投影平面，则投影平面的法向量为v₂，坐 标为(a,b,c)，设投影后的数据集为T′＝{t′_i|i＝1,2,...N_t}，t′_i表示第i个点t_i投影到XOZ 平面上的点，对于任意点t_i(t_i ^x,t_i ^y,t_i ^z)∈T，t′_i(t′_i ^x,t′_i ^y,t′_i ^z)∈T′，

因为

与v₂同向，令l为向量

的长度，则

因此

因为

且向量

在v₂方向的投影长度为向量

的长度 l，所以

即

结合式(4)和式(5)，得到点t_i在XOZ平面上的投影点t′_i的坐标为：

步骤2.2，对噪声点进行处理，剔除投影后的重复点，基于去噪的数据 提取局部极大值，根据局部极大值定位潜在树顶点，计算所有潜在峰值点之 间的欧式距离，将峰值点之间的距离按升序排序，重新计算峰值点对之间的 距离并进行距离评估，直到所有峰值点之间的距离均不小于距离阈值时合并 结束。

步骤2.2.1，如图5所示，对噪声点进行处理，剔除投影后的重复点：对 于投影轮廓上的点t′_i(t′_i ^x,t′_i ^y,t′_i ^z)和t′_j(t′_j ^x,t′_j ^y,t′_j ^z)，点t′_j(t′_j ^x,t′_j ^y,t′_j ^z)表示与 点t′_i(t′_i ^x,t′_i ^y,t′_i ^z)不重复的点，如果t′_i ^x＝t′_j ^x且t′_i ^z＝t′_j ^z，则只保留其中一个点； 对于t′_i ^x＝t′_j ^x的所有点，只保留z值最大的点，基于去噪的数据提取局部极 大值，从去噪后的场景获得的局部极大值集中于树冠的外部轮廓位置；

步骤2.2.2，根据局部极大值定位潜在树顶点：从X坐标的最小值到X坐 标的最大值为所有局部极大值指定顺序，表示为M＝{m_i|i＝1,2,...,N_m}，N_m为局 部极大值的点数，计算点m_i(m_i ^x,m_i ^y,m_i ^z)在Z轴的差值DM_i，

如图6所示为检测到潜在树顶点的正视图和俯视图，从图中可以看出， 检测到的树顶点个数为8个，而实际树木个数为6棵。因此，需要进一步对 潜在树顶点进行优化处理；

步骤2.2.3，寻找平滑差值的正负变化点来定位潜在树顶点，先利用符号 函数来判断点m_i平滑差值的正负,

式(8)中，sign(m_i)为1时代表点m_i差值为正数，sign(m_i)为0时代表点 m_i差值为零，sign(m_i)为-1时代表点m_i差值为负数，因此，当sign(m_i)＞sign(m_i+1) 时为交叉点，通过检测其中Z轴坐标超过高度阈值的点来找到潜在峰值点；

设定树顶点位于整个场景高度的1/2以上，以避免其他明显不属于树顶 点数据的干扰，设置高度阈值为，

h_th＝(z_max-z_min)×1/2+z_min (18)

式(9)中，z_max和z_min为所有数据在Z轴上的最高点和最低点，高度小 于h_th的所有交叉点都被剔除；

对于符合上述条件的点，找到其邻域内的两个点，将其中Z坐标最大的 点即为潜在树顶点；

步骤2.2.4，计算所有潜在树顶点之间的欧式距离，将潜在树顶点之间的 距离按升序排序，对于距离最近的一对潜在树顶点，判断它们之间的距离是 否小于距离阈值d_th，如果小于该阈值，则把这两个潜在树顶点用其中心点 进行替换，重新计算峰值点对之间的距离并进行距离评估，直到所有峰值点 之间的距离均不小于距离阈值时合并结束，树顶点优化结果如图7所示。

步骤3，确定树顶点的位置后，基于距离进行初步聚类来获得每个簇的 矩形和圆所在的区域，如图8所示，利用点云数据与圆的关系进行自顶向下 逐层分割实现单棵树木的提取；

步骤3.1，对于给定的树顶点G＝{G_i|i＝1,2,...,N_g}，根据其指定半径的邻域 信息进行初始聚类得到聚类簇Clu＝{Clu_i|i＝1,2,...,N_g}，N_g为树顶点的个数；

步骤3.2，根据初始聚类结果，求取簇clu_i中点的x坐标和y坐标的最大 值(x_i ^max,y_i ^max)和最小值(x_i ^min,y_i ^min)，组成边界集Bou_i(x_i ^max,y_i ^max,x_i ^min,y_i ^min)，所有簇 的边界集为Bou＝{Bou_i|i＝1,2,...,N_g}。然后根据边界值计算簇clu_i所在圆的半径R_i和 圆心O_i(O_i ^x,O_i ^y)：

R_i＝((x_i ^max-x_i ^min)+(y_i ^max-y_i ^min))/4, (19)

最终，每个树顶点所在圆的半径组成集合R＝{R_i|i＝1,2,...,N_g}，圆心组成集 合O＝{O_i|i＝1,2,...,N_g}

步骤3.3，根据高度将所有树点分为N_l层，划分的层数依赖于场景中树 木的分布情况和树种差异，N_l越大单棵树木提取结果越精确，但同时处理时 间也会增加，因为场景中的树木之间的遮挡程度不高，因此选定的层数N_l为15层，则每层高度h为：

h＝(z_max-z_min)/N_l, (21)

按照从高到低的顺序对所有分层段排序得到L＝{L_i|i＝1,2,...,N_l}，其中L_i层的点数为N_Li，依次开始处理每一层的数据，分两种情况分割树点，若L_i层中 的点

位于簇Clu_j边界点Bou_j范围内，且到圆心G_j的距离小于半径R_j， 则将点

分割到簇Clu_j中，如果不满足以上条件，则按照式(13)对树点进行 分割，

步骤3.4，令f＝1，基于初始半径和圆心对L_f层数据分割，转到步骤 3.2更新当前的圆心的半径；

步骤3.5，令f＝f+1，如果f＜N_l，则对L_f层的数据进行分割，转到 步骤3.2更新圆心和半径，继续执行步骤3.5；否则，直接对该层树点进行分 割，直至所有点处理完毕，实现单棵树木的提取。

通过上述方式可知，本发明一种基于半径扩展的层次级单颗树木提取方 法，通过结合深度学习与标记点控制的分水岭算法思想，以期提高单棵树木 提取准确性。主要工作包括：1)运用深度学习的方法检测树点和非树点；2) 提出了一种基于投影的树顶点定位方法；3)设计了一种基于半径扩展的单 棵树木提取方法。解决了提取结果受噪声干扰、场景复杂性以及数据不完整 性等因素影响的问题。

Claims (9)

1.一种基于半径扩展的层次级单棵树木提取方法，其特征在于，具体按照如下步骤实施：

步骤1，使用主成分分析算法PCA获取原始场景中每个点的维度特征，通过将这些特征输入到PointNet网络来对场景进行分类，提取出其中的树点并利用滤波算法进行去噪；

步骤2，对于去噪后的树点，通过建立局部坐标参考系来确定投影方向，从投影点中找到局部极大值平滑差值的正负交叉点，通过设置高度阈值找到潜在树顶点，通过欧式距离对潜在树顶点进行合并优化；

步骤3，确定树顶点的位置后，基于距离进行初步聚类来获得每个簇的矩形和圆所在的区域，利用点云数据与圆的关系进行自顶向下逐层分割实现单棵树木的提取。

2.根据权利要求1所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤1具体为：

步骤1.1，维度特征提取，将维度特征引入到语义分割中；

步骤1.2，基于语义分割的树木检测，在获得维度特征后，使用PointNet网络对原始场景的点云数据进行分割；

步骤1.3，分割结果优化，依语义分割结果提取出所有的树点，使用PCL中的直通滤波器和统计滤波器进行噪声点和离群点的剔除。

3.根据权利要求2所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤1.1具体为：

步骤1.1.1，给定点云数据集P＝{p₁,p₂,...p_N}，N为点云数据的个数，给定邻域点的范围[k_min,k_max]，然后将邻域k_i初始化为k_min；

步骤1.1.2，从三维点云中选取点p_i局部邻域点并构建协方差矩阵，依据PCA算法计算得到其特征值λ₁、λ₂、λ₃，λ₁≥λ₂≥λ₃；

步骤1.1.3，利用特征值计算三个不同维度特征的值，

式(1)中，pro_1D、pro_2D、以及pro_3D分别表示该点属于线状点、平面点和散乱点的概率；

步骤1.1.4，计算点p_i的香农熵值Entropy：

Entropy＝-pro_1Dln(pro_2D)-pro_2Dln(pro_2D)-pro_3Dln(pro_3D), (2)

记录此时的熵值，熵值越小，代表该点包含的信息越小，维度特征越单一；

步骤1.1.5，如果k_i＜k_max，转到步骤1.1.2，令k_i＝k_i+Δk，Δk为邻域k_i的增量，值为1，计算对应p_i维度特征值，否则，转到步骤1.1.6；

步骤1.1.6，统计不同邻域半径下的维度特征对应的熵值，将熵值最小时的维度特征值作为点p_i的最佳特征值。

4.根据权利要求2所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤1.2具体为：

步骤1.2.1，将点云数据的空间坐标和维度特征组合成六维数据(x，y，z，pro_1D，pro_2D，pro_3D)数据输入到网络中，经过一个共享多层感知机MLP和两个T-Net将N*6维的数据转变为利于分割N*64维的数据，得到点云数据的局部特征；

步骤1.2.2，通过一个共享MLP得到N*1024的高维数据，并利用maxpool 得到1*1024的全局特征；

步骤1.2.3，将64维的局部特征和1024维的全局特征进行拼接，得到N*1088的数据，再经过一个共享的MLP，得到表示每个点属于树点和非树点的概率的N*2维矩阵，完成场景语义分割，将点云数据分为树点和非树点。

5.根据权利要求2所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤1.3中直通滤波器和统计滤波器进行噪声点和离群点的剔除具体为：直通滤波器通过确定点云数据X轴,Y轴,Z轴的范围，快速去除离群点；统计滤波器通过计算每个点到临近点的平均距离，并与给定均值和方差对比剔除掉噪声点。

6.根据权利要求4所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤2具体为：

步骤2.1，根据局部坐标参考系确定投影方向，利用PCA算法计算得到投影点簇的最大的三个特征值和特征向量，选取X轴和Z轴所在平面作为投影方向；

步骤2.2，对噪声点进行处理，剔除投影后的重复点，基于去噪的数据提取局部极大值，根据局部极大值定位潜在树顶点，计算所有潜在峰值点之间的欧式距离，将峰值点之间的距离按升序排序，重新计算峰值点对之间的距离并进行距离评估，直到所有峰值点之间的距离均不小于距离阈值时合并结束。

7.根据权利要求6所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤2.1具体为：

步骤2.1.1，根据局部坐标参考系确定投影方向：假设检测到的树点集合T＝{t_i|i＝1,2,...,N_t}，N_t为树点个数；定义平面模型ax+by+cz+d＝0，其中a＝b＝d＝0，c＝1；将T中的点云数据进行水平面投影；

步骤2.1.2，利用PCA算法计算得到投影点簇的最大的三个特征值和特征向量，将最大特征值对应的特征向量v₁所在方向作为X轴，中间特征值对应的特征向量v₂所在方向作为Y轴，最小特征值对应的特征向量v₃所在方向作为Z轴来建立局部坐标参考系；

步骤2.1.3，选取X轴和Z轴所在平面作为投影方向，将场景数据沿着Y轴进行投影，对所有的树点计算质心，

式(3)中，n是所有投影点的质心点；

将质心n所在的XOZ平面作为投影平面，则投影平面的法向量为v₂，坐标为(a,b,c)，设投影后的数据集为T'＝{t_i'|i＝1,2,...N_t}，t_i'表示第i个点t_i投影到XOZ平面上的点，对于任意点t_i(t_i ^x,t_i ^y,t_i ^z)∈T，t_i'(t_i'^x,t_i'^y,t_i'^z)∈T'，

因为

与v₂同向，令l为向量

的长度，则

因此

因为

且向量

在v₂方向的投影长度为向量

的长度l，所以

即

结合式(4)和式(5)，得到点t_i在XOZ平面上的投影点t_i'的坐标为：

8.根据权利要求4所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤2.2具体为：

步骤2.2.1，对噪声点进行处理，剔除投影后的重复点：对于投影轮廓上的点t_i'(t_i ^x',t_i ^y,'t_i ^z)和t_j'(t_j'^x,t_j'^y,t_j'^z)，点t_j'(t_j ^x',t_j ^y,'t_j ^z)表示与点t_i'(t_i ^x',t_i ^y,'t_i ^z)不重复的点，如果t_i'^x＝t_j'^x且t_i'^z＝t_j'^z，则只保留其中一个点；对于t_i'^x＝t_j'^x的所有点，只保留z值最大的点，基于去噪的数据提取局部极大值，从去噪后的场景获得的局部极大值集中于树冠的外部轮廓位置；

步骤2.2.2，根据局部极大值定位潜在树顶点：从X坐标的最小值到X坐标的最大值为所有局部极大值指定顺序，表示为M＝{m_i|i＝1,2,...,N_m}，N_m为局部极大值的点数，计算点m_i(m_i ^x,m_i ^y,m_i ^z)在Z轴的差值DM_i，

步骤2.2.3，寻找平滑差值的正负变化点来定位潜在树顶点，先利用符号函数来判断点m_i平滑差值的正负，

式(8)中，sign(m_i)为1时代表点mi差值为正数，sign(m_i)为0时代表点mi差值为零，sign(m_i)为-1时代表点mi差值为负数，因此，当sign(m_i)＞sign(m_i+1) 时为交叉点，通过检测其中Z轴坐标超过高度阈值的点来找到潜在峰值点；

设定树顶点位于整个场景高度的1/2以上，以避免其他明显不属于树顶点数据的干扰，设置高度阈值为，

h_th＝(z_max-z_min)×1/2+z_min (9)

式(9)中，z_max和z_min为所有数据在Z轴上的最高点和最低点，高度小于h_th的所有交叉点都被剔除；

对于符合上述条件的点，找到其邻域内的两个点，将其中Z坐标最大的点即为潜在树顶点；

步骤2.2.4，计算所有潜在树顶点之间的欧式距离，将潜在树顶点之间的距离按升序排序，对于距离最近的一对潜在树顶点，判断它们之间的距离是否小于距离阈值d_th，如果小于该阈值，则把这两个潜在树顶点用其中心点进行替换，重新计算峰值点对之间的距离并进行距离评估，直到所有峰值点之间的距离均不小于距离阈值时合并结束。

9.根据权利要求1所述的一种基于半径扩展的层次级单棵树木提取方法，其特征在于，所述步骤3具体为：

步骤3.1，对于给定的树顶点G＝{G_i|i＝1,2,...,N_g}，根据其指定半径的邻域信息进行初始聚类得到聚类簇Clu＝{Clu_i|i＝1,2,...,N_g}，N_g为树顶点的个数；

步骤3.2，计算簇所在圆的初始半径R＝{R_i|i＝1,2,...N_g}和圆心O＝{O_i|i＝1,2,N._g.}.；

步骤3.3，根据高度将树点分为N_l层，每一层点组成集合L＝{L_i|i＝1,2,...N_l}，N_l表示将树点分为的层数；

步骤3.4，令f＝1，基于初始半径和圆心对L_f层数据分割，转到步骤3.2更新当前的圆心的半径；

步骤3.5，令f＝f+1，如果f＜N_l，则对L_f层的数据进行分割，转到步骤3.2更新圆心和半径，继续执行步骤3.5；否则，直接对该层树点进行分割，直至所有点处理完毕，实现单棵树木的提取。

Images