CN106931883B - 一种基于激光点云数据的树干材积获取方法 - Google Patents

Abstract

本发明公开提供了一种基于激光点云数据的树干材积获取方法，包括：步骤1，对树干点云垂直分段；步骤2，逐一计算各个垂直分段的体积；步骤3，计算所有垂直分段体积的累加和。本发明在获取树干断面轮廓的基础上，通过准确计算树干断面的断面积来准确计算树干每段的树干材积，具有较高的精确性，并为准确提取树干生物量、碳储量与树干量材优化等提供基础数据。

Description

一种基于激光点云数据的树干材积获取方法

技术领域

本发明涉及一种树干参数提取技术，尤其是涉及一种基于激光点云数据的树干材积获取方法。

背景技术

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

树干材积是测树学中的基本测树因子，是树木经济价值的直接体现，也是计算树木生物量与碳储量的基础。传统上将树干简化为一个绕中心轴旋转的旋转体，以此为基础估算树干材积。然而树干形状是不规则的，旋转体的简化结构过于简化了树干的几何形状，在此基础上计算得到的树干材积必定与树干的真实材积存在着差异，也正因如此，如何准确计算树干材积是林业工作者一直关注的问题。

发明内容

有鉴于此，本发明的目的是针对现有技术的不足，提供一种基于激光点云数据的树干材积获取方法，用以根据树干点云获取精确的树干材积。

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

一种基于激光点云数据的树干材积获取方法，包括如下步骤：

步骤1，对树干点云垂直分段；

步骤2，逐一计算各个垂直分段的体积；

步骤3，计算所有垂直分段体积的累加和。

优选的，所述步骤1，对树干点云垂直分段，包括：

以高度为标准，计算树干点云中的最低点p_lowest，其中z_lowest为最低点p_lowest的高度值，并以所述最低点p_lowest为基础计算树干的高度；

构建平面π₀与π₁，其方程分别为z＝z_lowest与z＝z_lowest+thickness，其中thickness是每一个垂直分段的高度；则位于平面π₀与π₁间的点云形成树干点云的一个垂直分段Q₀；继续向上构建平面π₂，其方程为z＝z_lowest+thickness*2，则位于平面π₁与π₂间的树干点云形成另一个垂直分段Q₁；则，平面π_i的方程可表示为z＝z_lowest+thickness*i，位于平面π_i与π_i+1间的树干点云形成点云垂直分段Q_i；其中，平面π_i与平面π_i+1分别称为垂直分段Q_i的下平面与上平面；按照上述垂直分段的构建方法直至树梢部位得到若干个树干点云的垂直分段。

优选的，所述步骤2，逐一计算各个垂直分段的体积，包括：

计算垂直分段在下平面的投影点集；

计算垂直分段的断面轮廓点；

构建断面轮廓曲线并计算树干断面的断面积；

计算垂直分段的体积。

优选的，所述计算垂直分段在下平面的投影点集，包括：

计算垂直分段Q_i在其下平面的投影点集Q_i′；将垂直分段Q_i中的树干点云投影至垂直分段Q_i的下平面π_i得到一个平面点集Q_i′；计算平面点集Q_i′的凸包多边形的质心点作为平面点集Q_i′的几何中心点c_i。

优选的，所述计算垂直分段的断面轮廓点，包括：

设点p_m为垂直分段Q_i中的一个点，p_m′表示点p_m在垂直分段Q_i下平面的投影，则几何中心点c_i与投影点p_m′构成的向量与X轴的夹角可表示为：

根据角度分区的角度值θ，将平面点集Q_i′中的点分为个角度分区，每个角度分区Q_ij中的点可表示为：

计算每一个角度分区Q_ij中点的重心点，所述重心点为所述角度分区中所有点的几何平均值；并将所述重心点作为所述角度分区的轮廓点，由此得到垂直分段Q_i的个轮廓点。

优选的，所述构建断面轮廓曲线并计算树干断面的断面积，包括：设垂直分段Q_i的断面轮廓点集合用R_profile表示，首先采用非有理三次B样条曲线构建一条插值于点集R_profile的闭合曲线，所述闭合曲线即为断面轮廓曲线；断面轮廓曲线所围区域的面积即为树干断面的断面积。

优选的，所述构建断面轮廓曲线并计算树干断面的断面积，包括：

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

使用弦长参数法计算插值点R的节点参数值

令d为总弦长，节点参数值定义为：

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

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

采用全局曲线插值的方法求解非有理三次B样条曲线的控制点，得到系数矩阵为(n+1)×(n+1)的线性方程组：

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

求解方程组得到n+1个控制点{P_i}，根据非有理三次B样条曲线的表达式

p＝3；可得到一条连续曲线，设置R_profile的第一个节点R₀为首节点，从首节点处开始构建，当曲线环绕一周到达首节点时停止曲线的构建，即可得到完整的闭合曲线，所述的闭合曲线即为断面轮廓曲线；

断面轮廓曲线是一条闭合曲线，根据格林公式将该闭合曲线所围区域的面积转换为第一类曲线积分的计算，从而得到树干断面的断面积S_i。

优选的，所述计算垂直分段的体积，包括：

计算垂直分段Q_i的体积V_i，将垂直分段近似为一个柱体，断面形状为断面轮廓曲线所围区域的形状，则根据柱体体积公式，垂直分段Q_i的体积V_i可表示为V_i＝S_i*thickness。

优选的，所述步骤3，计算所有垂直分段体积的累加和，包括：所有垂直分段体积的累加和即为树干材积，树干材积的计算公式可表示为V＝thickness*∑S_i，V为树干材积。

优选的，所述垂直分段的高度值thickness为0.5厘米。

本发明的有益效果是：

传统上使用区分求积法计算伐倒木的树干材积，其建立在将树干简化为树干横断面是圆，树干中心轴为直线的旋转体的基础上，该方法将树干区分成若干等长或不等长的区分段，使各区段的干形更接近于正几何体，再根据近似求积式计算各分段材积，合计各分段材积得到树干材积。计算立木的树干材积的传统方法是采用经验公式，将立木高度与直径的测量值代入材积经验公式得到立木的树干材积。由于材积经验公式建立在伐倒木材积计算的基础上，因此伐倒木与立木的材积计算都建立在树干横断面为圆，树干中心轴为直线的规则旋转体的几何理论基础之上。然而树干形状是不规则的，树干横断面也不是规则的圆形，因此在此简化基础上计算的材积与树干的真实材积存在着差异。

本发明建立在精确计算树干断面积的基础上，通过精确获取树干断面的轮廓点，进而构建一条闭合的断面轮廓曲线，闭合曲线所围面积即为断面面积，该段断面积计算方式建立在获取树干断面轮廓的基础上，贴近于树干断面的真实轮廓；而且使用的thickness值为0.5厘米，对于一段长度为0.5厘米的树干，其在各处的断面几乎相同，因此本发明计算的垂直分段的体积非常接近于垂直分段的真实体积。相对比传统上伐倒木和立木的不同材积计算方法，本发明在获取树干点云的基础上，无需区分伐倒木还是立木，均可使用本发明的方法计算树干材积。

本发明在精确计算树干断面断面积的基础上精确计算树干材积，相比于传统的将树干断面简化为圆来计算树干断面的断面积的方法，本发明在获取树干断面轮廓的基础上计算树干断面的断面积，具有较高的精确性，从而为准确提取树干生物量、碳储量与树干量材优化等与树干材积的相关研究与应用提供基础数据。

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

附图说明

图1是本发明一种基于激光点云数据的树干材积获取方法的流程图。

图2是本发明一种基于激光点云数据的树干材积获取方法步骤2的流程图。

图3是树干点云垂直分段的示例图。

图4是断面轮廓曲线的示例图。

具体实施方式

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

如图1所示，本发明提供一种基于激光点云数据的树干材积获取方法，包括如下步骤：

步骤1，对树干点云垂直分段；

步骤2，逐一计算各个垂直分段的体积；

步骤3，计算所有垂直分段体积的累加和。

步骤1，对树干点云垂直分段，包括：

根据z轴坐标，也即高度值，计算树干点云中的最低点p_lowest，设最低点p_lowest的Z轴坐标为z_lowest，并根据最低点为基础计算树干的高度；

构建平面π₀与π₁，其方程分别为z＝z_lowest与z＝z_lowest+thickness，其中thickness是每一个垂直分段的高度值。则位于π₀与π₁间的点云形成树干点云的一个垂直分段Q₀；继续向上构建平面π₂，其方程为z＝z_lowest+thickness*2，则位于平面π₁与π₂间的树干点云形成另一个垂直分段Q₁；则，平面π_i的方程可表示为z＝z_lowest+thickness*i，位于平面π_i与π_i+1间的树干点云形成点云垂直分段Q_i；其中，平面π_i与平面π_i+1分别称为垂直分段Q_i的下平面与上平面；按照上述垂直分段的构建方法直至树梢部位得到若干个树干点云的垂直分段；如图3所示，图中上平面1与上平面2之间是垂直分段Q_i。

thickness的值可根据实际需求预先设定，其值的大小可调，其值越小越好，但考虑到点云密度的影响，在一个实施例中，thickness的值为0.5厘米。thickness的值越小，计算得到的体积越精确，但受点云密度的影响，若thickness的值太小，则垂直分段中包括的点云太少，容易导致计算的断面轮廓点无法反映断面的真实轮廓，若thickness的值太大，树干在此段范围内干形有较大变化，则导致体积计算精度降低。

如图2所示，步骤2，逐一计算垂直分段Q_i的体积V_i，包括如下步骤：

步骤201，计算垂直分段在下平面的投影点集；

步骤202，计算垂直分段的断面轮廓点；

步骤203，构建断面轮廓曲线并计算树干断面的断面积；

步骤204，计算垂直分段的体积。

对于一个特定的垂直分段Q_i，将Q_i中的树干点云投影至垂直分段的下平面得到平面投影点集Q_i′，计算平面点集Q_i′构成的凸包多边形的质心点c_i，根据点c_i和Q_i′中点构成的向量与角度分区的角度值θ度对投影点集Q_i′进行角度分区，再计算每一个角度分区的重心点，个角度分区的重心点构成断面的轮廓点R_profile；再采用非有理三次B样条曲线构建一条过点集R_profile的闭合曲线，该曲线为断面轮廓曲线；然后再使用格林公式计算闭合曲线所围区域的面积的积分转换为第一类曲线积分，计算得到树干断面的断面积S_i；垂直分段Q_i的体积V_i可表示为V_i＝thickness*S_i平方厘米。

在一个实施例中，θ值为1度(θ的值越小，断面轮廓点越多，断面轮廓会越准确，但受点云密度的影响，当θ太小，容易导致角度分区中无点云，当θ值为1度时，可确保断面在360个方位上都有轮廓点，基本上能真实反映树干断面的轮廓特征)。该实施例中，采用平面点集凸包多边形的质心点作为平面点集的几何中心点；根据角度分区的角度值与平面点集的几何中心点对平面点集进行角度分区的方法。

步骤2具体包括如下步骤：

计算垂直分段Q_i在其下平面的投影点集Q_i′，将垂直分段Q_i中的树干点云投影至垂直分段Q_i的下平面π_i得到一个平面点集Q_i′；计算平面点集Q_i′的凸包多边形的质心点作为平面点集Q_i′的几何中心点c_i；

设点p_m为垂直分段Q_i中的一个点，p_m′表示点p_m在垂直分段Q_i下平面的投影，则几何中心点c_i与投影点p_m′构成的向量与X轴的夹角可表示为：

根据角度分区的角度值θ，将平面点集Q_i′中的点分为个角度分区，每个角度分区Q_ij中的点可表示为：

p_m′∈Q_i′}，计算每一个角度分区Q_ij中点的重心点(该角度分区中所有点的几何平均值)，并将该重心点作为该角度分区的轮廓点，由此得到垂直分段Q_i的个轮廓点；

构建断面轮廓曲线并计算树干断面的断面积，设垂直分段Q_i的断面轮廓点集合用R_profile表示，首先采用非有理三次B样条曲线构建一条插值于点集R_profile的闭合曲线，该闭合曲线即为断面轮廓曲线；断面轮廓曲线所围区域的面积即为树干断面的断面积。具体方法为：

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

使用弦长参数法计算插值点R的节点参数值

令d为总弦长，参数值定义为：

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

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

采用全局曲线插值的方法求解非有理三次B样条曲线的控制点，得到系数矩阵为(n+1)×(n+1)的线性方程组：

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

求解方程组得到n+1个控制点{P_i}，根据非有理三次B样条曲线的表达式p＝3；可得到一条连续曲线，设置R_pro_file的第一个节点R₀为首节点，从首节点处开始构建，当曲线环绕一周到达首节点时停止曲线的构建，即可得到完整的闭合曲线，所述的闭合曲线即为断面轮廓曲线；图4断面轮廓曲线的示例图；

断面轮廓曲线是一条闭合曲线，根据格林公式可将该闭合曲线所围区域的面积转换为第一类曲线积分的计算，从而得到树干断面的断面积S_i；

计算垂直分段Q_i的体积V_i，将垂直分段近似为一个柱体，断面形状为断面轮廓曲线所围区域的形状，则根据柱体体积公式，垂直分段Q_i的体积V_i可表示为V_i＝S_i*thickness。

步骤3，计算所有垂直分段体积的累积和，包括：

合计步骤2计算的每一个垂直分段的体积即得到树干的材积，其公式为V＝thickness*∑S_i。在一个实施例中，thickness的值为0.5厘米，则公式可表示为V＝0.5*∑S_i。

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

Claims (6)

1.一种基于激光点云数据的树干材积的计算方法，其特征在于：包括如下步骤：

步骤1，对树干点云垂直分段；

步骤2，逐一计算各个垂直分段的体积；

步骤3，计算所有垂直分段体积的累加和；

所述步骤1，对树干点云垂直分段，包括：

以高度为标准，计算树干点云中的最低点p_lowest，其中z_lowest为最低点p_lowest的高度值，并以所述最低点p_lowest为基础计算树干的高度；

构建平面π₀与π₁，其方程分别为z＝z_lowest与z＝z_lowest+thickness，其中thickness是每一个垂直分段的高度值；则位于平面π₀与π₁间的点云形成树干点云的一个垂直分段Q₀；继续向上构建平面π₂，其方程为z＝z_lowest+thickness*2，则位于平面π₁与π₂间的树干点云形成另一个垂直分段Q₁；则，平面π_i的方程可表示为z＝z_lowest+thickness*i，位于平面π_i与π_i+1间的树干点云形成点云垂直分段Q_i；其中，平面π_i与平面π_i+1分别称为垂直分段Q_i的下平面与上平面；按照上述垂直分段的构建方法直至树梢部位得到若干个树干点云的垂直分段；

所述步骤2，逐一计算各个垂直分段的体积，包括：

计算垂直分段在下平面的投影点集；

计算垂直分段的断面轮廓点；

构建断面轮廓曲线并计算树干断面的断面积；

计算垂直分段的体积；

所述计算垂直分段在下平面的投影点集，包括：

计算垂直分段Q_i在其下平面的投影点集Q_i′；将垂直分段Q_i中的树干点云投影至垂直分段Q_i的下平面π_i得到一个平面点集Q_i′；计算平面点集Q_i′的凸包多边形的质心点作为平面点集Q_i′的几何中心点c_i；

所述计算垂直分段的断面轮廓点，包括：

设点p_m为垂直分段Q_i中的一个点，p_m′表示点p_m在垂直分段Q_i下平面的投影，则几何中心点c_i与投影点p_m′构成的向量与X轴的夹角可表示为：根据角度分区的角度值θ，将平面点集Q_i′中的点分为360/θ个角度分区，每个角度分区Q_ij中的点可表示为：计算每一个角度分区Q_ij中点的重心点，所述重心点为所述角度分区中所有点的几何平均值；并将所述重心点作为所述角度分区的轮廓点，由此得到垂直分段Q_i的360/θ个轮廓点。

2.根据权利要求1所述的一种基于激光点云数据的树干材积的计算方法，其特征在于：构建断面轮廓曲线并计算树干断面的断面积，包括：设垂直分段Q_i的断面轮廓点集合用R_profile表示，首先采用非有理三次B样条曲线构建一条插值于点集R_profile的闭合曲线，所述闭合曲线即为断面轮廓曲线；断面轮廓曲线所围区域的面积即为树干断面的断面积。

3.根据权利要求1所述的一种基于激光点云数据的树干材积的计算方法，其特征在于：构建断面轮廓曲线并计算树干断面的断面积，包括：

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

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

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

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

采用全局曲线插值的方法求解非有理三次B样条曲线的控制点，得到系数矩阵为(n+1)×(n+1)的线性方程组：

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

求解方程组得到n+1个控制点{P_i}，根据非有理三次B样条曲线的表达式

可得到一条连续曲线，设置R_profile的第一个节点R₀为首节点，从首节点处开始构建，当曲线环绕一周到达首节点时停止曲线的构建，即可得到完整的闭合曲线，所述的闭合曲线即为断面轮廓曲线；

断面轮廓曲线是一条闭合曲线，根据格林公式将该闭合曲线所围区域的面积转换为第一类曲线积分的计算，从而得到树干断面的断面积S_i。

4.根据权利要求1所述的一种基于激光点云数据的树干材积的计算方法，其特征在于：计算垂直分段的体积，包括：

计算垂直分段Q_i的体积V_i，将垂直分段近似为一个柱体，断面形状为断面轮廓曲线所围区域的形状，则根据柱体体积公式，垂直分段Q_i的体积V_i可表示为V_i＝S_i*thickness。

5.根据权利要求1所述的一种基于激光点云数据的树干材积的计算方法，其特征在于：所述步骤3，计算所有垂直分段体积的累加和，包括：所有垂直分段体积的累加和即为树干材积，树干材积的计算公式可表示为V＝thickness*∑S_i，V为树干材积。

6.根据权利要求1所述的一种基于激光点云数据的树干材积的计算方法，其特征在于：所述用于垂直分段的值thickness为0.5厘米。