本发明涉及单木识别，具体为一种结合树干探测和树冠分割的单木识别方法。
背景技术：
1、精确的森林信息对于森林可持续规划、修复以及管理具有至关重要的作用。遥感技术现已成为获取精确森林信息的有效手段，这些信息包括：单株树木的位置、树高、冠幅、材积以及生物量。遥感技术中的机载激光雷达技术是一种主动遥感技术，它从森林冠层上方向下扫描获取数据，激光脉冲可以穿透冠层直至地面，因此可以获取森林的三维结构信息。
2、现有的基于无人机载点云数据的单木识别与分割算法大多数适用于针叶林，对于阔叶林，尤其是亚热带阔叶林，由于具有树种多样、树冠形状不规则，相邻树冠重叠度高，林下结构复杂等特点，现有的算法大都无法适用，不能获得较高的单木识别与分割精度，因此需要开发一种能够适用于阔叶林的基于无人机载点云数据的单木识别与分割算法。
技术实现思路
1、本发明的目的在于提供一种结合树干探测和树冠分割的单木识别方法，以解决上述背景技术提出的问题。
2、为实现上述目的，本发明提供如下技术方案：一种结合树干探测和树冠分割的单木识别方法，包括以下步骤：
3、s1、点云预处理：对lidar点云数据进行预处理，包括去噪、点云分类和点云归一化；
4、s2、基于tpdi的单木探测方法：利用落叶期采集的lidar数据，采用基于tpdi的单木探测方法进行树干探测，获得树干探测初始结果；
5、s3、树顶点探测和树冠分割：利用落叶期或有叶期采集的lidar数据进行树顶点探测和树冠分割，树顶点探测采用局部最大值法，树冠分割使用基于标记符控制的分水岭分割算法；
6、s4、优化单木探测初始结果：将树干探测结果与树冠分割结果叠加，识别出包含多个树干的树冠分割区域，然后基于树干和枝叶的三维结构信息去除错误的树干位置。
7、优选的，在步骤s1中，所述点云预处理的具体步骤为：
8、a1、去除点云中包含的高点和低点，这些噪声点是由测量过程中仪器或环境的异常造成的；
9、a2、使用csf算法将lidar点云分为地面点和非地面点，csf算法的总体思路是将原始lidar点云倒置，将布料模型放置在倒置的点云上，分析布料节点的位置约束，模拟出地面形状，从而确定地面位置；
10、a3、利用地面点生成dem，对非地面点进行归一化处理以消除地形影响，如果地形平坦，可将dem的分辨率设为较大的值。
11、优选的，在步骤s2中，tpdi是一个表达lidar点沿着树干的分布模式的指标，定义在一个三维的规则格网上，tpdi的定义如下：
12、
13、其中，i、j代表所处网格也就是体素的行、列号，m是在(i,j)这个位置竖直方向上体素的总个数，ni,j,k是在第k个体素内的lidar点个数；
14、为降低体素大小对tpdi结果的影响，按照下式将tpdi归一化处理：
15、
16、式中：ntpdii,j代表(i,j)位置处归一化后的tpdi值，tpdimin和tpdimax是最小和最大的tpdi值。
17、在生成ntpdi网格图之后，使用一个固定大小的方形窗口在每个网格上移动，识别窗口中心处网格的ntpdi在窗口内是否是最大值，如果具有最大值，则认为该网格可能包含了树干，并提取网格中心作为可能的树干位置。窗口大小设置为：
18、
19、式中：dmin是所有树干之间的最小距离，dbhmax是所有树干的最大胸径，dh是体素的水平尺寸。
20、以每个可能的树干位置为中心提取一定半径(dmin/2)范围内的lidar点，然后对这些点应用ransac拟合三维直线，以区分树干和林下植被。
21、优选的，所述ransac是一种随离散点拟合参数模型的迭代方法，包括以下步骤：
22、b1、在lidar点数据中随机选取两个点作为假设内部点；
23、b2、根据假设的内部点拟合一个三维直线模型；
24、b3、计算剩余的点与三维直线的垂直距离，按照距离的大小将剩余点分为内部点和外部点，并记录内部点的个数；
25、b4、重复b1-b3步骤n次，选取内部点个数最多的那一次建立的三维直线模型作为最终模型；
26、拟合的三维直线代表这些lidar点在空间的延伸方向，如果三维线与垂直方向之间的夹角小于一定的阈值，则认为这些lidar点基本沿着垂直方向分布，可能是树干点，否则，将该树干位置认为是错误的探测结果，予以删除。
27、优选的，在步骤s3中，所述树顶点探测和树冠分割的具体方法为：
28、利用归一化后的点云数据生成chm：将点云数据划分为一定分辨率的规则格网，然后在每个网格内搜寻最高点，将最高点的离地高度作为该网格的值，最终得到chm；
29、采用高斯平滑法对chm进行平滑，以消除冠层表面不规则的变化。
30、优选的，所述高斯平滑法具体为：将一个固定大小的窗口在chm网格上移动，提取窗口内所有的网格值，将网格值与二维高斯函数确定的模板进行卷积，得到平滑后的chm网格值；
31、基于平滑后的chm，使用一个固定大小的方形窗口在每个网格上移动，识别窗口中心处网格的chm值在窗口内是否是最大值，如果是最大值，则认为该网格可能包含了树顶点，并提取网格中心作为可能的树顶点位置，窗口大小设置为：
32、
33、式中：cs代表chm的分辨率。
34、优选的，基于原始chm应用标记符控制的分水岭分割算法获取单个树冠的分割区域，具体步骤为：首先对chm应用基于重建的开、闭操作，提取具有固定值的联通区域作为前景标记，利用sobel边缘检测算子计算chm的梯度幅值，然后修改梯度幅值图使其在前景标记处取得局部最小值，最后对梯度幅值图进行分水岭分割，获得单个树冠的分割区域。
35、优选的，在步骤s4中，优化单木探测初始结果的具体方法为：
36、d1、识别可能错误的树干位置：对于一个包含了多个树干位置的树冠分割区域，提取一定高度区间的lidar点，应用二维dbscan算法对点进行聚类，该高度区间与用于计算tpdi的高度区间相同，dbscan算法将点定义为三类：核心点、边缘点和噪声点，从任何一个点p开始，搜寻p周围一定半径ε内的点作为邻近点，若邻近点个数不小于阈值minpts，则p为核心点，边缘点与核心点相连，但没有足够数量的邻近点，剩余的点则是噪声点，相连的核心点和边缘点构成了点簇，计算每个树干位置与各个点簇之间的平面距离，当距离小于dbhmax/2，则认为该树干位置与点簇相关，当多个树干位置与同一个点簇相关，则认为这些树干位置有可能与灌木或者大的枝叶有关，将对其作进一步的分析；
37、d2、去除树冠点：以0.1米间隔在垂直方向上分层，计算每一层lidar点个数与总点数的比值，也就是频数，从而生成高度直方图，直方图中最大峰值为maxp，出现在高度hp处，假设所有点中的最高点的高度为maxh，如果hp位于maxh的四分之三以上，cbh近似估计为频数小于0.1×maxp的最高层的高度，如果估计的cbh低于maxh的二分之一，则重新计算cbh为hp-(maxh-hp)，cbh以上的点认为是树冠点并删除，如果hp位于maxh的四分之三以下，有可能存在被上层冠层遮挡的低矮小树，或者该树干位置并非真正树干而是由大树枝造成，需要在直方图中搜寻空值区间，如果最高的空值区间大于0.5米，则删除该空值区间以上的lidar点；
38、d3、去除树枝点：树干与树枝的lidar点在空间中的分布方向不同，因而，需要移除树枝点以便做进一步的判断，该过程包含五个步骤：
39、第一步，对移除了树冠点的剩余点实施pca计算点的最大延伸方向，该方向对应pca的第一分量；
40、第二步，将剩余点投影到包含了最大延伸方向向量的竖直面内；
41、第三步，对投影点在投影平面内进行dbscan聚类，邻域半径ε设置为0.1米，最小邻近点个数minpts设置为3；
42、第四步，对每一个点簇在三维空间内进行dbscan聚类，去除孤立点和邻近点个数不超过minpts的小点簇，邻域半径ε设置为dbhmax/2，最小邻近点个数minpts设置为3；
43、第五步，通过pca计算每一个点簇的最大延伸方向以及与垂直方向的夹角，若夹角超过45°，则认为该点簇为树枝点，从点云中删除；
44、d4、区分正确的与错误的树干位置：在去除了树枝点后，如果在正确的树干位置上，剩余点应大多为树干点，基本沿垂直方向分布，如果在错误的树干位置上，剩余点应杂乱分布，或分布方向与垂直方向的夹角较大，因为一些树的树干倾斜，因而将剩余点分为两层，分别进行ransac三维线拟合，第一层与计算tpdi指标的高度层相同，第二层从第一层的顶部开始，直至最高点的高度，计算第一层拟合直线的顶端端点与第二层拟合直线底端端点之间的水平距离，如果距离不超过一定的阈值td，并且至少一条直线与垂直方向的夹角不超过45°，则认为该树干位置为正确的树干位置，保留在最终的结果中，阈值td与dbhmax有关，td定义为dbhmax+0.2。
45、优选的，在步骤d1中，计算提取的树顶点与可能错误的树干位置之间的平面距离，如果距离小于一定阈值，认为该树干位置是正确的，不做进一步分析，阈值设置为dmin/2，dmin为所有树干之间的最小间隔，否则，提取该树干位置周围一定半径内的点，半径定义为：dmin-dbhmax/2。
46、优选的，在步骤d1中，点与点之间的距离为二维平面距离。minpts定义为3，ε与树木胸径有关，定义为dbhmax/2，dbhmax是所有树木的最大胸径。
47、与现有技术相比，本发明的有益效果是：
48、本发明中，通过同时利用树干和树冠信息，在对树冠重叠度高，树种组成复杂的阔叶混交林中的单木识别时能够获得较高的单木识别与分割精度，有利于森林可持续规划、修复以及管理。
1.一种结合树干探测和树冠分割的单木识别方法，其特征在于，包括以下步骤：
2.根据权利要求1所述的结合树干探测和树冠分割的单木识别方法，其特征在于：在步骤s1中，所述点云预处理的具体步骤为：
3.根据权利要求1所述的结合树干探测和树冠分割的单木识别方法，其特征在于：在步骤s2中，tpdi是一个表达lidar点沿着树干的分布模式的指标，定义在一个三维的规则格网上，tpdi的定义如下：
4.根据权利要求3所述的结合树干探测和树冠分割的单木识别方法，其特征在于：所述ransac是一种随离散点拟合参数模型的迭代方法，包括以下步骤：
5.根据权利要求1所述的结合树干探测和树冠分割的单木识别方法，其特征在于：在步骤s3中，所述树顶点探测和树冠分割的具体方法为：
6.根据权利要求5所述的结合树干探测和树冠分割的单木识别方法，其特征在于：所述高斯平滑法具体为：将一个固定大小的窗口在chm网格上移动，提取窗口内所有的网格值，将网格值与二维高斯函数确定的模板进行卷积，得到平滑后的chm网格值；
7.根据权利要求6所述的结合树干探测和树冠分割的单木识别方法，其特征在于：基于原始chm应用标记符控制的分水岭分割算法获取单个树冠的分割区域，具体步骤为：首先对chm应用基于重建的开、闭操作，提取具有固定值的联通区域作为前景标记，利用sobel边缘检测算子计算chm的梯度幅值，然后修改梯度幅值图使其在前景标记处取得局部最小值，最后对梯度幅值图进行分水岭分割，获得单个树冠的分割区域。
8.根据权利要求3所述的结合树干探测和树冠分割的单木识别方法，其特征在于：在步骤s4中，优化单木探测初始结果的具体方法为：
9.根据权利要求8所述的结合树干探测和树冠分割的单木识别方法，其特征在于：在步骤d1中，计算提取的树顶点与可能错误的树干位置之间的平面距离，如果距离小于一定阈值，认为该树干位置是正确的，不做进一步分析，阈值设置为dmin/2，dmin为所有树干之间的最小间隔，否则，提取该树干位置周围一定半径内的点，半径定义为：dmin-dbhmax/2。
10.根据权利要求8所述的结合树干探测和树冠分割的单木识别方法，其特征在于：在步骤d1中，点与点之间的距离为二维平面距离。minpts定义为3，ε与树木胸径有关，定义为dbhmax/2，dbhmax是所有树木的最大胸径。
