CN106408608A - 一种从地面激光雷达点云数据中提取树干直径的方法 - Google Patents

Abstract

本发明公开了提供一种从地面激光雷达点云数据中提取树干直径的方法，步骤1，从地面三维激光扫描仪获取的树木点云中提取树干点云；步骤2，计算树干在高度h处树干横断面的定位点；步骤3，计算树干在高度h处的生长方向与树干横断面；步骤4，获取用于计算树干直径的树干点云及其凸包点；步骤5，构建经过凸包点的一条闭合三次B样条曲线以模拟围尺轨迹，提取树干直径。本发明能显著提高从树干点云中提取树干直径的精度，能满足当前林业测量精度的要求，为地面三维激光扫描仪在林业调查中的实际应用提供了技术支撑。

Description

一种从地面激光雷达点云数据中提取树干直径的方法

技术领域

本发明涉及地面激光扫描点云数据处理和信息提取技术领域，尤其是涉及一种从地面激光雷达点云数据中提取树干直径的方法。

背景技术

地面三维激光扫描技术(Terrestrial Laser Scanning Technology，简称TLS)是一种新型的测绘技术，产生于20世纪90年代，其可以快速准确的获得目标物体表面的点云数据信息，从而为构建目标物体的三维模型提供了基础数据。经过近20多年来的发展，三维激光扫描仪已经可以连续快速地对被观测物体进行非接触式测量，其通过获取物体表面至扫描仪的距离和发射强度获取大量物体表面的三维点云数据。使用TLS提取林业相关参数是近10年来的研究热点。

森林是由单棵树木组成的，在林业研究与生产中，林业工作者需要提取单棵树木(简称单木)的直径与树高等数据，为计算单棵树木的材积与生物量、构建树木生长方程等研究与应用提供基础数据。因此，使用TLS提取单木的参数是使用TLS提取林业其他参数的基础。树干直径作为单木的直接测定因子之一，是林业生产与研究中基础参数之一。目前使用TLS提取树干材积与生物量的方法多建立在使用TLS提取树干直径的基础上，因此使用TLS精确地提取树干直径是精确提取树干材积与生物量及林分相关参数的基础，也是TLS实际应用于林业测量的技术基础。

目前，使用TLS提取树干直径的方法主要有以下三种方法：拟合圆柱体、拟合圆与凸包折线。这三类方法都是根据高度在树干的某一个高度范围处(如1.2米至1.4米)选定用于提取树干直径的树干点云，以此部分选定的树干点云为基础提取树干直径；拟合圆柱体的方法是直接将选定的树干点云拟合为圆柱体，圆柱体的直径即为树干直径；拟合圆的方法则将选定的树干点云投影至水平面后在二维平面上对投影点云拟合圆，圆的直径即为树干直径；凸包折线方法将选定的树干点云投影至水平面后计算投影点集在二维平面上的凸包点，相邻凸包点相连构成凸多边形，凸多边形周长除以圆周率Ω得到树干直径。

根据树干直径的定义，树干直径是树干横断面上的直径。因此，准确定位树干的横断面是准确提取树干直径的前提。树干横断面是树干生长方向垂直的一个断面。树干是一个不规则的几何体，存在着不同程度的弯曲与扭曲，因此树干横断面不一定与水平面平行，特定高度处的树干横断面只与此高度处的树干生长方向相关。外业调查时，调查人员根据围尺所围轨迹与树干生长方向是否垂直调整围尺的测量轨迹以确保直径测量的正确性。围尺在人为拉力的作用下将紧贴着树干横断面的凸出部分，此时围尺轨迹可看作一条过横断面凸包部位的闭合连续曲线。

根据树干直径与树干横断面的定义，现有使用TLS提取树干直径的三种方法存在着如下缺点：(1)用于提取树干直径的点云并未建立在定位树干横断面的基础上，直接根据树干高度选定用于提取树干直径的点云必然导致提取的树干直径值不准确；(2)提取树干直径的拟合方法构建的树干横断面几何形状与围尺测量时构成的几何形状不符，树干横断面是一个不规则的几何体，拟合圆柱体与拟合圆构建的是规则的圆形横断面，凸包折线构建的凸多边形不能反映围尺测量时围尺轨迹的连续变化特征；现有方法构建的树干横断面形状与围尺测量直径时形成的几何形状存在着偏差，这也必然导致现有方法提取的树干直径值不准确。

发明内容

有鉴于此，本发明的目的是针对现有技术的不足，提供一种从地面激光雷达点云数据中提取树干直径的方法，用以实现是从地面三维激光扫描仪扫描的树干点云中提取给定高度h处的树干直径。

为达到上述目的，本发明采用以下技术方案：

一种从地面激光雷达点云数据中提取树干直径的方法，包括如下步骤：

步骤1，从地面三维激光扫描仪获取的树木点云中提取树干点云；

步骤2，计算树干在高度h处树干横断面的定位点；

计算树干点云中的最低点，并以所述最低点为基础计算树干的高度；

以最低点为基准，以向量为平面的法向量，在最低点的上方h处与h+预设距离处分别构建的两个平面；位于两个平面中间的树干点云构成了一个树干点云块，将树干点云块点云投影到树干高度h处的平面上得到一个平面点云集，计算平面点云集构成的凸包多边形的质心点，所述质心点为平面点云集的几何中心点与树干在高度h处树干横断面的定位点；

步骤3，计算树干在高度h处的生长方向与树干横断面；

树干在高度h处的生长方向即是树干在高度h处的树干横断面的法向量；

采用迭代方法计算树干在高度h处的生长方向；第k次迭代计算为：构建过定位点且法向量为的一个平面，所述平面为基准平面；

计算树干在高度h处的生长方向其中，协方差矩阵的最大特征值对应的特征向量为树干在高度h处的生长方向平面法向量赋值为进行第k+1次迭代；初次迭代时法向量为迭代结束条件为：相邻两次迭代得到的两个树干生长方向与的夹角θ_k小于0.5度或者θ_k与θ_k+1的角度差小于0.5度；迭代结束时的树干生长方向即为树干在高度h处的生长方向；以树干在高度h处的生长方向为平面的法向量，通过定位点构建树干在高度h处的树干横断面；

步骤4，获取用于计算树干直径的树干点云及其凸包点；

步骤5，构建经过凸包点的一条闭合三次B样条曲线以模拟围尺轨迹，提取树干直径。

优选的，所述步骤1，从地面三维激光扫描仪获取的树木点云中提取树干点云，包括：

使用地面三维激光扫描仪配备的点云处理软件对地面三维激光扫描仪在多站扫描的点云进行配准得到单木的点云数据，经移除地表点云、移除树枝点云、移除离散噪声点后得到树干点云。

优选的，所述步骤2中，遍历树干点云集得到z坐标值最小的一个点，并将得到的所述点设为树干点云的最低点。

优选的，所述预设距离为0.5厘米。

优选的，所述步骤3中，计算树干在高度h处的生长方向包括：

以基准平面为基础，分别构建5个与此平行的平面，其中3个平面位于基准平面上方，2个平面位于基准平面下方，且相邻平面间的距离为0.5厘米；

得到5个树干点云块与5个几何中心点；

对5个几何中心点采用主成分分析方法计算树干在高度h处的生长方向

优选的，所述步骤4，获取用于计算树干直径的树干点云及其凸包点，包括：

在树干横断面上方构建一个与树干横断面平行且距离为1厘米的平面，位于树干横断面与平面间的树干点云即为用于计算树干直径的树干点云，将此部分点云投影至树干横断面上得到一个三维空间上的平面点云集，采用几何变换中的旋转操作，将此平面点云集旋转为与水平面平行的一个二维空间上的平面点云集，根据凸包算法计算此平面点云集的凸包点集Q_Dhull。

优选的，所述步骤5，构建经过凸包点的一条闭合三次B样条曲线以模拟围尺轨迹，提取树干直径，包括：

将Q_Dhull首部的3个点依次加入到Q_Dhull的尾部，同时也将尾部的3个点依次加入到Q_Dhull的首部构成插值点Q＝{Q_k},k＝0,1,...,n；

使用向心参数法计算插值点Q的节点参数值令d为总弦长，参数值定义为：

再根据节点参数值使用平均值法计算节点向量U＝{u₀,u₁,...,u_m}，有

其中，p＝3表示样条曲线的次数，n＝|Q|表示插值点的个数，m＝n+p-1；然后采用全局曲线插值的方法求解非有理三次B样条曲线的控制点，根据等式关系：

其中，N_i,p(u)表示为：

得到一个系数矩阵为(n+1)×(n+1)的线性方程组，求解方程组得到n+1个控制点{P_i}，根据非有理三次B样条曲线的表达式：

得到一条连续的曲线，该曲线会在Q_Dhull的第一个节点Q₀附近相交；设置Q_Dhull的第一个节点Q₀为首节点，从首节点处开始构建，当曲线环绕一周到达首节点时停止曲线的构建即可得到完整的闭合曲线，此条闭合曲线就是模拟围尺测量轨迹的曲线；非有理三次B样条曲线是一条参数曲线，可表示为：

根据弧长公式，对于给定的参数区间u∈α,β，曲线在此区间的弧长s计算公式可表示为：

使用数值积分中的复合辛普森求积公式方法计算得到闭合曲线的长度L，则树干直径可表示为L/π。

本发明的有益效果是：

本发明根据树干点云定位在树干高度h处树干横断面，获取用于提取树干直径的树干点云及其凸包点集的方法，在定位树干横断面基础上，结合围尺测量直径时的场景，通过构建一条闭合的连续曲线模拟围尺测量轨迹以提取树干的直径；本发明提取的直径相较于现有技术，其误差小、准确性更高。

本发明首先在给定高度处准确定位树干横断面，再根据围尺宽度选择用于提取树干直径的树干点云，然后将选定的点云投影至树干的横断面上，在计算投影点云凸包点的基础上，使用三次B样条方法构建一条过凸包点的闭合曲线以模拟围尺的测径轨迹，从而提取树干直径。实验表明，本发明能较好地模拟围尺的测量轨迹，将从TLS提取树干直径的精度提高到了毫米级，满足了林业测量的精度要求。

本发明的其它特征和优点将在随后的说明书中阐述，并且，部分地从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点可通过在所写的说明书、权利要求书、以及附图中所特别指出的结构来实现和获得。

附图说明

图1是本发明的方法流程图。

图2是定位树干横断面的示意图。

图3是获取用于计算树干直径的点云示例图。

图4是模拟的围尺轨迹示例图。

具体实施方式

下面结合附图和实施例对本发明作进一步描述。

如图1所示，本发明提供一种从地面激光雷达点云数据中提取树干直径的方法，包括如下步骤：

步骤1，从地面三维激光扫描仪获取的树木点云中提取树干点云；

步骤2，计算树干在高度h处树干横断面的定位点；

计算树干点云中的最低点，并以所述最低点为基础计算树干的高度；

以最低点为基准，以向量为平面的法向量，在最低点的上方h处与h+预设距离处分别构建的两个平面；位于两个平面中的树干点云构成了一个树干点云块，将树干点云块点云投影到树干高度h处的平面上得到一个平面点云集，计算平面点云集构成的凸包多边形的质心点，所述质心点为平面点云集的几何中心点与树干在高度h处树干横断面的定位点；

步骤3，计算树干在高度h处的生长方向与树干横断面；

树干在高度h处的生长方向即是树干在高度h处的树干横断面的法向量；

采用迭代方法计算树干在高度h处的生长方向；第k次迭代计算为：构建过定位点且法向量为的一个平面，所述平面为基准平面；

计算树干在高度h处的生长方向其中，协方差矩阵的最大特征值对应的特征向量为树干在高度h处的生长方向平面法向量赋值为进行第k+1次迭代；初次迭代时法向量为迭代结束条件为：相邻两次迭代得到的两个树干生长方向与的夹角θ_k小于0.5度或者θ_k与θ_k+1的角度差小于0.5度；迭代结束时的树干生长方向即为树干在高度h处的生长方向；以树干在高度h处的生长方向为平面的法向量，通过定位点构建树干在高度h处的树干横断面；

步骤4，获取用于计算树干直径的树干点云及其凸包点；

步骤5，构建经过凸包点的一条闭合三次B样条曲线以模拟围尺轨迹，提取树干直径。

在一个实施例中，本发明提供一种从地面激光雷达点云数据中提取树干直径的方法，包括：

步骤1，从地面三维激光扫描仪获取的树木点云中提取树干点云。

使用地面三维激光扫描仪配备的点云处理软件对地面三维激光扫描仪在多站扫描的点云进行配准得到单木的点云数据，经移除地表点云、移除树枝点云、移除离散噪声点后得到树干点云。

步骤2，计算树干在高度h处树干横断面的定位点。

根据z轴的坐标值，计算树干点云中的最低点，并以此点为基础计算树干的高度；以最低点为基准，以向量为平面的法向量，在最低点的上方h处与h+0.5厘米处分别构建的两个平面，树干高度h处的平面如图2所示。图2中，虚线2为计算定位点时的下平面，点1为定位点，虚线3为树干在高度h处的的横断面，箭头所指方向为树干横断面的法向量。图3是获取用于计算树干直径的点云示例图。图3中，虚线2是树干的横断面，虚线1是与虚线2平行且到虚线2所在平面距离为1厘米的一个平面，位于虚线1与虚线2中间的树干点云为即为用于计算树干直径的点云。位于两个平面中的树干点云构成了一个树干点云块，将此部分点云投影到树干高度h处的平面上得到一个平面点云集，计算此平面点云集构成的凸包多边形的质心点，质心点即为此平面点云集的几何中心点，同时也是树干在高度h处树干横断面的定位点。

步骤3，计算树干在高度h处的生长方向与树干横断面。

树干在高度h处的生长方向即是树干在高度h处的树干横断面的法向量。采用迭代方法计算树干在高度h处的生长方向，第k次迭代计算为：构建过定位点且法向量为的一个平面，称此平面为基准平面。

以基准平面为基础，分别构建5个与此平行的平面，其中3个平面位于基准平面上方，2个平面位于基准平面下方，且相邻平面间的距离为0.5厘米；由此得到5个树干点云块与5个几何中心点，对5个几何中心点采用主成分分析方法计算树干在高度h处的生长方向其中，协方差矩阵的最大特征值对应的特征向量为树干在高度h处的生长方向平面法向量赋值为进行第k+1次迭代；初次迭代时法向量为迭代结束条件为：相邻两次迭代得到的两个树干生长方向与的夹角θ_k小于0.5度或者θ_k与θ_k+1的角度差小于0.5度；迭代结束时的树干生长方向即为树干在高度h处的生长方向；以树干在高度h处的生长方向为平面的法向量，过定位点构建树干在高度h处的树干横断面。

步骤4，获取用于计算树干直径的树干点云及其凸包点。

在得到树干在高度h处的树干横断面的基础上，在此树干横断面上方构建一个与树干横断面平行且距离为1厘米(围尺的宽度)的平面，位于树干横断面与平面间的树干点云即为用于计算树干直径的树干点云，将此部分点云投影至树干横断面上得到一个三维空间上的平面点云集，采用几何变换中的旋转操作，将此平面点云集旋转为与水平面平行的一个二维空间上的平面点云集，根据凸包算法计算此平面点云集的凸包点集Q_Dhull。

步骤5，构建过凸包点的一条闭合三次B样条曲线以模拟围尺轨迹，从而提取树干直径。

将Q_Dhull首部的3个点依次加入到Q_Dhull的尾部，同时也将尾部的3个点依次加入到Q_Dhull的首部构成插值点Q＝{Q_k},k＝0,1,...,n；使用向心参数法计算插值点Q的节点参数值令d为总弦长，参数值定义为

再根据节点参数值使用平均值法计算节点向量U＝{u₀,u₁,...,u_m}，有

其中，p＝3表示样条曲线的次数，n＝|Q|表示插值点的个数，m＝n+p-1；

然后采用全局曲线插值的方法求解非有理三次B样条曲线的控制点，根据等式关系：

其中，N_i,p(u)表示为

得到一个系数矩阵为(n+1)×(n+1)的线性方程组，求解方程组得到n+1个控制点{P_i}，根据非有理三次B样条曲线的表达式：

得到一条连续的曲线，该曲线会在Q_Dhull的第一个节点Q₀附近相交，为此可设置Q_Dhull的第一个节点Q₀为首节点，从首节点处开始构建，当曲线环绕一周到达首节点时停止曲线的构建即可得到完整的闭合曲线，此条闭合曲线就是模拟围尺测量轨迹的曲线。图4模拟的围尺轨迹示例图，(a)与(c)分别是树干横断面点云旋转后得到的平面点云集，(b)与(d)是对应的模拟围尺测径轨迹的闭合曲线。

非有理三次B样条曲线是一条参数曲线，可表示为：

其中，u是参数，(x_i,y_i,z_i)是曲线s(u)上的点，根据弧长公式，对于给定的参数区间u∈α,β，曲线在此区间的弧长s计算公式可表示为：

本发明中，f(u)的原函数难以使用初等函数明确表达，因此使用数值积分中的复合辛普森求积公式方法计算得到闭合曲线的长度L，则树干直径可表示为L/π。

以外业调查时获取的31条树干直径数据作为真实值，采用本发明中的获取用于计算树干直径的树干点云的方法获取对应位置处的树干点云，再分别使用拟合圆、凸包折线与本发明中的模拟围尺测径的方法提取树干直径，与外业调查的实测直径相比，3种方法提取树干直径的精度比较如下表所示。

表1不同方法实验结果对照比

	MAPE/％	RMSE/cm	R2
				拟合圆	2.6659	1.0315	0.9853
凸包折线	0.3448	0.1628	0.9868
				本发明	0.3499	0.1556	0.9885

表中评测指标的含义及计算公式如下：平均绝对百分比误差MAPE(mean absolutepercentage error)、均方根误差RMSE(root mean square error)和决定系数R²；对应计算公式如下：

n为样本单元，D是直径的实测值，其平均值用表示，表示从树干点云提取的直径值。上述评价指标中MAPE与RMSE越小越好(R²则越大越好)，这说明从树干点云中提取直径的准确性越高。

本发明通过获取单木的树干点云数据，定位树干在指定高度处的树干横断面，获取用于提取树干直径的树干点云并计算得到用于模拟围尺测量轨迹的凸包插值点，根据凸包插值点构建一条闭合的非有理三次B样条曲线模拟围尺测径轨迹，最后根据闭合曲线长度计算树干直径。与现有方法相比，本发明能显著提高从树干点云中提取树干直径的精度，能满足当前林业测量精度的要求，为地面三维激光扫描仪在林业调查中的实际应用提供了技术支撑。

需要说明的是，以上实施例仅用以说明本发明的技术方案而非限制，本领域普通技术人员对本发明的技术方案所做的其他修改或者等同替换，只要不脱离本发明技术方案的精神和范围，均应涵盖在本发明的权利要求范围当中。

Claims (7)

1.一种从地面激光雷达点云数据中提取树干直径的方法，其特征在于：包括如下步骤：

步骤1，从地面三维激光扫描仪获取的树木点云中提取树干点云；

步骤2，计算树干在高度h处树干横断面的定位点；

计算树干点云中的最低点，并以所述最低点为基础计算树干的高度；

以最低点为基准，以向量为平面的法向量，在最低点的上方h处与h+预设距离处分别构建的两个平面；位于两个平面中间的树干点云构成了一个树干点云块，将树干点云块点云投影到树干高度h处的平面上得到一个平面点云集，计算平面点云集构成的凸包多边形的质心点，所述质心点为平面点云集的几何中心点与树干在高度h处树干横断面的定位点；

步骤3，计算树干在高度h处的生长方向与树干横断面；

树干在高度h处的生长方向即是树干在高度h处的树干横断面的法向量；

采用迭代方法计算树干在高度h处的生长方向；第k次迭代计算为：构建过定位点且法向量为的一个平面，所述平面为基准平面；

计算树干在高度h处的生长方向其中，协方差矩阵的最大特征值对应的特征向量为树干在高度h处的生长方向平面法向量赋值为进行第k+1次迭代；初次迭代时法向量为迭代结束条件为：相邻两次迭代得到的两个树干生长方向与的夹角θ_k小于0.5度或者θ_k与θ_k+1的角度差小于0.5度；迭代结束时的树干生长方向即为树干在高度h处的生长方向；以树干在高度h处的生长方向为平面的法向量，通过定位点构建树干在高度h处的树干横断面；

步骤4，获取用于计算树干直径的树干点云及其凸包点；

步骤5，构建经过凸包点的一条闭合三次B样条曲线以模拟围尺轨迹，提取树干直径。

2.根据权利要求1所述的一种从地面激光雷达点云数据中提取树干直径的方法，其特征在于：所述步骤1，从地面三维激光扫描仪获取的树木点云中提取树干点云，包括：

使用地面三维激光扫描仪配备的点云处理软件对地面三维激光扫描仪在多站扫描的点云进行配准得到单木的点云数据，经移除地表点云、移除树枝点云、移除离散噪声点后得到树干点云。

3.根据权利要求1所述的一种从地面激光雷达点云数据中提取树干直径的方法，其特征在于：所述步骤2中，遍历树干点云集得到z坐标值最小的一个点，并将得到的所述点设为树干点云的最低点。

4.根据权利要求1或3所述的一种从地面激光雷达点云数据中提取树干直径的方法，其特征在于：所述预设距离为0.5厘米。

5.根据权利要求1所述的一种从地面激光雷达点云数据中提取树干直径的方法，其特征在于：所述步骤3中，计算树干在高度h处的生长方向包括：

以基准平面为基础，分别构建5个与此平行的平面，其中3个平面位于基准平面上方，2个平面位于基准平面下方，且相邻平面间的距离为0.5厘米；

得到5个树干点云块与5个几何中心点；

对5个几何中心点采用主成分分析方法计算树干在高度h处的生长方向

6.根据权利要求1所述的一种从地面激光雷达点云数据中提取树干直径的方法，其特征在于：所述步骤4，获取用于计算树干直径的树干点云及其凸包点，包括：

在树干横断面上方构建一个与树干横断面平行且距离为1厘米的平面，位于树干横断面与平面间的树干点云即为用于计算树干直径的树干点云，将此部分点云投影至树干横断面上得到一个三维空间上的平面点云集，采用几何变换中的旋转操作，将此平面点云集旋转为与水平面平行的一个二维空间上的平面点云集，根据凸包算法计算此平面点云集的凸包点集Q_Dhull。

7.根据权利要求6所述的一种从地面激光雷达点云数据中提取树干直径的方法，其特征在于：所述步骤5，构建经过凸包点的一条闭合三次B样条曲线以模拟围尺轨迹，提取树干直径，包括：

将Q_Dhull首部的3个点依次加入到Q_Dhull的尾部，同时也将尾部的3个点依次加入到Q_Dhull的首部构成插值点Q＝{Q_k},k＝0,1,...,n；

使用向心参数法计算插值点Q的节点参数值令d为总弦长，参数值定义为：

${\stackrel{\‾}{u}}_0=0,{\stackrel{\‾}{u}}_n=1,{\stackrel{\‾}{u}}_k={\stackrel{\‾}{u}}_{k-1}+\sqrt{|Q_k-Q_{k-1}|}/d,k=1,2,...,n-1;$

再根据节点参数值使用平均值法计算节点向量U＝{u₀,u₁,...,u_m}，有

$u_0=...=u_p=0,u_{m-p}=...=u_m=1,u_{j+p}=\frac{1}{p}\underset{i=j}{\overset{j+p-1}{\Σ}}\stackrel{\‾}{u_i},j=1,2,...,n-p;$

其中，p＝3表示样条曲线的次数，n＝|Q|表示插值点的个数，m＝n+p-1；然后采用全局曲线插值的方法求解非有理三次B样条曲线的控制点，根据等式关系：

其中，N_i,p(u)表示为：

得到一个系数矩阵为(n+1)×(n+1)的线性方程组，求解方程组得到n+1个控制点{P_i}，根据非有理三次B样条曲线的表达式：

$C\left(u\right)=\underset{i=0}{\overset{n}{\Σ}}{N}_{i,p}\left(u\right)P_i,p=3;$

得到一条连续的曲线，该曲线会在Q_Dhull的第一个节点Q₀附近相交；设置Q_Dhull的第一个节点Q₀为首节点，从首节点处开始构建，当曲线环绕一周到达首节点时停止曲线的构建即可得到完整的闭合曲线，此条闭合曲线就是模拟围尺测量轨迹的曲线；非有理三次B样条曲线是一条参数曲线，可表示为：

$s\left(u\right)=\left\{\begin{array}{c}{x}_i=x\left(u\right)\\ y_i=y\left(u\right)\\ z_i=z\left(u\right)\end{array}\right.;$

根据弧长公式，对于给定的参数区间u∈α,β，曲线在此区间的弧长s计算公式可表示为：

$s=f_{\mathrm{\α}}^{\mathrm{\β}}\left(u\right)du=\sqrt{x{\left(u\right)}^2+y{\left(u\right)}^2+z{\left(u\right)}^2}_{\mathrm{\α}}^{\mathrm{\β}}du;$

使用数值积分中的复合辛普森求积公式方法计算得到闭合曲线的长度L，则树干直径可表示为L/π。