周宇飞 王振师 钟映霞 吴泽鹏魏书精 李 强 李小川
(广东省森林培育与保护利用重点实验室/广东省林业科学研究院, 广东 广州 510520)
桉树Eucalyptus spp.作为世界著名的速生树种,凭借其适应性强、用途广泛、速生高产等特点,是华南地区人工速生林的优选树种[1-2]。桉树人工林蓄积调查是林业生产经营的重要工作[2]。通常的蓄积调查方法包括机械角规样地调查法、标准样地调查法、全林每木检尺法以及样圆调查法等[3],这些传统的蓄积量调查方法均需要耗费大量人工,且数据质量难以保障。无人机可以获取较大面积林地的林木几何参数,在林业应用中具有成本低、效率高等优势[4-5]。无人机机载激光雷达测树技术是近年迅速发展起来的遙感应用技术[6],机载激光雷达可以自动获取较大面积林地的激光点云数据,经过专业软件处理后可得到株数、单木树高及冠幅等结构信息的客观数据,为胸径和蓄积的计算提供基础,从而作为桉树人工林蓄积量估测的依据[6-9]。
国外随着无人机机载激光雷达应用的成熟,已经在估算蓄积量[10]、生物量[11]、树木高度[12-13]等多方面得到了应用。利用无人机机载激光雷达估测林木的蓄积量已经成为林业经营管理应用的趋势。
试验区位于广东省肇庆市大南山林场,位于北 纬23°21′15″ 至23°21′32″, 东 经112°45′26″ 至112°45′44″之间,总面积约17.5 hm2。该试验区属亚热带季风气候,年平均气温 21.3℃。极端最高气温38.5℃,极端最低气温-1.2℃,最低气温多出现在1月或2 月份,年平均降雨量约1 803.6 mm,年均日照1 702.3 h,年均雷暴日89 d。林地海拔50~200 m,坡度10~30°。土壤为砖红壤,土层厚0.5~1.2 m,适合桉树生长。详细测量的试验样地位于有道路等标志分割的一处独立林班,面积0.81 hm2,是尾巨桉Eucalyptus urophylla×E. grandis 人工林,林龄约6 a的成熟林,种植株行距2.5 m ×3 m,林木长势良好,无明显风折风倒木、病木和死木。
图1 试验区和样地位置卫星影像Fig.1 Satellite imagery of test area and plot location
无人机上搭载RIEGL VUX-1 机载激光雷达系统,于2019 年3 月18 日对调查区域进行机载激光雷达数据采集,最终获得las 数据格式的点云数据。无人机搭载激光系统的飞行技术性能见表1。
无人机载激光雷达测量完成后,于2019 年7—8 月对试验样地内的林分进行采伐,详细记录试验样地内所有桉树单木的胸径、树高、株数。采伐前使用胸径尺测量单株胸径,并在伐倒后使用皮尺量取树高。
单木蓄积量采用广东省桉树类二元立木材积表的计算公式[8]:
式(1)中V 为单木蓄积量、D 为单木胸径、H 为单木树高。
使用LiDAR360 软件处理激光雷达点云数据。首先将采集获取的激光雷达点云原始数据进行裁剪、去噪、地面点自动分类形成数字高程模型(Digital elevation model,DEM)等操作;然后利用生成的DEM 对点云数据进行归一化处理,消弱地形起伏的影响;最后使用LiDAR360 软件的单木分割功能自动进行单木分割。
表1 RIEGL VUX-1 机载激光雷达系统技术性能Table 1 Technical performance of airborne LiDAR system RIEGL VUX-1
表2 实测树高和胸径的统计分析Table 2 Statistical analysis of measured tree height and DBH
表3 实测胸径和各模型反演胸径的欧氏距离Table 3 The euclidean distance between the measured DBH and inverted DBH of various models
单木分割功能可以自动获得测量区域内的林木株数,以及每株林木的树高、冠幅等[14-16]。由于激光雷达软件处理无法获取单木的胸径数据,故可利用样地调查所获取的单木数据建立胸径对树高的回归模型,并利用该模型为激光雷达测量的单木胸径进行赋值,单木的蓄积也可以利用式(1)的方法计算,由此估测整个试验样地的蓄积量。
经调查,该试验样地内桉树共有491 株,单木胸径和树高的测量值如表2 所示,所有单木按胸径和树高的分布如图2 所示。
图2 实测桉树单木树高和胸径分布Fig.2 Measured data distribution of eucalyptus’ singletree height and DBH
由于所有单木的胸径和树高之间的相关系数r值为0.943 2,因此,这些树木的胸径与树高之间呈强相关关系。利用树高,采用常用的线性函数、二次函数、幂函数和指数函数分别对单木的胸径进行回归,以选择最合适的胸径对树高的回归模型,这些模型的表示如下所示:
其中,式(2)为线性函数,式(3)为二次函数,式(4)为幂函数,式(5)为指数函数。这些式中,X 是树高,Y 为胸径。通过实测树高,利用各个模型反演胸径,再与实测胸径进行对比,以判断最为适合的回归模型。对于相同数目的两组数组,可以采用欧氏距离来判断反演胸径和实测胸径之间的相似性。欧氏距离越小,则两组数组之间的相似性越高。各模型反演获得的胸径和实测胸径欧氏距离见表3。
由此可见,采用指数函数建立的回归模型,即式(5)来表示样地内桉树林树高与胸径的关系最为合适。
利用式(1)计算桉树单株的材积并求和,可得到该试验样地内桉树总蓄积量为90.485 m3。
激光雷达点云数据处理后的单木分割如图3所示。
图3 激光雷达点云数据的单木分割Fig.3 Tree partitioning diagram of LiDAR point-cloud data
通过Lidar360 软件处理后得到试验样地内桉树株数共有358 株,激光雷达自动获取的单木树高及通过式(5)回归获得的单木胸径如表4 所示。
利用式(1)计算桉树单株的材积并求和,可得到激光雷达测量的该试验样地的桉树总蓄积量为73.757 m3。
激光雷达测量值与实际测量值的差异见表5。
总体看来,激光雷达测量的株数比实际株数偏少27.09%,平均树高偏高8.67%,平均胸径偏大8.25%,蓄积量偏小22.68%。究其原因,是由于LiDAR360 采用的是基于种子点的单木分割算法,在无人为干预的情况下,点云数据经单木分割后所获取的林木数量,是基于种子点的冠层顶部的峰值数量,而低于林分主冠层的林木,受上层树冠的覆盖和包涵,通常会被忽略,因此激光雷达获取的总株数偏少。
由于激光雷达数据忽略低于主冠层的林木属于普遍现象,导致最后的总蓄积量比实测偏少。在大面积桉树成熟林的蓄积计算中,可以用22.68%的偏差对激光雷达测量的总蓄积量进行修正。其总蓄积量可表达为:
其中V 为总蓄积量,V0为激光雷达测量的总蓄积量。
3.1 随着机载激光雷达技术在森林资源调查中的应用,将改变传统的森林调查模式,是提高调查数据准确性和调查效率的重要技术手段。机载激光雷达经过处理,可自动获取较大范围桉树人工林的株数、树高、冠幅等信息,利用桉树树高与胸径的强相关性可建立二者之间的回归模型,即利用式(5)间接获取胸径,获得了株数、树高和胸径,就可估算测量区域的桉树蓄积量。
表4 激光雷达获取的树高和胸径的统计分析Table 4 Statistical analysis of tree height and DBH obtained by LiDAR
表5 激光雷达测量值与实际测量值的差异比较Table 5 Difference comparison between the value obtained via LiDAR and field measured value
在无人为干预处理激光雷达点云单木分割的情况下,比较实地测量与激光雷达测量数据,激光雷达测量的株数比实际株数偏少27.09%,平均树高偏高8.67%,平均胸径偏大8.25%,总蓄积量偏小22.68%。这主要是由于在桉树成熟林中,受上层树冠的阻挡,激光雷达点云的单木分割算法难以将低于主冠层的林木纳入计算。由于不计算低于主冠层林木在激光雷达单木分割算法中属于普遍现象,桉树成熟林的总蓄积量可以采用式(6)的方法来对激光雷达测量的总蓄积量进行修正获得。
3.2 由于受郁闭度和地面灌草植被覆盖的影响,激光雷达点云难以直接打到地面,故激光雷达获取的树高较实际树高要低。本次试验区实测最高树木为27.0 m,而激光雷达测得的最高树高为25.84 m,其高差相差约1.1 m,这与其它研究发现的LiDAR 技术测量树高比实际树高低0.3~3.0 m相符[17]。由于本次试验难以将激光雷达获取的单株林木与试验区的单株林木一一对应,需进一步研究激光雷达树高测量误差的原因。
此外,由于机载激光雷达测量无法准确获取单株林木的胸径,胸径是通过对试验区实际测量树高与胸径建立的回归模型,再与激光雷达测量的树高回归计算得出,仅能代表该试验区的树高与胸径关系,今后可进一步扩大测试区域,调整该回归模型使之更具普适性。