基于2期地基激光数据的杨树干形变化监测

蒋昳萱，温小荣*，蒋佳文，朱 硕，孙 圆，翁卫松，徐 达

(1.南京林业大学 南方现代林业协同创新中心，江苏 南京 210037；2.南京林业大学 林学院，江苏 南京 210037；3.浙江省森林资源监测中心，浙江 杭州 310020)

干形用以描述树木外部形状，与树木材积和林分蓄积量息息相关[1]。传统测量通过伐倒木分段测量来获取树木任意位置处直径，计算树木干形，并基于伐倒木数据，采用统计的方法推算大范围树木干形特征[2]。受到森林采伐限额条例限制，伐倒大量树木用于科学研究已然不大现实。近年来，随着激光雷达技术的迅速发展，研究人员将其引入林业调查之中[3-6]。

地基激光雷达(terrestrial laser scanning,TLS)通过自下而上的扫描获取点云数据，可以实现林木参数的精准提取实现对树干的详细描述[7]。目前TLS数据在立木干形研究中的应用主要分为2个方面[8-9]。一是建立干曲线方程，描述树干形状。有学者探讨了不同树种的干曲线分形维数和生长趋势，提出利用TLS数据进行干形研究[10]。随后大量学者尝试利用TLS数据绘制树干曲线，分析树干形状变化。另一种则是采用干形指标描述树木干形[11-14]，其中TAP是描述树干削度的较理想的干形指标[15]。

目前，科研人员不满足于对单期点云数据的分析，着眼于树木干形的动态分析，利用多期数据分析一定期间内树木干形变化[16-17]。本研究在前人研究的基础上，提出试验形率，并以此建立试验形数形率模型，用于研究区杨树形数的推导和估算；根据杨树生长特性，对树干削度指标进行调整，同时计算研究期内各干形指标变化情况，分析了树干形状变化情况。为杨树人工林科学经营提供参考，充分发挥地面激光雷达在林业动态监测中的优势。

1 材料与方法

1.1 研究区概况

研究区位于江苏省宿迁市泗洪县陈圩林场马浪湖分场，全区面积800 km2，活立木蓄积2.3万m3，气候适宜，水热资源充沛，植物种质资源丰富。培养了美洲黑杨(Populus deltoides)中南林95杨、南林895杨、南林797杨等优良品种，是国家首批杨树良种基地、南京林业大学杨树研究开发中心试验基地(图1)。

1.2 数据来源

分别于2017年和2019年对研究区进行了2次测量。第1次测量于2017年冬季进行(T1)，并于2019年进行重复调查(T2)，2次调查遵循相同的原则进行，采用相同的扫描设备与设置。使用RIEGL VZ 400i仪器在试验区确定的杨树标准地测站点上进行扫描，每块样地的扫描不少于5站。在每块样地上选取中心点，采集由中心点发出的胸径(DBH)>5 cm、树高(h)>2 m，生长匀称、没有偏冠的点云数据进行分析(图1)。除地基激光数据外，在实地测量时，还利用直径卷尺和测高器测量了样地内树木的胸径和树高。

图1 试验区位置Fig.1 location map of test area

本次研究以南林速生系列中797杨为对象，选取2017-2019年一直存活的树木进行测算。最终选出4块不同密度样地，共计178棵树进行分析(图2)。样地中树木DBH为15.1～35.3 cm，h为17.1～27.9 m(表1)。

图2 2019年第1样地第1排示意Fig.2 Schematic diagram of row 1 of plot 1 in 2019

表1 2017年(T1)样地概况Table 1 A survey of sample plots in 2017 (T1)

1.3 数据预处理与测量方法

由于被测物表面形状、纹理，以及气候条件等因素的影响，不可避免产生噪声点和离散点，需要对数据进行相应预处理[18]：1)通过采用人机交互的方式进行点云去噪和光顺，提升点云质量。2)借助RiSCAN Pro软件实现对多站点采集的局部重叠点云自动拼接。3)通过操作软件实现地形过滤，精简点云数据以便上部直径的提取。

选用迭代K均值聚类法对林分活立木进行定位和分割，每隔1 m对相应树高处的点云进行切片，做降维处理，采用最小二乘法提取树木直径(图3)。树高则利用数字表面模型(DSM)和数字地形模型(DEM)做差，生成冠层高度模型(CHM)提取。

图3 2019年第1样地第1排第3株立木切片示意Fig.3 Slice diagram of the third tree in row 1 of plot 1 in 2019

1.4 杨树干形指标分析

经过预处理后，点云数据内包括林分内所有活立木信息，选取部分测树因子描述林木干形，分析生长变化。由于TLS数据可以提供准确的林木直径信息[19-20]，使得树木干形指标的测算变得便捷容易。参考试验形数数值稳定，探索(h+3)/2处直径在计算干形指标的作用。

通过定义试验形率(式1)，根据书斯托夫提出的模型(式2)，利用2019年点云数据建立胸高形数与试验形率关系模型计算材积。选取一定干形指标对杨树进行树木干形分析(表2)。

表2 干形指标汇总Table 2 Summary table ofstem form indicators

(1)

(2)

式中，q实为试验形率；d(h+3)/2为(h+3)/2处直径；d1.3为胸径 ；f1.3为胸高形数。

其中，重新设计单位面积TAP累积量，用以评价不同种植密度下树干削度的变化情况。单位面积TAP累积量愈大，表明TAP负生长的树木数量愈多，干形生长更加通直圆满的树木比例愈多。

对于所有干形指标中，通过测算样本树之间的变化量来评估T1与T2之间的变化，计算T1与T2之间的差值来确定它们的绝对变化；计算绝对变化量与T1时期的原始属性值的比值，确定它们的相对变化。

2 结果与分析

2.1 TLS数据提取树高胸径精度分析

研究利用胸径、树高实测值检验TLS点云数据的提取精度。将实测值与提取值进行线性拟合，利用决定系数(R2)和均方根误差(RMSE)反映拟合效果，结果显示2期TLS数据提取精度较高，拟合线性方程的决定系数(R2)均>0.9，表明提取的树高胸径与实测值之间存在线性关系,RMSE值较小，拟合效果良好，可将提取值用作实测值，进行后续分析(表3，图4)。

表3 树高和胸径线性拟合效果Table 3 Linear fitting effect of tree height and DBH

注：a：2017年树高；b：2017年胸径；c:2019年树高；d:2019年胸径。图4 树高胸径提取值与实测值回归分析Fig.4 Regression analysis of tree height and DBH extraction value and measured value

2.2 形数形率模型研建

随机选取样本木中80%的数据作为建模样本，剩余20%作为检验样本，根据材积三要素计算材积，并采用江苏省通用的杨树二元材积模型进行检验，结果如下:

(3)

采用决定系数，平均绝对误差(MAE)，均方根误差，相对均方根误差以及乖离率(Bias)对模型结果进行检验和评价。其中，R2越大越好；RMSE、MAE、rRMSE、Bias越小越好(表4)。

表4 模型计算结果分析Table 4 Analysis table of model calculation results

由表4模型检验结果来看，形数形率模型和材积估计值的决定系数(R2)分别为0.773和0.996，平均绝对误差(MAE)为0.25%和0.34%，均方根误差(RMSE)为 0.37%和0.63%，表明基于TLS计算的材积具有较高的准确性和客观真实性；相对均方根误差(rRMSE)表达了标准差相对平均值的变动程度，形数形率模型略高于材积模型。乖离率(Bias)较小且接近于0，表明预测值与实测值并未出现较大偏离[22-23]。结果表明通过定义试验形率拟合形数形率模型效果较好，可以根据调查需求合理使用该模型，为研究杨树干形参数胸高形数的测算提供参考。

2.3 基于2期TLS数据干形分析

2.3.1 基于2期TLS数据的不同造林密度下TAP累积量分析 根据V.Luomaetal[15]，树干削度为胸径与6 m处直径之差，本次研究分别选用胸径与6、7.6、9.6 m，0.382h以及(h+3)/2处直径进行对比，观察用以描述样本木树干削度的最优选择。采用最小二乘法提取树干上部直径，计算树木削度，采用单位面积TAP累积量对结果进行分析评价。

由表5可知，以6 m的传统定义计算树木削度，6 m×6 m、4.5 m×8 m、5 m×5 m、3 m×8 m这4种密度的单位面积TAP累积量分别为50%、45.5%、72.8%、65.7%。表明种植密度与树木削度之间无明显联系，不符合树木生长规律。利用传统定义的6 m或7.6 m、9.6 m，0.382 h处直径计算树木削度，均无法发现不同种植密度下树木干形变化的规律。说明杨树9.6 m以下直径较胸径生长速度较慢，变化较小，难以准确描述树木干形变化，因此无法用9.6 m以下的任意一段干形变化代表样本木整体的树木干形变化。以(h+3)/2处直径作为上部直径代表，可以描述杨树1.3 m与(h+3)/2之间树木的削度变化，最大限度描述了杨树干形变化，可以有效代表样本木整体的干形变化。

表5 不同直径下单位面积TAP累积量分析Table 5 Analysis table ofaccumulation of TAP per unit area under different diameters

2.3.2 基于2期TLS数据不同造林密度下材积与高径比分析 不同造林密度下测树因子变化分析结果见表6：造林密度对树木材积生长和高径比变化有显著影响。在同一无性系中，随造林密度的增加，材积增长量逐渐降低，表现为低密度的杨树材积增长量大于高密度杨树。分析不同样地的树木高径比变化可以发现，HDR的相对变化率均为负，且随造林密度的增加而增大。高密度株行距配置的3 m×8 m的杨树在3 m方向上的生长竞争加大，树木生长空间不足，严重影响了树木材积生长和高径比变化，表现为△VOL最小，为5.89%；而△HDR最大，为-3.44%。

表6 不同密度样地材积与高径比变化Table 6 VOL and DBH in different density plots

2.3.3 基于2期TLS数据干形变化分析 通过计算T1和T2时期的测树因子可知：TAP的变化有正有负，负增长居多，样本树的平均TAP变化量为负，T2时期的削度要

表7 干形指标变化Table 7 The changes in stem form indicators

3 结论与讨论

分析TLS数据提取林分树高胸径的精度情况，2期数据提取树高、胸径拟合线性方程的R2均>0.9，RMSE较小。显然TLS数据提取值可以代替实测值，用以后续分析。此外，基于TLS数据通过定义试验形率，探索了基于TLS数据提取林木参数的杨树形数形率模型，其中R2为 0.773，MAE为0.25%，RMSE 0.37%，rRMSE为0.477，Bias为0.000 14。显然试验形率可以较好地代替胸高形率，应用于林业调查中胸高形数的预测，为分析杨树干形参数提供参考。

在对比各高度直径后发现选用胸径与(h+3)/2处直径作为最具代表的下部直径与上部直径，计算TAP和TAP累积量最为合理，可以准确代表样本树的树木干形变化。分密度对林木干形进行分析显示，TAP累积量和材积变化量随造林密度增加而降低，高径比变化量随造林密度增加而增加，树木干形趋于尖削。

通过对研究区样木干形进行分析，发现样本树的平均TAP变化量为负，胸径与(h+3)/2处直径大小趋于一致，树木生长趋向通直圆满；样本树的材积有所增长，但增长率变化不大；q实和f1.3在研究期内变化均不明显，变化的范围也较为稳定。此外，HDR平均变化呈现负增长，可能与样本树的生长周期有关，树木生长中后期高生长趋于平缓，以直径粗生长为主，引起HDR的下降。

相较于传统干曲线严谨复杂的计算方法，TAP更加方便易算，用以代替连续的曲线方程，近似描述树木干形变化。但需注意的是，对于TAP和TAP累积量，其精确提取是树木干形变化分析的关键。研究中存在少数例外的观察结果。这可能是树木位置或周围树木对其生长产生了影响，例如外力导致树梢折断或生长阶段不同对树高、胸径增量的影响等。基于本试验结果，需要进一步探究不同无性系、不同发育阶段对树木材积增长和干形变化的影响，以便详细了解这一问题。