基于激光点云数据的单木骨架三维重构

赵永辉，刘雪妍，吕勇，万晓玉，窦胡元，刘淑玉

摘要：针对树木三维重构过程中面临的处理速度慢、重构精度低等问题，提出一种采用激光点云数据的单木骨架三维重构方法。首先，根据点云数据类型确定组合滤波方式，以去除离群点和地面点；其次，采用一种基于内部形态描述子（ISS）和相干点漂移算法（CPD）的混合配准算法（Intrinsic Shape-Coherent Point Drift， IS-CPD），以获取单棵树木的完整点云数据；最后，采用Laplace收缩点集和拓扑细化相结合的方法提取骨架，并通过柱体构建枝干模型，实现骨架三维重构。试验结果表明，相比传统CPD算法，研究设计的配准方案精度和执行速度分别提高50%和95.8%，最终重构误差不超过2.48%。研究结果证明可有效地重构单棵树木的三维骨架，效果接近树木原型，为构建林木数字孪生环境和林业资源管理提供参考。

关键词：激光雷达；树木点云；关键点提取；树木骨架；几何模型

中圖分类号：S792.95；TN958.98文献标识码：A文章编号：1006-8023（2024）01-0128-07

3D Reconstruction of Single Wood Skeleton Based on Laser Point Cloud Data

ZHAO Yonghui， LIU Xueyan， LYU Yong， WAN Xiaoyu， DOU Huyuan， LIU Shuyu*

（College of Computer and Control Engineering， Northeast Forestry University， Harbin 150040， China）

Abstract：In response to the slow processing speed and low reconstruction accuracy encountered during the 3D reconstruction of trees， a method for 3D reconstruction of single-tree skeletons using laser point cloud data is proposed. Firstly， a combination filtering method is determined based on the point cloud data type to remove outliers and ground points. Secondly， a hybrid registration algorithm based on ISS （Intrinsic Shape Descriptor） and CPD （Coherent Point Drift algorithm）， called IS-CPD （Intrinsic Shape-Coherent Point Drift）， is employed to obtain complete point cloud data for individual trees. Finally， a method combining Laplace contraction of point sets and topological refinement is used to obtain the skeleton， and branch models are constructed using cylinders to achieve 3D skeleton reconstruction. Experimental results show that compared to traditional CPD algorithm， the proposed registration scheme improves accuracy and execution speed by 50% and 95.8% respectively， with a final reconstruction error of no more than 2.48%. The research demonstrates the effective reconstruction of the 3D skeleton of individual trees， with results close to the original trees， providing a reference for building digital twin environments of forest trees and forestry resource management.

Keywords：LiDAR; tree point cloud; key point extraction; tree skeleton; geometry model

0引言

激光雷达可用于获取目标稠密点云数据，是实现自动驾驶和三维重建的重要手段。使用机载或地基激光雷达可以获取树高、胸径和冠层等量化信息，用于树木的三维重建，为推断树木的生态结构参数和碳储量反演提供依据，也可为林业数字孪生提供数据支撑。

主流的点云数据去噪方法主要有基于密度、基于聚类和基于统计3种[1]。分离地面点和非地面点是点云数据处理的第一步，学者提出多种算法用于地面点分离。然而，即使是最先进的滤波算法，也需要设置许多复杂的参数才能实现。Zhang等[2]提出了一种新颖的布料模拟滤波算法（Cloth Simulation Filter， CSF），该算法只需调整几个参数即可实现地面点的过滤，但该算法对于点云噪声非常敏感。在点云配准方面，经典的算法是Besl等[3]提出的迭代最近点算法（Iterative Closest Point， ICP），但易出现局部最优解，从而限制了该算法的应用。因此，许多学者采用概率统计方法进行点云配准，典型的方法是相干点漂移算法（Coherent Point Drift， CPD）[4-5]等，但该方法存在运行时间长和计算复杂的问题。石珣等[6]结合曲率特征与CPD提出了一种快速配准方法，速度大大提高，但细节精确度有所下降。陆军等[7]、夏坎强[8]、史丰博等[9]对基于关键点特征匹配的点云配准方法进行了深入研究。三维树木几何重建从传统的基于规则、草图和影像重建，发展到如今借助激光雷达技术，可以构建拓扑正确的三维树木几何形态。翟晓晓等[10]以点云数据进行树木重建，由于受激光雷达视场角的约束，难以获得树冠结构的信息，因此仅重建了树干。Lin等[11]、You等[12]涉及点云骨架提取的研究，构建了树的几何和拓扑结构，但重构模型的真实感不够强。Cao等[13]使用基于Laplace算子的建模方法提取主要枝干的几何信息，拓扑连接正确，并保留了部分细枝。曹伟等[14]对点云树木建模的发展和前景进行了综述，但在结合点云数据提取骨架并重建等方面研究不足。

本研究提出一种基于骨架的方法，旨在准确地从单木的点云数据中重建三维模型。原始点云数据经过CSF算法和K维树（Kd-Tree）近邻搜素法的组合滤波后，提取了准确的单木数据。同时，基于树木特征点云的混合配准算法（Intrinsic Shape-Coherent Point Drift， IS-CPD），可显著提高配准效率。最后，通过提取单棵树木的骨架点，构造连接性，并用圆柱拟合枝干，实现了单木的三维建模。

1数据采集及预处理

1.1数据获取

数据采集自山东省潍坊市奎文区植物园内一株高约8.5 m、树龄约20 a的银杏树。使用RoboSense雷达从2个不同角度进行点云数据采集，雷达高为1.5 m，与树木水平距离约为10 m。通过对来自树木正东方向和正北方向的2组点云数据进行采集，如图1所示。

1.2点云预处理

为了提高后续处理点云数据的准确性和时效性，需要对数据进行预处理。首先，利用CSF滤波算法去除冗余的地面背景信息，该算法参数较少，分离速度快。通过使用落在重力下的布来获取地形的物理表示，单木点云可以被分离出来。

由于扫描环境和激光雷达硬件误差的影响，可能会出现离群点。因此，采用Kd-Tree算法对提取的点云进行降噪处理，提高单个树木数据的精度，以备在后续的算法使用中得到更准确的结果。通过搜索待滤波点云pi（xi，yi，zi）中每个点的空间邻近点pj（xj，yj，zj），计算之间的平均距离（di）、全局均值（μ）以及标准差（σ）。筛选符合范围（μ－α×σ≤di≤μ+α×σ）的点并过滤掉离群值（α为决定点云空间分布的参数），di、μ、σ的计算公式如下。

di=∑kj=1‖xi－yj‖kμ=∑ni=1dinσ=∑ni=1（di－μ）2n 。（1）

式中：k为决定点云密集度的参数；n为点云数量。

通过试验发现，最终选定参数k=20，α=1.2时，对点云数据进行处理结果最优，滤噪结果如图2所示，基本去除了离群噪声点和地面点同时又确保对点云模型轮廓的保护。

2单木骨架重构方法

单木骨架重构方法的过程主要包括以下几个步骤，如图3所示。首先，对预处理的2组点云数据进行特征提取，并进行精确的配准；其次，对点云进行几何收缩，获取零体积点集，并通过拓扑细化将点集细化成一维曲线，得到与点云模型基本吻合的骨架线；最后，基于骨架线对树木枝干进行圆柱拟合，以构建枝干的三维模型。

2.1三维点云配准

CPD配準是一种基于概率的点集配准算法，在对点集进行配准时，一组点集作为高斯混合模型（Gaussian Mixture Model， GMM）的质心，假设模板点集坐标为XM×D=（y1，y2，…，yM）T，另一组点集作为混合高斯模型的数据集，假设目标点集坐标为XN×D=（x1，x2，…，xN）T，N、M分别代表2组点的数目，D为Z组的维度，T为矩阵转置。通过GMM的最大后验概率得到点集之间的匹配对应关系。GMM概率密度函数如下。

p（x）=ω1N+（1－ω）∑Mm=11Mp（xm） 。（2）

式中：px|m=1（2πσ2）D2exp （－x－ym22σ2），； p（x）是概率密度函数；ω（0≤ω≤1）为溢出点的权重参数；m为1—M中的任何一个数。

GMM质心的位置通过调整变换参数（θ）的值进行改变，而变换参数的值可以通过最小化-log函数来求解。

Eθ，σ2=－∑Nn－1log∑Mm－1p（m）p（xn|m） 。（3）

式中，xn与ym之间的匹配关系可以由GMM质心的后验概率p（mxn）=p（m）p（xnm）来定义。

采用期望最大值算法进行迭代循环，从而对最大似然估计进行优化，当收敛时迭代停止。得到θ和σ2的解，即完成模板网格点集向目标网格点集的配准。

扫描设备采集的点云数据通常数量庞大，因此并非所有点云信息都对配准有效。此外，CPD算法的计算复杂度较高，匹配速度较慢。因此，本研究采用ISS（Intrinsic Shape Signaturs）算法[15]提取关键点，以降低几何信息不显著点的数量。通过对这些特征点进行精确配准，可以提高点云配准的效率。图4给出了IS-CPD配准过程。

IS-CPD点云配准算法流程如下。

（1）选择2个视角点云重叠区域。

（2）采用ISS算法提取特征点集。设点云数据有n个点，（xi，yi，zi），i=0，1，…，n-1。记Pi=（xi，yi，zi）。

①针对输入点云的每个点构建一个半径为r的球形邻域，并根据式（4）计算每个点的权重。

Wij=1||pi-pj||，|pi-pj|<r 。（4）

②根据式（5）计算各点的协方差矩阵cov及其特征值{λ1i，λ2i，λ3i}，并按从小到大的次序进行排列。

cov（pi）=∑|pi－pj|<rwij（Pi－Pj）（Pi－Pj）T∑|Pi－Pj|<rwij 。（5）

③设置阈值ε1与ε2，满足λ1iλ2iε1、λ2iλ3iε2的点即为关键点。

（3）初始化CPD算法參数。

（4）求出相关概率矩阵与后验概率p（m|xn）。

（5）利用最小负对数似然函数求出各参数的值。

（6）判断p的收敛性，若不收敛，则重复步骤（4）直到收敛。

（7）在点集数据中，利用所得到的转换矩阵，完成配准。

2.2点云枝干重建

传统的构建枝干的方法是直接在点云表面上进行重构，这种方法会导致大量畸变结构。因此，本研究先提取单木骨架线，再通过拟合圆柱来构建几何模型。图5为骨架提取并重建枝干的过程。

为精确提取树干和树枝，采用Laplace收缩法提取骨架。首先，对点云模型进行顶点邻域三角化，得到顶点的单环邻域关系。然后，计算相应的余切形式的拉普拉斯矩阵，并以此为依据收缩点云，直至模型收缩比例占初始体积的1%，再通过拓扑细化将点集细化成一维曲线。采用最远距离点球对收缩点进行采样，利用一环邻域相关性将采样点连接成初始骨架，折叠不必要的边，直到不存在三角形，得到与点云模型基本吻合的骨架线。

为准确地模拟树枝的几何形状，采用圆柱拟合方法。在树基区域，使用优化方法来获得主干的几何结构[16]。由于靠近树冠和树枝尖端的小树枝的点云数据较为杂乱，使用树木异速生长理论来控制枝干半径。最终，拟合圆柱体来得到树木点云的3D几何模型[17]，原理如图6所示。以粗度R为半径，以上端点M和下端点N为圆心生成多个圆截面，并沿着骨架线连接圆周点绘制出圆柱体，以此代表每个树枝，最终完成整棵树的枝干的绘制。

3试验结果与分析

3.1点云配准结果与分析

为验证IS-CPD配准算法的有效性，对滤波后的点云进行试验，比较该算法与原始CPD算法及石珣等[6]提出的方法在同一数据下的运行时间及均方根误差（RMSE，式中记为RMSE），其表达式见式（6），值越小表示配准效果越精确。图7及表1给出了3种配准算法的对比结果。

RMSE=∑ni－1（xi－x︿i）2n 。（6）

式中：n为点云数量；xi和x︿i分别为配准前后对应点之间欧氏距离。

经过配准结果图7和表1的分析，石珣等[5]算法虽提高了配准速度，但其细节精度下降，配准结果不佳。相比之下，CPD和IS-CPD算法均能成功地融合2个不同角度的点云，达到毫米级的精度，

2种方法可视为效果近乎一致。相比之下，本研究算法的时间复杂度要小得多。此外，由表2可知，配准时间缩短至10.77 s，平均配准精度相较CPD提高了约50%。

3.2点云枝干重建结果与分析

在几何重建部分（图8），采用基于Laplace收缩的骨架提取方法，仅需不到5次迭代，就可以将点收缩到较好的位置，如图8（b）所示。对收缩后的零体

积点集进行拓扑细化，得到与点云模型基本吻合的骨架线，如图8（c）所示。随后，对枝干进行圆柱拟合。至此，树木点云重建工作全部完成。图8（d）为树木骨架几何重建的最终结果。

本研究使用单棵树木的树高和胸径作为重建模型的精度评价指标。首先，采用树干点拟合圆柱的方法来将点云投影至圆柱轴向方向，通过求取该轴向投影的最大值和最小值来获取树高信息。同时，在Pitkanen等[18]研究方法的基础上，对树干点云进行分层切片处理，将二维平面上的分层点云进行投影，再通过圆拟合方法得到更为精确的胸径尺寸。

为验证该算法重建模型的准确性，进行20次试验，并将其与Nurunnabi等[16]的重建方法进行了比较。表2为2种方法分别获得的树高和胸径的平均值，并将其与真实测量值进行了对比。结果表明，该算法相较于Nurunnabi等[16]的重建方法具有更高的精度，胸径平均误差仅为2.48%，树高平均误差仅为1.64%。

4结论

本研究讨论了激光雷达重建单棵树木的流程，分析并改进了关键问题。充分发挥CSF滤波和Kd-Tree算法的优势，从而精准地分离出了单棵树木的数据，提高了处理速度。提出IS-CPD配准算法，可将点云配准的效率提高约95.8%。通过精确配准后的点云数据，成功提取骨架树，最终重构误差控制在2.48%以内。试验结果表明，研究方法在树木点云数据滤波、配准和骨架提取方面具有可行性，树木枝干结构重建效果良好，且重构模型可为评估农林作物、森林生态结构健康等提供支持。

【参考文献】

[1]鲁冬冬，邹进贵.三维激光点云的降噪算法对比研究[J].测绘通报，2019（S2）：102-105.

LU D D， ZOU J G. Comparative research on denoising algorithms of 3D laser point cloud[J]. Survey and Mapping Bulletin， 2019（S2）：102-105.

[2]ZHANG W， QI J， WAN P， et al. An easy-to-use airborne LiDAR data filtering method based on cloth simulation[J]. Remote Sensing， 2016， 8（6）：501.

[3]BESL P J，  MCKAY H D. A method for registration of 3-D shapes[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence， 1992， 14（2）：239-256.

[4]MYRONENKO A，SONG X.Point set registration： coherent point drift[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence， 2010， 32（12）：2262-2275.

[5]王愛丽， 张宇枭， 吴海滨， 等. 基于集成卷积神经网络的LiDAR数据分类[J]. 哈尔滨理工大学学报， 2021， 26（4）：138-145

WANG A L， ZHANG Y X， WU H B， et al. LiDAR data classification based on ensembled convolutional neural networks[J]. Journal of Harbin University of Science and Technology， 2021， 26（4）：138-145.

[6]石珣，任洁，任小康.等.基于曲率特征的漂移配准方法[J].激光与光电子学进展，2018，55（8）：248-254.

SHI X， REN J， REN X K， et al. Drift registration based on curvature characteristics[J]. Laser & Optoelectronics Progress， 2018， 55（8）：248-254.

[7]陆军，邵红旭，王伟.等.基于关键点特征匹配的点云配准方法[J].北京理工大学学报，2020，40（4）：409-415.

LU J， SHAO H X， WANG W， et al. Point cloud registration method based on key point extraction with small overlap[J]. Transactions of Beijing Institute of Technology， 2020， 40（4）：409-415.

[8]夏坎强.基于ISS特征点和改进描述子的点云配准算法研究[J].软件工程，2022，25（1）：1-5.

XIA K Q. Research on point cloud algorithm based on ISS feature points and improved descriptor[J]. Software Engineering， 2022， 25（1）：1-5.

[9]史丰博，曹琴，魏军.等.基于特征点的曲面点云配准方法[J].北京测绘，2022，36（10）：1345-1349.

SHI F B， CAO Q， WEI J， et al. Surface point cloud registration method based on feature points[J]. Beijing Surveying and Mapping， 2022， 36（10）：1345-1349.

[10]翟晓晓，邵杰，张吴明.等.基于移动LiDAR点云的树木三维重建[J].中国农业信息，2019，31（5）：84-89.

ZHAI X X， SHAO J， ZHANG W M， et al. Three-dimensional reconstruction of trees using mobile laser scanning point cloud[J]. China Agricultural Information， 2019， 31（5）：84-89.

[11]LIN G， TANG Y， ZOU X， et al. Three-dimensional reconstruction of guava fruits and branches using instance segmentation and geometry analysis[J]. Computers and Electronics in Agriculture， 2021， 184：106107.

[12]YOU A， GRIMM C， SILWAL A， et al. Semantics-guided skeletonization of upright fruiting offshoot trees for robotic pruning[J]. Computers and Electronics in Agriculture， 2022， 192：106622.

[13]CAO J J， TAGLIASACCJI A， OLSON M， et al. Point cloud skeletons via Laplacian based contraction[C]//Proceedings of the Shape Modeling International Conference. Los Alamitos： IEEE Computer Society Press， 2010：187-197.

[14]曹伟，陈动，史玉峰.等.激光雷达点云树木建模研究进展与展望[J].武汉大学学报（信息科学版），2021，46（2）：203-220.

CAO W， CHEN D， SHI Y F， et al. Progress and prospect of LiDAR point clouds to 3D tree models[J]. Geomatics and Information Science of Wuhan University， 2021， 46（2）： 203-220.

[15]YU Z. Intrinsic shape signatures： A shape descriptor for 3D object recognition[C]//IEEE International Conference on Computer Vision Workshops. IEEE， 2010.

[16]NURUNNABI A， SADAHIRO Y， LINDENBERGH R， et al. Robust cylinder fitting in laser scanning point cloud data[J]. Measurement， 2019， 138：632-651.

[17]GUO J W， XU S B， YAN D M， et al. Realistic procedural plant modeling from multiple view images[J]. IEEE Transactions on Visualization and Computer Graphics， 2020， 26（2）：1372-1384.

[18]PITKANEN T P， RAUMONEN P， KANGAS A. Measuring stem diameters with TLS in boreal forests by complementary fitting procedure[J]. ISPRS Journal of Photogrammetry and Remote Sensing， 2019， 147：294-306.