基于TLS数据构建落叶松树干削度方程

康庆江 　李丹丹 　贾炜玮

摘要：  利用地基激光雷达获取数据，构建TLS数据代替伐倒木实测数据的树干削度模型，采用地基激光雷达扫描林口林业局6块落叶松样地和18株解析木。首先对TLS数据提取的解析木不同高度处的直径进行精度分析，再分别采用解析木TLS数据与实测数据建立树干削度方程，利用R软件拟合5个基础模型，采用Bias、RMSE、R2、P%对削度方程进行评价。结果表明：地基激光雷达获取的树干直径精度达到98%以上;TLS数据与实测数据构建的削度方程拟合结果基本相同，最优模型都为Kozak （2002）-II方程。采用TLS获取数据具有较高的精度，0.7H为直径提取最适高度，模型拟合效果较好，可以利用TLS获取树干直径代替伐倒解析木量取数据。

关键词：  TLS;  相对高度处直径;  削度方程

中图分类号：   S 758               文献标识码：   A                文章编号：1001 - 9499（2022）03 - 0041 - 06

削度方程在樹干材积估算、树干三维空间模型重建、出材率估计以及森林规划、仿真与优化造材中有着广泛的应用[ 1 ]。利用削度方程能够估计树干上任意高度处的直径及任意直径处的高度、树干总材积、不同规格的商品材积[ 2 ]。使用简化的Kozak可变指数削度模型，通过混合效应方法添加辅助直径，能够提高模型的预测能力[ 3 ]。Doyog等人对兼容的Max and Burkhart （1976）削度方程的研究表明：树干底部直径的估计值高于真实值，而中部和上部的直径则被低估;特定直径下的整个树干高度和体积被高估[ 4 ]。使用削度方程与总材积方程和材积销售方程组成的兼容系统进行建模，最优模型的拟合精度能达到98%[ 5 - 7 ]。姜立春等采用非线性混合模型拟合树冠特征的兴安落叶松树立削度方程关系，结果表明，树干削度变化与冠长率有关，冠长率越大，干形质量越差[ 8 ]。蔡坚等[ 9 ]提出，中度间伐（伐去35.7%的植株） 后的林木干形相对较饱满。

激光雷达是采用激光发射机对目标物进行主动扫描，根据发射的激光束来探测目标物的位置、方向、形状、高度、速度等参数的雷达扫描系统，广泛应用于环境监测、土地利用、农业开发、资源勘探、水利工程、交通通讯、城市规划、防震减灾等方面。地基激光雷达于2000年开始引入[ 11 ]，与传统的测绘技术相比，TLS能够获取相对完整的目标三维坐标信息，这些信息能够以点云的形式进行重建，从而真实还原目标的整体结构和形态特征[ 12 ]。TLS获取的点云数据密度大、精度高，并具有自动处理数据的潜力，且地基激光雷达扫描林木冠层的分辨率高，对林木没有破坏性，能够自动获取林木参数，重建真实的森林结构三维场景[ 13 ]。

地面激光扫描技术通过非破坏性方法准确快速地获取立木的三维数据，具有很高发展的潜力，对于森林的可持续发展具有非常大的研究意义。已有研究表明，利用地基激光雷达提取树木的胸径、树高都达到较高的精度，但对于树干上部的直径提取效果却鲜见报道，地基激光雷达测树技术在大尺度的数据获取中仍然存在数据精度和结果验证的问题[ 15 ]。为了检验地基激光雷达在整个树干上的提取精度，本文利用提取直径数据建立削度方程，旨在为估算林木材积提供新的数据获取方式，同时使地基激光雷达在林业中发挥更大的作用。

1 数据来源

2019年4月，在黑龙江省牡丹江市林口林业局根据不同立地条件、不同密度、不同年龄、不同抚育间伐强度选择6块落叶松固定样地，测定样地内每株树的胸径、树高、冠幅及相对坐标。在6块样地附近找3株平均木进行树干解析，共计选取18株解析木。根据实测总树高计算解析木各相对高度0 m、0.1H、0.2H、0.3H、0.4H、0.5H、0.6H、0.7H、0.8H、0.9H，测量各相对高处的直径。

TLS数据在伐倒解析木前获取，采用Trimble TX8窄红外激光束进行扫描。每块样地扫描5站，每株解析木扫描3站，每站扫描时间为3 min。数据预处理后6块样地共获得419株落叶松样木，提取2 784个不同高度处的直径;18株解析木共提取到270个直径。将所提取的直径按75%和25%的比例分为建模样本和检验样本，样地测树因子统计量和实测解析木因子统计量分别见表1和表2。

2 研究方法

2. 1 点云数据处理

扫描后的原始数据为TZF格式，采用Trimble RealWorks11.1和LiDAR360软件对原始数据进行预处理;再根据树高使用最小二乘法拟合不同高度处的直径。为了减少树干高度对提取精度的影响，本文用相对高度代替绝对高度，在提取直径时，切片厚度为10 cm，如提取1.3 m处的直径则选择1.25～1.35 m树干高度的点云进行拟合。

选用伐倒的18株解析木各相对高度处的实测直径数据进行验证，也为构建削度方程所需数据的获取提供依据。将18株解析的各相对高度处直径作为实测值，地基激光雷达提取的同株样木相同位置处的直径为预测值，采用确定系数、均方根误差、绝对误差、提取精度评价每株树的提取精度。评价指标如下：

式中，n为每棵树提取的直径数量;xi为提取单木参数;xi为与提取出的单木对应的实测单木参数。

点云数据在建立数字高程模型时对地面与树基的识别能力有限，且实测树基的直径时差异也比较大。因此，单独将0 m处的直径与实测作比较，并对18株解析木的胸径、树高及所有提取直径进行分析。解析木树高在0.7H以上时，树干上部分点云由于树枝树叶的遮挡，数据有缺失，无法提取到直径。

由表3可知，10个不同相对高度处的直径提取的R2大部分都在0.92以上， R2最大值为0.995 4（0.1H处）、最小值为0.781 2（0.9H处）。均方根误差（RMSE）最小值为0.31 cm（0.1H），最大值为0.64 cm（0.8H）。绝对误差（Bias）为0.25～0.50 cm。胸径和树高的提取效果都较好，胸径提取精度为98.18%，树高提取精度达到99.36%。枝下高的精度为94.60%，冠幅的提取效果不是很好。

18株解析木的10个不同相对高位置直径提取精度都在85%以上，最高达98%以上。精度最高的是0.1H（98.43%），最低的是0.9H（85.51%）。直径提取精度的变化趋势为：在0.1H时最大值;随着高度的增加，提取精度逐渐减小。在树冠处，由于地基激光雷达扫描时，树干上部晃动，点云数据处理时多站数据不能完全匹配，在直径提取时造成较大视觉误差;到0.9H时，点云数据缺失、密度减小，直径提取误差较大。

2. 2 削度方程备选模型

参考最新版《测树学》（第四版）教材及国内外相关文献，本文选用5个常用的不同类型的落叶松树干削度方程作为备选模型，分别为“Kozak 等（1969）”和“Schumacher（1973）”;一个带有两个拐点的分段削度方程“Max and Burkhart （1976）”;2个可变指数削度方程，分别为“曾伟生，廖志云（1997）”和“Kozak （2002）-II”。

分段削度方程采用姜立春等人[ 12 ]研究兴安落叶松使用的Max and Burkhart （1976）削度方程，因为所建模型的树种同为落叶松，因此树干上两个拐点的相对高也采用姜立春等人拟合的兴安落叶松的拐点相对高度值，这样能消除树干干形的影响，则下部和上部拐点处的相对高度值分别为0.09H和0.77H。模型形式如下：

（1）简单削度方程：

（模型1） Kozak 等（1969）：

d2/D2=a1+a2（h/H）+a3（h/H）2（5）

（模型2） Schumacher（1973）：

d=aD（h-H）/H（6）

（2）分段削度方程：

（模型3） Max and Burkhart （1976）：

d2=D2[a1（q-1）+a2（q2-1）+a3（b1-q）2I1+a4（b2-q）2I2]

I1=1，if q≤b1;0 otehrwise

I2=1，if q≤b2;0 otehrwise（7）

（3）可變指数削度方程：

（模型4） 曾伟生，廖志云（1997）：

d/D=

K=a+a（h/H）+a（h/H）+a（D/H）（8）

（模型5） Kozak （2002）-II：

d=aDH·

 （9）

式中，D为树干胸径（cm）;d为树干h高处的直径（cm）;ai为模型参数;i=1，2，3，4，5，6，7，8，9;H为全树高（m）;h为从地面起算的高度或至某上部直径限或利用长度限处的高度（m）;q=h/H;t=1.3/H;b1、b2为树干下部拐点和上部拐点处的相对高度，当q≤b1时，I1=1;当q>b1时，I1=0;当q≤b2时，I2=1;当q>b2时，I2=0。

2. 3 模型评价与检验

用R软件对模型进行拟合，得到模型参数，模型拟合与检验结果采用R2、RMSE、Bias、预估精（P%）进行评价，R2、RMSE、Bias见式（1）、（2）、（3），P%公式如下：

预估精度 （P%）：

式中，n为样本数; p为模型参数的数量;xi为提取单木直径;xi为与提取出的单木对应的实测单木直径。

从各评价指标可以看出（表4），两种数据都是模型（5）的拟合精度最高;模型（1）、（2）的点云数据拟合效果优于实测数据;模型（3）、（4）、（5）的实测数据拟合效果优于点云数据。总体而言，点云数据的拟合效果与实测数据相差极小，R2都在0.94以上，最优模型的R2仅相差0.002 4，Bias、RMSE也非常接近。

由表5可知，两种数据检验的R2都在0.92以上，除了模型（4），其余模型点云数据的R2都大于实测数据;点云数据削度模型的预估精度达到99%以上，5个模型的检验效果与拟合效果大小一致。实测数据与TLS数据中都是模型（5）的R2最大，Bias、RMSE最小，且精度最高，因此模型（5）为两种数据的最优削度模型。

为了更清楚地描述TLS数据建立削度方程的效果，使用样地单木所提取的胸径、树高平均值作为模拟树木的胸径（16.98 cm）和树高（18.44 m），采用模型（5）绘制两种数据的干形变化趋势图（图1）。TLS削度模型与实测削度模型基本重合。根据解析木的拟合结果与检验结果可知，TLS数据与实测数据的结果一致、统计指标相差极微小，且精度较高，两种数据模拟的树干干形变化趋势基本重合，因此可以采用TLS数据代替伐倒木获取数据建立削度方程。

2. 4 不同相对高度削度模型的拟合与检验

由于地基激光雷达的扫描效果易受树木枝叶的遮挡，因此树干上部的结构信息提取精度相对偏低。直径提取到0.5H时的精度为96.78%，越往上精度越低，点云数据的精度会影响削度模型的拟合效果，因此本研究根据前文的提取精度，分别采用0 m到相对高度为0.5H、0.6H、0.7H、0.8H及0.9H的直径构建削度方程。备选模型采用前文拟合的最优削度方程（模型5：Kozak （2002）-II），采用同样的方法进行模型的拟合与检验。

由表6可知，采用0 m～0.6H的直径建立的削度模型精度最高（99.16%），比使用全部点云数据的建模精度高0.15%。0 m～0.5H时，直径提取的精度虽高，但表示树干上部干形的数据较少;0 m～0.7H、0.8H、0.9H时，模型的精度与直径提取的精度变化趋势一致，都是随着树干高度的增加精度逐渐降低。

结合前文中直径的提取精度与削度模型的拟合精度，直径提取到0.7H为地基激光雷达获取树干直径数据的最佳高度。0.7H时的直径提取精度为94.80%，与0.6H的精度相差不大，而比0.8H高出很多;0.7H～0.9H时，数据量越多模型精度越低，0.7H时的削度模型精度达到99.12%;也减少了大量的点云数据处理工作。

2. 5 样地削度模型的构建

本研究扫描6块样地，共正确分割出419株落叶松单木，采用样地单木直径构建的削度方程更能代表大尺度林分内树干的干形变化，也能较准确的预测树干上部直径、估算树干材积。根据前文的研究结果，本研究使用模型5（Kozak （2002）-II）、0 m～0.7H处的直径建立样地削度方程，模型参数及统计量如表7。t检验表明，模型参数估计都显著（P<0.000 1），拟合统计量（R2=0.965 4）表明该模型能很好地描述落叶松树干干形的变化趋势。

本研究采用未参与建模的检验数据进行模型的独立性检验。检验的R2、RMSE、Bias和精度分别为0.958 9、1.223 4 cm、0.802 8 cm和99.71%。说明使用地基激光雷达获取数据建立的削度模型满足精度要求，采用地基激光雷达获取的数据在林业中具有重要的实用价值，可以在林业上大力推广使用。

3 讨 论

3. 1 传统的构建削度方程都是通过砍伐树木的方式获取数据，这样就会造成对森林破坏和资源浪费，且费时费力。本研究提供一种新的非破坏性的方法，即激光雷达扫描法获取数据，既有实测解析木数据作为验证，又有大尺度的样地数据。本研究提取相对高度的直径，既减少了树干高度对提取精度的影响，又能通过削度模型预测树干上部的直径。

3. 2 根据解析木不同相对高的直径提取精度得知，整株树干从下往上，直径提取精度逐渐减小。研究发现，在地基扫描时因为风吹树干的影响，或者在测站配准时存在一定的误差，不同测站上部树干存在夾角或重影，在0.7H以上的直径提取时易造成较大的视觉误差;而0.8H时，多站树干交叉角较大，彼此分开，视觉误差相对降低;到0.9H时，点云密度减小，直径模型不易判读，提取误差增大。此外，林分密度、郁闭度、树高、枝生长状况以及扫描的测站数及测站之间的距离都会影响地基激光雷达提取信息的精度。

3. 3 国际上对基于TLS点云数据的树干检测率、树干位置、胸径、树高、树干体积、生物量等森林测量过程已经达到了高度的自动化水平[ 17]。Cabo等采用算法自动获取的胸径和树高与传统的树木测量方法相比，胸径差异小于1 cm的占85%， 树高差异小于0.5m的占92%[ 18 ]。本研究采用手动提取的方式，除0 m以外，所有的直径误差都小于1 cm，从数据质量上来说，本研究的精度更高，但数据处理时间却较长，工作量也更大。本文采用不同高度的直径拟合基础削度方程，得出对于TLS数据最适合的获取高度为0.7H，这个位置的数据具有较高的精度又能减少冗杂数据的提取时间。

3. 4 本研究另一特点是完全利用TLS技术获取树木的直径、树高、活枝高和冠幅数据建立削度方程模型，这一方法无需砍伐树木又能提高模型质量。非破坏性的测树方法还可采用RD1000激光测树仪，其测量直径在相对高度约0.64～0.8时的突变最大，与本文中采用TLS数据得到的结果相似（0.7）[ 4 ]。Sun等利用改进的Schumacher方程拟合地基激光雷达扫描的杨树削度方程，其拟合R2=0.96，与本研究基础模型的R2基本相同[ 19 ]。也有研究表明，作为计算材积的变量，使用树干横截面积比直径的预测效果更好[ 20 ]。

4 结 论

直径提取的最适高度到0.7H，且此时的精度满足林业调查的要求。若提取的高度过高，树干上部直径精度过低，则整体精度变低，因此会得到大量冗余数据;若提取的高度过低，则获取的数据量过少，会影响树干的整体信息获取效果。并且，获取的数据量和精度都会影响削度模型的拟合，进而影响树干干形的模拟和材积估算精度。

本研究采用标准解析木的点云数据与实测数据对比，结果表明TLS数据的精度满足林业调查要求，使用两种数据建立的削度方程模型拟合效果较好、精度较高。采用0 m～0.7H处的直径数据拟合模型的精度最高。采用激光雷达数据建立削度方程具有典型的代表性和超强的实用性，在今后的外业测量工作中，可通过本方式无须砍伐树木就能进行高精度的森林调查。

参考文献

[1] 孟宪宇.  削度方程和出材率表的研究[J].  南京林业大学学报：自然科学版， 1982（1）： 122-126.

[2] 姜立春，  刘瑞龙.  基于非线性混合模型的落叶松树干削度模型[J].  林业科学， 2011， 47（4）： 101 - 106.

[3] Geneviève Lejeune， Ung C H ， Fortin M ， et al. A simple stem taper model with mixed effects for boreal black spruce[J]. european journal of forest research， 2009， 128（5）：505-513.

[4] Doyog N D， Lee Y J， Sun J L， et al. Compatible taper and stem volume equations for Larix kaempferi（Japanese larch） species of South Korea[J]. 山地科学学报（英文）， 2017， （7）： 109 - 117.

[5] Menéndez-Miguélez M， Canga E， álvarez-álvarez P， et al. Stem taper function for sweet chestnut （Castanea sativaMill.） coppice stands in northwest Spain[J]. Annals of Forest Science， 2014， 71（7）： 761 - 770.

[6] Sharma M， Parton J. Modeling Stand Density Effects on Taper for Jack Pine and Black Spruce Plantations Using Dimensional Analysis[J]. Forest Science， 2009， 55（3）： 268 - 282.

[7] 唐誠.  西南桦人工林生长模拟及立地质量评价[D]. 中国林业科学研究院， 2017.

[8] 姜立春， 蒋雨航. 利用混合模型模拟树冠特征对兴安落叶松树干干形的影响[J].  北京林业大学学报， 2014， 36（2）： 10 -14.

[9] 蔡坚，  潘文，  王保华，  等.  林分密度对湿地松林木干形影响的研究[J]. 林业与环境科学， 2006， 22（2）： 6 - 10.

[10] Arp H， Griesbach J C， Burns J P. Mapping in tropical forests： a new approach using the Laser APR[J]. Photogrammetric Engine-

ering & Remote Sensing， 1982， 48（1）： 91 - 100.

[11] 王俞明.  基于地基激光雷达的杉木参数提取与材积估测[D]. 长沙： 中南林业科技大学， 2019.

[12] Henning J G. Detailed stem measurements of standing trees from ground-based scanning lidar[J]. Forest Science， 2006， 52（1）： 67 - 80（14）.

[13] 刘鲁霞，  庞勇.  机载激光雷达和地基激光雷达林业应用现状[J].  世界林业研究， 2014， 27（1）： 49 - 56.

[14] 刘鲁霞，  庞勇，  李增元，  等. 用地基激光雷达提取单木结构参数——以白皮松为例[J]. 遥感学报， 2014， 18（02）： 365 - 377.

[15] 邢万里.  林地TLS点云数据自动配准算法研究[D]. 哈尔滨： 东北林业大学， 2018.

[16] 姜立春，  李凤日，  刘瑞龙.  兴安落叶松树干削度和材积相容模型[J]. 北京林业大学学报， 2011， 33（5）： 1 - 7.

[17] Yuan S， Xinlian L， Ziyu L， et al. Deriving Merchantable Volume in Poplar through a Localized Tapering Function from Non-  Destructive Terrestrial Laser Scanning[J]. Forests （19994907）， 2016， 7（4）： 87.

[18] Cabo C， Ordóez， Celestino， López-Sánchez， Carlos A， et al. Automatic dendrometry： Tree detection， tree height and diameter estimation using terrestrial laser scanning[J]. International Journal of Applied Earth Observation and Geoinformation， 2018， 69： 164 - 174.

[19] Marchi M， Scotti R， Rinaldini G， et al. Taper Function for Pinus nigra in Central Italy： Is a More Complex Computational System Required[J]. Forests， 2020， 11（4）： 405.

[20] Gregoire， T.G.; Schabenberger， O.; Kong， F. Prediction from an integrated regression equation： A Forestry application[J]. Biome-  trics 2000， 56， 414–419.

第1作者简介：  康庆江（1964-），  男，  高级工程师，研究方向：  森林经理。

通讯作者：  贾炜玮（1978-），  男，  教授，  研究方向：  林分生长与收获模型。

收稿日期： 2022 - 02 -  20