Abstract

精确测量树木的结构特征，如高度、直径、宽度和锥度，是森林资源调查的重要组成部分。目前，地面和空中激光雷达都被用来产生点云数据，通过这些数据可以确定清单指标。陆地/地面扫描通常提供每平方米数千个点的点云分辨率，从中可以观察到树干并直接进行清单测量，而航空扫描的典型分辨率(每平方米数十个点)需要使用从地面收集的库存参考数据从激光雷达变量回归清单度量。微型LiDAR传感器的最新发展使低空飞行的飞机能够以高分辨率（每平方米数百个点）从空中捕捉点云，从中开始直接看到树干信息，从而使不需要进入地面就能进行小区规模的清查成为可能。在本文中，我们使用基于深度监督机器学习的算法开发了自动树木检测、分割和树干重建的新方法，这些算法设计用于空中获取的高分辨率激光雷达点云。我们的方法能够隔离单棵树，确定树干点，并进一步建立主树干的分段模型，包括树高、直径、锥度和扫描。通过使用深度学习模型，我们的方法能够适应点云密度的变化和部分遮挡，这在从空中捕捉数据时尤为普遍。我们使用从澳大利亚新南威尔士州的两个辐射松森林上空的直升机上捕获的高分辨率激光雷达点云来呈现我们的算法的结果。

关键词: 机载激光扫描；森林清查；深度学习

1. Introduction

准确的森林清单是有效森林管理的一个重要组成部分，有助于评估受管理的商业种植园的潜在价值、评估潜在的火灾危险以及监测病虫害[1]。森林调查的重要指标包括结构指标，如树高、胸径(DBH)、树干断面积和材积，以及树干形状，如锥度和弯曲度。传统上，这些属性是在大规模取样地块的单个树木水平上使用实地工作/地面人工测量来测量的[2]。随着陆地激光扫描(TLS)传感器和系统的出现，树木的库存指标和其他相关结构属性可以被捕获并在现场自动处理[3–13]。TLS产生每平方米数千个点的激光雷达点云密度，可用于直接测量树干；然后可以开发自动化过程来分割树干点并拟合圆形或圆柱形模型，然后用于提取树级结构度量。

对于通过地面难以进入的大面积森林、偏远地区的森林或其他危险地区的调查，必须通过空中手段进行调查和调查。机载激光扫描(ALS)涉及与TLS类似的操作原理，但扫描仪安装在高空飞行的飞机上，产生密度低得多的点云。由于分辨率较低且通常缺乏直接的树干测量，过去大多数关于森林调查ALS的工作都集中在从点云信息[14–20]而不是从直接测量的单个树干回归调查指标的方法上。这些方法依赖于清单指标的实地测量样本来构建可以处理低分辨率点云的回归模型。最近，小型高分辨率激光雷达系统(例如[21])能够从低空飞行的飞机上捕捉点云密度，这些飞机能够直接测量树干的结构和形状(即密集的ALS，密度为每平方米数百个点)，但密度和覆盖范围仍然比TLS低得多。从密集的ALS中直接测量与清单相关的树干特性的能力将使不依赖于地面测量的精确航空清单成为可能，这有可能增加清单活动的区域覆盖面、效率和安全性。

在本文中，我们开发了使用密集的航空激光雷达点云来检测单个树木、分割树干点和重建树干几何形状(即半径、扫描、锥度)的新方法。与TLS不同，在TLS中，通常可以在每棵树的整个长度和圆周覆盖范围内直接观察到树干[22]，由于从空中扫描和穿过森林树冠时发生的遮挡，密集的航空点云仍然在水平和垂直方向上仍然包含树干上的许多缺失点。我们的算法可以通过利用围绕3D点云的深度学习构建的方法，适应由于航空扫描的性质而导致的缺失数据和点分辨率的变化。3D点云学习是一个新兴的研究领域[22–24]，由自动驾驶汽车[25–28]遥感、机载激光雷达和林业[13]等应用驱动。从我们的方法中为单个树木提取的几何模型可用于直接测量森林清单度量(例如，蓄积量、材积、DBH等)。从而有可能在不需要从地面收集信息的情况下进行空中森林调查。

1.1. Related Work

由于扫描角度、点密度不足和林冠遮挡来自树干的脉冲等因素，许多有价值的结构度量通常不能从ALS数据中直接测量。相反，方法通常确定可以从ALS中精确导出的变量(例如，冠层高度)如何凭经验与实地测量的结构度量相关。为了找到相关性并对这种相关性进行建模，方法包括拟合线性模型[14，15]，使用非参数回归[16，17]和copulas[29]以及其他技术。使用基于区域的方法[18–20]在林分水平估计结构度量，或者使用单棵树检测[16，17，29]或多源单树调查(MS-STI)技术[30，31]在树水平估计结构度量。通常存在一种折衷，即简单方法通用性好，但在感兴趣的树木与实地测量的树木存在显著差异的极端情况下表现不佳，对于复杂方法则相反。

最近，森林调查的航空方法转向探测和描绘单个树木。低分辨率ALS中的单棵树检测通常通过使用树冠的3D结构来完成。方法可以大致分为几类，最常见的方法是将树冠作为树冠高度模型(CHM)中的局部最大值进行检测，CHM是一种2D栅格，其中的像元值表示树冠和地面上其他植被的高度。一旦检测到树木，就使用基于标记的分水岭[32–34]或种子区域生长[14，35]等算法来描绘树木点，其中引导描绘的标记或种子就是树木检测。其他方法包括聚类[36，37]，对3D数据的2D投影进行形态学运算[38]，基于图形的方法[39，40]和一系列其他技术[41–45]。树木检测允许推断出诸如树木数量、密度和位置图等库存指标。

TLS能够从地面感知树冠下的树木。由于没有树冠的阻碍，可以直接从数据中测量树干，由于其接近目标，通常比ALS的分辨率高得多[3–13]。不是找到与LiDAR变量的相关性，而是通常围绕自动分割树干和/或对其拟合模型(例如，圆形、圆柱形)来构建方法，其然后被用于提取树级结构度量。这些方法工作得很好，因为树干在数据中被精确地重建，因此结构度量可以比ALS数据通常使用的统计估计方法更精确地提取。然而，TLS不像ALS那样适用于估计大型森林资源存量，因为它的可行区域覆盖范围小得多，并且需要基于地面的访问/实地工作，这在偏远地区、崎岖地形或有大量林下植物的森林中可能是困难的。因此，TLS可用于从较小的样本地块中自动提取结构指标，以及非破坏性地开发和更新异速生长模型[46]。

在自动驾驶汽车等应用的需求推动下，计算机视觉和机器人领域的高分辨率3D点云数据处理取得了进展。深度学习算法包括Voxnet[22]和Shapenet[23]，其中过滤器在3D中卷积，pointnet[24]及其扩展[47，48]，其中网络直接在点上训练，以及Voxelnet[49]结合了前两者，用于点云数据中的目标检测、分类和分割等任务。对于自动驾驶汽车应用，稳健的算法旨在处理在具有挑战性的非结构化户外环境中捕获的激光雷达点云数据[25–28]。有了收集高分辨率ALS数据的能力，就有可能从这些领域的最新研究成果中吸取经验，并将其应用于林业。遥感领域逐步采用了深度学习技术。它们用于处理从卫星图像[50–52]、超光谱图像[53，54]和激光雷达[55–57]收集的数据。深度学习方法在最近的工作中被纳入激光雷达林业方法，最突出的是在点密度高的TLS应用中。使用从3D数据转换为其他形式的点云进行树种分类，如波形[58]、2D投影[59]和深度图像[60]。最近，Xi等人在TLS数据中使用3D-CNN对树干进行了标记[13]。Hamraz等人将2D-CNN应用于ALS数据的2D表示，以对树木是针叶树还是落叶树进行分类[61]，Ayrey和Hayes使用3D-CNN将体素化ALS图映射到图级别的森林指标(地上生物量、树木计数和针叶百分比)[62]。

1.2. Contributions of This Work

在之前的工作中，作者开发了一种方法，使用基于区域的CNN目标检测框架和使用体素化点云表示的3D-CNN的组合来检测树木和分割树木点[21]。我们当前的工作基于并扩展了本文中的工作，方法是开发和评估新的数据表示和深度学习架构，将这些与传统的非机器学习方法进行比较，并使用来自多个密集航空点云数据集的附加注释来评估结果。我们还开发了一种通过非线性最小二乘模型拟合进行树干重建的新算法，以完成单棵树检测和描绘、树干分割和树干模型拟合的统一管道，这些管道可用于估计对林业工作者有价值的树级结构属性。

本文的具体贡献是:

先前工作[21]中提出的树检测方法通过将3D点云数据编码的新表示添加到2D光栅化摘要中进行检测而得到扩展。使用多个航空数据集进行评估。
通过将包括LiDAR返回强度的体素表示合并到学习表示中，扩展了先前工作[21]中提出的树干分割方法，并且我们开发了一种新的基于点的深度学习架构(基于Pointnet[24])，用于树点云分割。使用多个航空数据集进行评估，并比较不同的分割架构。
我们使用随机样本一致性(RANSAC)和非线性最小二乘法开发了一种新的树干重建技术，该技术将树木主干的灵活几何模型拟合到分割的树干点，以计算沿着树干长度的多个点处的树木中心线位置和树干半径。然后，该模型可用于提取库存指标，如高度、直径等。

2. Materials

2.1. Study Areas

点云数据是从一架直升机上安装的Riegl VUX-1激光雷达扫描仪获得的，该直升机飞越澳大利亚新南威尔士州Tumut(2016年11月收集)和澳大利亚新南威尔士州cara bost(2018年2月收集)的商业松树人工林。这两个森林主要由23年生(Carabost)和26年生(Tumut)辐射松树种组成，密度分别为每公顷400(Tumut)和600(Carabost)株。每个地点的地形相对平坦，稀疏的灌木丛约1-2米高，Carabost地点略高。

2.2. Data Collection

在激光雷达采集期间，直升机飞行距离地面约60-90米，导致点云密度约为每平方米300-700点，分辨率远高于从更高飞行的有人驾驶飞机收集的常规ALS(通常为每平方米5-80点)，但仍低于常规TLS，后者可能提供每平方米数千个点的扫描。生成的数据集显示了激光雷达在树干上的命中，但由于数据中的遮挡和周围森林树冠和林下植物的高度“杂乱”点，沿树干的命中经常会出现缺失部分。

3. Methods

这种方法的目的是在ALS点云中定位单棵树，并估计它们的结构属性，如胸径(DBH)、树冠宽度、树高和树干体积。为了实现这一点，提出了一个多阶段流水线，包括地面表征和移除、单棵树的描绘、将树点分割成树干和树叶以及模型拟合树干点。每棵树的结构属性可以从它的拟合模型中推断出来。

3.1. Overview

图1给出了拟议管道的高级概述。从森林上空获得的高分辨率ALS数据中，用数字地形模型（DTM）来描述地面特征，并从点云中去除地面点。然后，检测单个树木并划定其树点，以便推断属性，如树干图、树数和小区密度，以及树冠宽度和树高等树级结构指标。接下来，从各个树木的点云中分割出茎点，然后通过模型拟合过程将小圆柱体与每个检测到的树木的分割茎点进行拟合（即茎重建）。从这些圆柱体中，可以得到每棵树的结构指标，如DBH、树干轮廓和树干体积。图2显示了该管道的图形概述。

图 1. 获取库存森林属性的管道的高级概述。绿色框表示高级流程，蓝色框表示可以获取的清单属性。

图 2. 获取库存森林属性的管道图形概览。有关词干分割模块的更详细的图形描述，请参见图5。

3.2. Ground Characterisation and Removal

接收到的许多ALS脉冲都是从地面返回的。这些点对于使用DTM对地面进行表征很重要，这本身就是清单的有用组成部分。DTM可用于估计树高等属性。还需要知道哪些点是从地面返回的，以便可以将它们从点云中移除。这简化了后续的检测、分割和重建步骤。本文使用的方法产生了一个相对粗糙的DTM。如果需要更精细的分辨率（例如，用于精确测量树高），则应使用其他工具，例如Open Topography或Simple Morphological Filter(https: //pdal.io/)来自PDAL库。

为了为地面构建DTM，整个点云在xy平面上被离散为2 × 2米的bin。存储每个bin中高度最小的点（x,y,z坐标）。接下来，创建一个跨越点云大小的4 × 4米单元格的网格。K-D树[63]用于从先前存储的子集到每个网格单元的中心位置找到最近的四个点（根据xy平面中的欧几里德距离）。这四个点的平均值，由它们到网格单元中心的距离加权，计算为网格中给定xy位置的地面高度。以这种方式计算每个网格单元的高度，然后使用delauney三角剖分[64]进行网格划分，从而输出平滑的DTM。

一旦估计了DTM，就可以通过丢弃高于DTM(沿xy平面变化)的高度阈值内的所有点来移除地面点。这个阈值应该针对不同的任务进行调整。例如，在进行树木检测时，将其设置为两米是有用的，因为这将移除许多属于地面植被的点，这些点可能会被注册为假阳性树木。然而，当进行树干分割和重建时，阈值应该更低，因为DBH通常在大约1.3m处测量(尽管树干模型在没有这些点的情况下仍然可以拟合)。在实验中，2m的阈值用于树木检测，1m的阈值用于树干分割。

3.3. Individual Tree Detection

一旦从点云中移除了地面点，就将个体树木检测为点的子集(即，树木被描绘)。描绘的任务比简单地检测树冠更困难，但是它允许提取诸如树冠宽度的清单度量，并且对于树干分割是必要的。在本研究中，使用点云中物体检测的常用方法[25，26]来完成描绘，通过点云的鸟瞰投影(BEV)(即，投影到xy平面，见图3)生成2D 3通道图像，用于更快的R-CNN物体检测器[65]，这是一种专为2D图像设计的方法，用于进行边界框检测。通过将2D边界框检测(图4a)从xy平面投影到3D立方体来描绘各个树，在这种情况下，这些点属于每棵树(图4b)。

图 3.Carabost图的两种不同BEV表示的比较。(a) Carabost绘图点云。(b)绘制到喷射色标(vd)上的垂直密度。(c )垂直密度、最大高度和平均回报(vd/mh/mr)。

图 4. 在xy平面（来自BEV）中的2D栅格上预测的边界框检测示例，并将其投影到描绘树点的3D长方体。(a)2D边界框检测；(b)3D长方体检测。

BEV Faster R-CNN目标检测方法主要用于车辆感知应用[25,26]。然而，有理由证明其在林业中的应用是合理的，主要是树木在整个森林中相当均匀，具有垂直圆柱形状，可以方便地投射到BEV图像中而没有太多遮挡。[66,67]研究了BEV栅格作为ALS森林数据表示的有效性。但是，这些方法使用更简单的方法来描绘栅格中的树木（例如，分水岭）。

3.3.1. Individual Tree Detection: BEV Representations

BEV是在点云中对3D信息进行编码的三通道图像。在车辆感知任务中，BEV通道包括垂直密度、最大高度和平均返回，是通过将点云离散为体素获得的[25,26]。垂直密度是通过将每个xy网格位置的z轴上占用的体素数相加并除以占用和未占用体素的总数来计算的。最大高度计算为每个xy网格位置的最高占用单元的z位置。通过对每个xy网格位置沿z轴取所有网格中所有点的返回值的平均值来找到平均返回值。

这项工作通过测试两种不同的表示来探索BEV表示对树检测的有效性。第一种表示是转换为三通道彩色图像的垂直密度，对垂直轴上的点密度信息进行编码——类似于Windrim等人2019年提出的树检测方法[21]。在这种表示法中，树木具有特征性的外观，因为它们的茎在较低密度的周围树叶的中心附近形成了高密度的圆形。在这项工作中，网格大小为600 × 600 × 1000 m，分辨率为0.2 m。这产生了一个2D矩阵，使用“喷射”色标将其映射到3通道彩色图像，产生600 × 600 × 3彩色图像。第二个表示将垂直密度、最大高度和平均回报连接起来，形成一个三通道图像。图3显示了同一曲线图的两种不同BEV表示的示例。

3.3.2. Individual Tree Detection: Training

为了训练物体检测器，来自地面移除的点云(x和y方向的尺寸分别小于120×120m，以适合分辨率为0.2m的600 × 600网格)的几个3D裁剪物被转换成BEV图像。每个图像中的树都用轴对齐的边界框标注。灌木和部分树木背景类也被注释以减少误报树木检测的数量。部分树类包含位于作物边界的树。带有边界框标签的BEV图像用于训练更快的RCNN目标检测器。具有Resnet-101后端[68]的fast-RCNN网络使用随机梯度下降(具有动量)训练10，000次迭代，学习率为0.003，动量为0.9，批量为1。使用MS COCO数据集上训练的权重对网络进行初始化[69]。

3.3.3. Individual Tree Detection: Inference

为了推断，一个窗口在去除地面的点云的xy平面中滑动，并且在每个位置，其x和y值包含在窗口边界内的点被裁剪掉并转换为BEV。使用经过训练的Faster-RCNN模型通过边界框和类名检测BEV中的树木、灌木和部分树木。灌木和部分树边界框被丢弃。请注意，窗口会重叠滑动，这样部分树将在另一个窗口中显示为完整的树。点云中所有树的边界框被累积在一起，它们的坐标从局部BEV坐标更改为更大点云中的位置。最后，将2D边界框投影到3D长方体中，使得x和y坐标落在边界框内的所有点都被检测为树。结果是与点云中的每棵树相关联的一组x、y、z点。

3.4. Stem Segmentation

一旦每棵树被检测为一组3D点，这些点就可以进一步分割成茎和叶成分。为此，测试了两种用于点云分割的深度学习方法。第一个使用体素表示，将给定树的点集转换为3D占用网格，并使用3D卷积架构推断每个体素的语义标签，然后将其转换回点。第二种方法称为pointnet，它使用数据的原始点表示来训练神经网络来推断标签。

为了训练分割网络，使用CloudCompare软件工具[70]手动标记单个树点云。对于每棵树，点被标记为树冠植被（在手稿中称为树叶，统称树枝、树枝和树叶）、茎或杂物类（杂物是除茎或树叶之外的任何东西，例如地面植被）。

3.4.1. Stem Segmentation: 3D-FCN Architecture for Voxel Segmentation

用于树分割的体素方法，如图5a所示，借鉴了体积深度学习中的两个开创性工作。第一个是VoxNet[22]，这是一种3D-CNN架构，用于对城市环境中捕获的目标（例如汽车、行人、自行车）的LiDAR扫描进行分类。第二种方法是V-net[71]，它是一种3D全卷积编码器-解码器神经网络 (3D-FCN)，设计用于分割3D医学图像（例如MRI）。在这项工作中为树分割设计的网络采用了V-net的3D-FCN结构，但具有与VoxNet相似的卷积层（在反卷积层中镜像）。

图 5. 体素和基于点的分割方法的图形比较。 (a)基于体素的3D-FCN方法。 (b) Pointnet方法。

对于输入，具有图6a中给出的架构的网络为每棵大小为150 × 150 × 100体素的树在x、y和z轴上分别采用分辨率为0.1 × 0.1 × 0.4m的占用网格。两种类型的占用网格是实验：二进制类型（图7a），只有1和0表示占用和未知占用/空闲空间，以及非二进制类型（图7b），它对网格中的每个占用体素都有一个脉冲返回值（和未知占用/可用空间为零）。二进制类型类似于Windrim等人2019[21]中提出的词干分割方法。该网络由一个编码器和解码器组成，经过训练，可以重建目标语义网格，该语义网格由三个类中的每一个的一个二进制占用网格组成——空白空间、树叶和茎（见图5a）。为了适应目标标签的“one-hot”样式表示，对于跨三个输出网格的每组对应体素，只占用一个体素（一个），使得每个体素只有一个类标签。对于没有点的位置，自由空间网格体素将被一个占据。标记为杂波的点在输入中表示为占用的单元格，但在输出中它们占用空白空间网格。这是为了训练网络从其输出中去除杂乱。如图6a所示，网络包含五个完全卷积层，其中两个在解码器中是反卷积的。相应的编码器和解码器层之间存在跳跃连接，以在上采样（即反卷积）时保持体素的分辨率。

图 6. 分割网络架构。蓝色块代表网络的层，红色块显示数据在网络中移动时不断变化的形状。对于基于体素的3D-FCN，3D卷积发生在最后三个维度。对于Pointnet 架构，2D卷积发生在最后两个维度（尽管滤波器没有空间在最后一个维度进行实际卷积，其大小为 1）。请注意，Pointnet通常包含 MLP层，但在这里它是使用等效的2D卷积层实现的。 (a) 基于体素的3D-FCN方法。 (b)PointNet方法。

为了训练体素网络，由网络输出的三个150 × 150 × 100占用网格以类似于独热向量的方式被处理用于分类:在交叉熵损失函数将预测与标记的目标占用网格进行比较之前，在相应的单元上应用softmax非线性。这些网络是用一批批单独的树点云来训练的。每批点云(输入和标记的目标)都被动态转换为其占用网格表示。这是为了使批次可以通过围绕z轴的随机旋转和翻转来增加，这必须在体素化之前发生在点云上。所有的体素网络都是从头开始训练的。训练网络直到收敛，学习率为0.001，在500次迭代后衰减到0.0001。使用Adam优化器[72]，并通过对总损失的贡献进行加权来平衡每批的等级。批量大小和数据扩充取决于训练样本的总数，并受到GPU内存(11GB)的限制。批量大小为六棵树，每个样本有四个额外的增量(每批总共30棵树)。

图 7. 输入3D-FCN（分割网络）的数据的体素3D表示的示例比较。(a)二进制。 (b)LiDAR脉冲返回。

对于推理，使用第3.3节的方法检测为树的点云作物被转换为二进制占用网格（即体素）并通过训练有素的体素网络。网络产生三个二进制占用网格（每个类别一个），它们被转换到标记的3D点云。由于体素化过程的下采样，生成的点云通常具有较低的分辨率。为了恢复原始分辨率，从网络输出的低分辨率标记点云构建K-D树，并使用原始高分辨率点云中的点来查询K-D树，以找到最近的标记点（在阈值距离）。结果是树的高分辨率标记点云。理想情况下，低分辨率点云中没有杂波点输出，因此原始高分辨率点云中的杂波点将超过距离阈值并且没有标签。

3.4.2. Stem Segmentation: Pointnet Architecture for Point Segmentation

与体素架构不同，点网架构使用基于点的数据表示来训练网络（图5b）。这项工作中采用的方法基于流行的Pointnet[24]架构，这是一种为点设计的神经网络风格，并且足够灵活，可用于分类和分割。在这项工作中，它用于将单个树木的点映射到树叶、茎和杂波类标签（即分割）。

所用架构(图6b)与[24]类似，但由于绕z轴旋转对称，没有T-net元件。此外，为了计算效率，MLP层被实现为2D CNN。与体素架构类似，开发了两种类型的网络，一种是用x、y和z元素简单表示的点，另一种是用x、y、z和激光雷达脉冲回波表示的点。后者的架构与图6b类似，但输入尺寸为bs × 4 × n × 1，以纳入回波强度。因此，第一层中的滤波器大小为4 × 1 × 1，而不是3 × 1 × 1。两种类型网络的架构构成了一个全局特征，它是通过在跨越n个点的轴上汇集而获得的，并且在学习更多层的特征之前被复制和附加到一个较低级的特征表示中。网络输出三类中每一类的分数，即所有分数。

为了训练网络，分数通过softmax层，并使用交叉熵损失函数与单热点分类标签进行比较。虽然网络架构可以支持具有不同点数的点云，但计算约束要求在给定批次的点云中具有相同的点数。因此，为了避免显著的下采样，batchsize保持较低(2)，并且每个点云与相似大小的点云配对，以便最小化下采样。方便的是，较低的批量产生的网络泛化能力更好[73]。如在[24]中，树点云 $\in \mathbb{R}^{\mathrm{N} \times 3}$ ，其中N是点的数量，被归一化到单位球中:

$Pij‾=Pij−1N∑i=1NPi,j(1)\overline{P_{i j}}=P_{i j}-\frac{1}{N} \sum_{i=1}^N P_{i, j} \quad\quad\quad\quad(1)$

$P˙ij=Pij‾max⁡∑j=13(Pij‾)2(2)\dot{P}_{i j}=\frac{\overline{P_{i j}}}{\max \sqrt{\sum_{j=1}^3\left(\overline{P_{i j}}\right)^2}} \quad\quad\quad\quad(2)$
网络以0.00001(无衰减)的学习率从头开始训练，并使用Adam优化器训练直到收敛。与体素网络一样，使用了批次上的类平衡。

为了进行推断，被检测为树的点云通过训练好的点网体系结构一批一批地被前馈，从而没有下采样。因此，不需要任何上采样，而基于体素的网络则需要上采样(3.4.1节)。

3.5. Stem Reconstruction/Model Fitting

一旦沿着每棵树的树干的点被识别，它们被传递到处理的最后阶段，用于基于对点云的最小二乘拟合来估计树干的几何模型。在主干点上计算主成分分析(PCA ),以确定主方差的向量方向，对应于树的中心线的粗略方向(图8b)。然后根据该向量计算旋转坐标变换，使得所得的z轴沿着中心主干的方向对齐，并且所有主干点都被变换到该坐标系中。识别树的顶点(具有最大z轴值的点)。

沿着z方向，每隔1米将点分成组，对于每个bin，使用随机样本一致性(RANSAC)算法基于与点的最佳拟合来估计2D/水平圆。RANSAC算法从容器中的所有茎点中随机选择三个点，并计算与这三个点相交的圆的参数 $x_c, y_c$ 水平中心位置)和 $r$ (半径)。位于圆半径5cm以内的剩余点被计数并用于给圆评分。用不同的随机选择的点重复该过程N = 200次迭代，并且具有最大数量的拟合共有点的圆被保留作为在该bin高度的树干中心线和半径的初始估计。对于每个bin，位于圆半径5cm以内的点(内点)被保留用于进一步处理。

每个柱上每个圆的半径和圆心的x-y坐标为树的主干轮廓创建了一个初始估计值(图8c)。这些参数在处理的最后阶段被细化，该处理使用成本函数的非线性最小平方最小化，该成本函数结合了与LiDAR内点的拟合、树的顶部的定位以及沿着树干的曲率约束。由树段组成的状态向量被构造为 $x=[xc1,yc2,r1,xc2,yc2,r2,…,xcM,ycM,rM]\mathbf{x}=\left[x_{c 1}, y_{c 2}, r_1, x_{c 2}, y_{c 2}, r_2, \ldots, x_{c M}, y_{c M}, r_M\right]$ 用于 $M$ 个垂直段。然后构建非线性最小平方成本函数:

$f(x)=∑j=1M∑i=1Npj(rj−r^j,i)2+Ecurve+Etip(3)f(\mathbf{x})=\sum_{j=1}^M \sum_{i=1}^{N_{p j}}\left(r_j-\hat{r}_{j, i}\right)^2+E_{c u r v e}+E_{t i p} \quad\quad\quad\quad(3)$

其中 $N_{p j}$ 是茎节 j 的内点 LiDAR 点的数量， $r^j,i\hat{r}_{j, i}$ 是茎节 $j$ 中的内点 $i$ 与当前估计的茎节 $j$ 中心线的水平距离 $j(r^j,i=(xij−xcj)2+(yij−ycj)2j\left(\hat{r}_{j, i}=\right. \sqrt{\left(x_{i j}-x_{c j}\right)^2+\left(y_{i j}-y_{c j}\right)^2}$ 。术语 $E_{\text {curve }}$ 和 $E_{t i p}$ 代表代价术语，其惩罚每个树干部分之间的树曲率的急剧变化以及树的顶部的最终树干部分与树的顶点的偏差:
$=Wxy∑j=1M−2(xcj+2−2xcj+1+xcj+ycj+2−2ycj+1+ycj)2+Wr∑j=1M−2(rj+2−2rj+1+rj)2(4)E_{\text {curve }}=W_{x y} \sum_{j=1}^{M-2}\left(x_{c j+2}-2 x_{c j+1}+x_{c j}+y_{c j+2}-2 y_{c j+1}+y_{c j}\right)^2+W_r \sum_{j=1}^{M-2}\left(r_{j+2}-2 r_{j+1}+r_j\right)^2 \quad\quad\quad\quad(4)$

$Etip=Wtip((xcM−xtip)2+(ycM−ytipj)2)(5)E_{t i p}=W_{t i p}\left(\left(x_{c M}-x_{t i p}\right)^2+\left(y_{c M}-y_{t i p j}\right)^2\right) \quad\quad\quad\quad(5)$

其中， $W_{x y}, W_r$ 和 $W_{t i p}$ 是权重项，用于控制拟合激光雷达点、定位树梢和最小化曲率(x-y方向上的树木扫描以及树木锥度(半径)的变化)之间的平衡。对于我们的数据集，发现 $W_{x y}=1.0, W_r=10.0$ 和 $=1.0W_{\text {tip }}=1.0$ 的值对观察到的树产生良好的拟合。等式(4)中的项提供了对沿着树干长度的 $x_c, y_c$ 和 $r$ ( $x - y$ 中的曲率，称为树扫描，r中的曲率，称为树锥度)的二阶导数的离散近似。

从最初的RANSAC圆拟合开始，使用Levenburg-Marquardt优化来产生最小化成本函数 $f(x)f(\mathbf{x})$ 的精确估计 $x^=argmin⁡(f(x))\hat{\mathbf{x}}=\operatorname{argmin}(f(\mathbf{x}))$ 。代表茎几何形状的最终模型是通过将茎段连接成一个锥形的扫掠圆柱体而构建的(图8d，e)。

图 8. 树干重建过程：（a）分割点云（树干点以红色显示），（b）树干点和使用PCA估计的树干粗略主方向，（c）每个树干部分的RANSAC圆和（d）基于稳健的最小二乘拟合过程的细化茎节估计。(e) 显示了最终拟合的茎模型的示例。

4. Experimental Setup

在实验中，从Tumut地点的点云中提取了17个不重叠的地块，包括188棵树（图9），从Carabost地点的点云中提取了8个不重叠的地块，包括259棵树。为了训练R-CNN检测器，这些地块被转换为光栅图像，所有的树木都被注上了边界框（这个过程在第3.3.2节给出）。在检测和分割实验中，从这些地块（Tumut有75棵，Carabost有81棵）中取样和提取单个树木，并在点水平上进行注释。因此，每个地块都有围绕所有树木的边界框注释，以及对树木样本的点标签注释。

图 9. 从Tumut站点收集的点云数据集，实验中使用的17个提取的地块用颜色标出。带有黑色边界的地块被用于测试。图中显示了其中一个测试地块的例子，其右边是相应的点级语义注释，左边是带有边界盒注释的光栅图像。

为了评估检测和分割方法，对每个站点使用了三重交叉验证。对于每个折叠，使用单个图进行测试（该折叠独有）。测试图要么对每棵树都有点级注释（请注意，检测是在3D中评估的，而不是使用2D边界框），或者它被裁剪为每棵树都有点级注释。所有其他注释的剩余注释给定站点的图用于训练该折叠的检测和分割模型。对于Tumut站点，三个试验区有12、8和11棵树。对于Carabost站点，三个试验区有 9、17和13棵树。表1中给出了每个折叠的训练/验证/测试拆分。

表 1. 用于训练、验证和测试交叉验证每一折的分割和检测模型的树的数量。测试拆分在折叠之间没有重叠。请注意，光栅边界框注释用于训练和验证检测模型，单个树点云注释用于测试检测模型，点级注释用于训练、验证和测试分割模型。

第3.3节中提出的树木检测方法与ALS树木检测方法进行了比较。一种技术为试验地块找到一个CHM [74]，并在CHM中检测树冠作为局部最大值。然后使用标记控制的分水岭分割来描绘属于单棵树的点，类似于[32–34]中采用的方法。注意，这与基于种子的区域生长分割方法原理相似(见[14，75])。另一种技术使用DBSCAN来聚类点云，如[36]中所示，并丢弃点低于阈值的聚类。每个剩余的点簇被认为是一棵树。

将第3.4节中提出的树木分割方法与在树木点云中分割树干的其他方法进行了比较。为高分辨率TLS森林数据开发的一种方法使用点的本征特征表示来训练分类器[76]。另一种方法，也是为TLS数据开发的，将茎点识别为非垂直圆柱体内的那些点，使用线性回归进行3D线拟合[3]。

为了评估重建的树干模型，从重建的树干中提取胸高直径(DBH)的测量值，并与在 LiDAR捕获时使用传统/手动库存技术在现场进行的DBH测量值进行比较。通过计算重建轮廓在1.3m高度处的直径，通过在相邻树干段的直径之间进行线性插值，提取了DBH的估计值。DBH值与使用传统清单在现场测量的值进行了比较。

计算两组结果。首先，为了测试树干重建的准确性，独立于之前处理步骤的点分割准确性(第3.4节)，我们使用每棵树的ground truth树干点标签(gt DBH)产生树干重建和DBH估计。然后，我们使用完整处理流水线(pred DBH)提供的茎点生成茎重建。

所有实验都使用了一台64位台式计算机，配有英特尔酷睿i7-7700K四核CPU@4.20GHz处理器和Nvidia GeForce GTX 1080Ti显卡。使用TensorFlow在python中实现了神经网络。检测和分割模型的训练时间大约为几个小时。用单个图进行检测的推断时间为秒量级，在单棵树上进行分割也是如此。因此，根据森林的大小和点密度，使用流水线可以在几分钟到几小时内完全处理一大片森林(即几块地)。

7. Conclusions

本文使用基于深度监督机器学习的算法，从高分辨率航空激光雷达点云中开发了自动树木检测、分割和树干重建的新方法。我们的方法能够确定树干点，并进一步建立主树干的分段模型，该模型包括树高、直径、锥度和扫描，这在诸如森林调查的应用中可能是有用的。通过使用深度学习模型，我们的方法能够适应点云密度的变化和部分遮挡，这在从空中捕捉数据时尤为普遍。我们的组合处理流水线允许确定单个树干的特征，即使由于通过树冠成像时的遮挡而导致沿着或穿过树干的直接观测点缺失。这里开发的方法在森林科学、生态学、管理和商业林业经营的森林调查中具有潜在的应用，通过允许从空中而不是从地面/田间测量单个树木的调查度量(包括树干体积、树干形状和潜在的产品级信息的直接测量)。

未来的研究将集中在验证其他森林类型的方法。我们将探索此处开发的技术在测量由不同物种组成的森林中的树木方面的潜力，在这种情况下，树木的高度、树干形状、树冠形状和其他结构属性在林分中变化更大。我们将继续开发和扩展我们的主干重建模型，探索从分割的主干点直接回归库存指标的方法。一种可以提高这里提出的树干重构算法的精度的潜在方法是通过使用取决于树种的树木锥度的特定参数模型，以便提高直径估计，特别是沿着包含许多缺失的激光雷达点的树干部分。特定于物种的锥形模型可以与一种使用混交林内树冠中的点云信息进行物种分类的方法相结合，这种方法可以依赖于这里使用的3D点云深度学习模型的类型。

论文链接：https://mdpi-res.com/d_attachment/remotesensing/remotesensing-12-01469/article_deploy/remotesensing-12-01469-v2.pdf?version=1589272871