朴津欣 于颖 范文义
(东北林业大学,哈尔滨,150040)
叶面积指数(LAI)是陆表关键的特征参量之一,已被广泛应用于气候、生态、水文、生物地球化学等模型研究中。测量叶面积指数有直接法与间接法。直接法主要包括破坏性采样法、落叶收集法和异速生长方程法等。落叶收集法通过比叶面积估计叶面积指数,不具有破坏性[1-2]。破坏性采样法相对简单可靠,但过程具有破坏性且费时费力。与直接法相比,间接法更加快捷,但是基于光学仪器的叶面积指数可靠性仍然需要使用直接法进行验证,常用的光学仪器主要有LAI-2200、TRAC和半球摄影。LAI-2200的理论基础是贝尔定律,贝尔定律是在连续植被水平及垂直分层的假设下,得到的间隙率、消光系数、聚集度指数、叶倾角和叶面积指数之间的关系,应用贝尔定律进行估算叶面积指数的理论方法目前在全世界广泛应用[3-5]。
激光雷达(LiDAR)是一种以激光作为光源对目标物进行探测的主动遥感技术[6],在获取森林垂直结构参数方面有很大优势,近年来,机载激光雷达(ALS)和地基激光雷达(TLS)已经被广泛应用于森林的树高、冠层郁闭度、叶面积指数和地上生物量的估算中[7-8]。通过贝尔定律计算有效叶面积指数需要得到冠层的间隙率,激光雷达具有一定的穿透性,激光穿透指数可以用来定量描述冠层的间隙率,进而估算有效叶面积指数[9]。机载激光雷达一般可记录第一次回波、最后一次回波或者多次回波信息,然后利用这些回波信息计算激光穿透指数。目前已有的研究表明,利用间隙率反演叶面积指数没有出现饱和的情况[10-11]。无人机(UAV)平台与机载平台相比具有较高的空间分辨率和较低的成本,近年来利用无人机平台搭载激光发射装置以及摄像机获取激光雷达数据并估算森林生物量的研究越来越多[12-15]。然而由于受到森林冠层遮挡以及激光点云密度的影响,冠层内部及林下的一些信息无法被精确获取[16]。
地基激光雷达能够非破坏性地获取高分辨率的冠层详细三维信息,实现各种冠层参数的获取[17]。相较于以往使用光学遥感数据估算叶面积指数,地基激光雷达具有不受天气光照等条件影响的优点。目前已有很多研究基于地基激光雷达数据估算叶面积指数[18-24]。Hosoi et al.[18]利用TLS数据建立了单木的三维体素模型,通过分层体素切片法计算了叶面积密度与叶面积指数。但该方法没有具体区分出叶片与木质部分从而造成了叶面积指数的高估。Li et al.[19]在此基础上进一步评估了去噪、木质部分离和体素大小对叶面积指数估算的影响,并提出了一种新的相较于原有的方法计算叶面积指数精度更高的体素测定方法。由于地基激光雷达是从地面对整个林地进行扫描,在林下环境比较复杂的情况下,扫描可能受到较高灌木的遮挡、郁闭度较高的冠层遮挡以及激光点云密度偏小等因素的影响,导致林分上层信息无法准确获取,从而造成叶面积指数的低估。
虽然激光雷达可以收到多次回波,但在密度较高的林分下受树叶与枝干的遮挡,地基激光雷达对树冠上层点云数据的获取可能出现缺失的情况、无人机载激光雷达对林地下部点云数据的获取可能出现缺失的情况[15],所以将同一时期的地基激光雷达点云数据和无人机载激光雷达点云数据通过坐标配准、重采样等步骤进行融合,进一步研究融合后数据对估算叶面积指数的影响。
本研究应用地基融合无人机载激光雷达点云数据对叶面积指数进行估算,并且在此基础上对融合数据进行抽稀,分析抽稀后数据对估算叶面积指数精度的影响;同时使用融合数据提取树高、胸径、冠幅、冠面积、冠体积5个参数,应用这些参数通过多元线性回归法构建叶面积指数模型反演真实叶面积指数,并使用落叶收集实测数据对估算的真实叶面积指数进行验证,使用LAI-2200实测数据对估算的有效叶面积指数结果进行验证。
本研究以黑龙江省尚志市东北林业大学帽儿山林场为研究区域,隶属于尚志国有林场管理局,地理坐标为(45°2′20″~45°18′16″N,127°18′0″~127°41′6″E),平均海拔400 m,研究区域属大陆性季风气候,年平均气温2.0 ℃,年均降水量676 mm,主要乔木树种有红松(Pinuskoraiensis)、云杉(Piceaspp.)、樟子松(Pinussylvestrisvar.mongholica)、落叶松(Larixgmelinii)、水曲柳(Fraxinusmandschurica)、胡桃楸(Juglansmandshurica)、柞树(Querusmongolica)、椴树(Tiliaamurenisis)、色树(Acermono)、榆树(Ulmuspumila)、白桦(Betulaplatyphylla)、杨树(Populusdavidiana)等10多种。本研究选取为一块100 m×100 m的针叶纯林样地和一块100 m×100 m的针阔混交林样地。针叶纯林样地包含树种有樟子松、红松、落叶松。针阔混交林包含阔叶树种有白桦、水曲柳、柞树、榆树;针叶树种有红松、落叶松。在样地内每10 m×10 m节点处设置凋落物收集装置,如图1所示。样地基本信息如表1所示。
图1 帽儿山样地位置、LAI-2200测量点布设
表1 样地实测基本信息
由Trimble公司的Trimble TX6对研究区域进行多站扫描将样地完全覆盖得到样地的地基激光雷达数据。扫面原理是基于水平旋转基础上的竖直转镜,扫描速度为每秒储存500个数据点,架设仪器高为1.5 m,扫面范围是0.6~80.0 m,视场角为360°×317°,并集成相机,分别在两块样地内架设10个测站进行扫描,确保样地覆盖完全。扫面完成后进行内业处理,应用Trimble RealWorks软件对采集到的原始数据采用单点配准的方法进行拼接配准。两块样地地基激光雷达点云图如图2和图3。
图2 针阔混交林地基激光雷达点云图
图3 针叶纯林地基激光雷达点云图
由D200无人机飞行平台搭载RIEGL mini VUX-IUAV对研究区域进行扫描得到无人机载激光雷达数据,其相对高度为80 m,平均飞行速度为5 m/s,平均点云密度为每平方米10 000个点,激光脉冲波长为905 nm,扫描当天天气状况良好。数据存储格式为标准LiDAR储存格式las 1.2,包含了坐标信息、回波强度信息、回波次数信息、高程信息等。扫描完成后进行内业处理将航带进行拼接,确保扫描范围将样地完全覆盖。两块样地无人机载激光雷达点云图如图4和图5。
图4 针阔混交林无人机载激光雷达点云图
图5 针叶纯林无人机载激光雷达点云图
首先将地基激光雷达数据和无人机载激光雷达数据转换为统一的坐标系,通过外业测得的样地4个角点以及样地内4个控制点GPS坐标进行点云数据粗拼接,但由于冠层遮挡林下GPS信号相对较弱,单纯根据GPS坐标点进行拼接误差较大。进一步应用邻近点迭代算法[25](ICP)对粗拼接的点云数据进行校正,完成两者的数据融合。邻近点迭代算法是根据最小二乘法对临近点进行最优匹配的方法。数据融合完成后,通过计算地基激光雷达数据与无人机载激光雷达数据之间的平均距离进行融合结果检验,检验结果显示融合数据的平均距离小于无人机载激光雷达点云平均间距,融合数据点云分层小于6 cm,融合数据可以使用[26]。融合后激光雷达点云图如图6。
图6 针阔混交林融合激光雷达点云图
本研究对两块研究区进行每半月一次的外业落叶收集,根据收集的落叶进行样地真实叶面积指数的计算。计算方法如下:首先将收集的叶子进行内业处理,按树种分类。随机抽样计算不同树种叶面积(Sa),并计算比叶面积(ASL)[27]。
ASL=Sa/W。
(1)
式中:Sa是测量不同树种的平均鲜叶面积;W是Sa对应树种烘干之后叶子的干质量。
每次收集完毕分别称出每个收集框中不同树种的总湿质量进行记录,记为M,按树种将整块样地同一树种的叶子混合均匀后,随机抽取一部分叶子测其湿质量,记为M1,然后将抽取出的叶子进行24 h烘干处理,测其干质量进行记录,记为M0,计算不同树种的干湿比,记为a。根据比叶面积、干湿比计算该收集时间段内样地各树种凋落物对应的凋落物叶面积指数值,记为ILA。计算公式为(2~3)。
a=M0/M1,
(2)
ILA=aMASL。
(3)
对于落叶树种理论上,全年凋落物之和即为全年最茂盛时期的ILA值;对于常绿针叶树种,通过周转率计算出全年最大ILA,不同时期对应的叶面积指数值即为最大值减掉凋落物ILA值。计算得到两块样地真实ILA数据如表2所示。
表2 落叶收集法计算叶面积指数结果
本研究使用LAI-2200冠层分析仪对样地进行有效叶面积指数(LE)的测量。在样地垂直于太阳光照的方向上每隔10 m设置一条测线,将LAI-2200置于离地面高1.5 m的位置进行测量,测量使用90°遮盖帽(得到A值),A值的测量设置在样地外的空地,将另一只探杆带入样地,每条样线上隔10 m设置一个测点进行采样(得到B值),整块样地采样完毕,将A值和B值进行匹配计算有效叶面积指数。测量在阴天或清晨(03:00—05:30)或黄昏(16:30—18:00)进行,LAI-2200测量点布设如图1所示。
数据预处理主要包括点云数据配准、去除噪点、地面点分类以及点云数据归一化。点云去噪算法是根据目标点云邻域10个以内的点,计算该点与相邻点的平均距离值,计算所有目标点平均距离值的中值以及标准差,并根据中值、标准差计算最大范围值,计算公式如(4),如果该点的平均距离值超出计算得到的最大范围值则被认为是噪点,进行去除。
最大范围值=中值+标准差倍数×标准差。
(4)
点云数据归一化主要是为了消除地形对点云数据产生的影响,处理方法是将每一处点云的高程与数字高程模型做差,数字高程模型由去噪后的激光雷达数据提取出地面点通过线性插值生成。
地基、无人机载以及融合激光雷达数据均采用体素法[18]对有效叶面积指数进行估算。首先分别由X、Y、Z的最小值和最大值来确定点云数据的域,然后根据以下公式进行体素化。
(5)
(6)
(7)
式中:i、j、k是体素坐标组;Int是取整函数X、Y、Z经过预处理后的激光雷达点云数据的点坐标;Xmin、Ymin、Zmin是X、Y、Z的最小值;Δi、Δj、Δk是体素单元大小;i×j×k是体素单元的数量。
叶面积密度(DLA)是单位群落体积的总植物叶面积[19]。DLA在垂直方向上的累加即为叶面积指数。判断体素单元里是否有激光点,有激光点的体素单元被赋值为1,没有激光点的单元被赋值为0,其中没有点的单元被认为是冠层中存在间隙。高度h和高度h+ΔH内的DLA计算式如下。
(8)
(9)
其中,N(k)为第k层激光点频率,ΔH是体素水平层的厚度,mh和mh+ΔH是垂直轴上的体素坐标,与垂直坐标上的高度h和高度h+ΔH等效(h=Δk×ΔH),nI(k)是第k层标志为1的体素单元数量,nT(k)是第k层的总体素单元数量。估算有效叶面积指数的公式如下。
(10)
(11)
式中:G(θ)是与消光系数相关的投影函数,α(θ)为叶片倾斜度改正因子[27]。
对融合数据进行重采样处理,即点云抽稀,得到不同密度的点云数据。为了充分研究不同抽稀程度的融合数据对叶面积指数估算结果的影响,基于抽稀后点云密度方面的考虑,以0.5以下的倍数对融合数据进行点云抽稀后点云密度小于地基点云密度,故本文选择以0.9、0.8、0.7、0.6、0.5的倍数对整体融合数据进行抽稀。根据不同采样率抽稀后的融合数据通过体素法估算有效叶面积指数。
对地基-无人机载激光雷达点云融合数据进行单木分割[28]。根据分割结果提取树高、胸径、冠幅、冠面积、冠体积,将这些参数通过多元线性回归法建立方程,反演真实叶面积指数。
y=a0+a1x1+a2x2+…+anxn+ε。
(12)
式中:y为真实叶面积指数;x1、…、xn分别为提取的树高、胸径、冠幅、冠面积、冠体积信息;a为常数;ε为误差。采用决定系数(R2)和均方根误差(ERMS)评价回归模型的精度。
(13)
(14)
根据无人机载点云数据通过贝尔定律对有效叶面积指数进行计算,贝尔定律公式如下[29-30]。
(15)
式中:LE是有效叶面积指数;ang是平均扫描角;P(θ)是间隙率;k是消光系数。
P(θ)=nground/n。
(16)
式中:nground是提取地面点的数量;n是总体点的数量。
地基激光雷达数据体素法估算有效叶面积指数结果。根据地基激光雷达体素法估算两块样地的有效叶面积指数结果与LAI-2200实测数据对比图如图7~8所示。由图可知,两块样地估算的有效叶面积指数结果与LAI-2200实测数据有效叶面积指数相关系数均大于0.7,显著相关。针叶纯林叶面积指数的估算结果精度为R2为0.79,均方根误差(ERMS)为0.69;针阔混交林叶面积指数的估算结果精度为R2为0.74,ERMS为0.63。说明应用地基激光雷达数据基于体素法可以较好的估算有效叶面积指数,不论是针叶纯林还是针阔混交林在条件允许的情况下,都可以应用地基激光雷达数据对有效叶面积指数进行估算,结果比较可靠。但就不同林分类型而言,使用地基激光雷达数据估算的针阔混交林有效叶面积指数结果精度没有针叶纯林有效叶面积指数结果精度高。
图7 针叶纯林TLS数据-LAI-2200有效LAI比较
图8 针阔混交林TLS数据-LAI-2200有效LAI比较
无人机载激光雷达数据体素法估算有效叶面积指数结果。根据无人机载激光雷达体素法估算两块样地的叶面积指数结果与LAI-2200实测数据对比图如图9~10所示。由图可知,两块样地估算的有效叶面积指数结果与LAI-2200实测叶面积指数相关系数均为0.76。虽然相关系数均大于0.7,但ERMS分别为0.83、0.97,均方根误差都偏大,总体精度较低,于是补充贝尔定律估算有效叶面积指数进行对比验证。
图9 针叶纯林UAV数据-LAI-2200有效LAI比较
图10 针阔混交林UAV数据-LAI-2200有效LAI比较
融合数据体素法估算有效叶面积指数结果。应用融合数据基于体素法估算针阔混交林有效叶面积指数与LAI-2200实测数据对比结果如图11所示,结果表明两者的相关系数为0.84,说明显著相关,ERMS为0.54。相比于该样地地基激光雷达数据估算有效叶面积指数的结果(R2为0.74,ERMS为0.63),无人机载激光雷达数据估算有效叶面积指数的结果(R2为0.75,ERMS为0.66),融合数据估算的结果精度更高。图12为融合数据估算有效叶面积指数结果与地基激光雷达数据估算有效叶面积指数结果的比较,R2为0.90,ERMS为0.46,表明两者相关性很显著,原因是地基激光雷达的点云密度相比于无人机载激光雷达的点云密度高的多,因此其对整体融合数据的影响较大。
图11 针阔混交林融合数据-LAI-2200有效LAI比较
图12 针阔混交林融合数据-TLS数据估算有效LAI比较
分别以0.9、0.8、0.7、0.6、0.5采样率抽稀的融合数据估算有效叶面积指数的结果如图13所示,结果精度如表3所示。由图13可分析得到未采样抽稀的融合数据估算叶面积指数的结果与LAI-2200实测的叶面积指数相比,未抽稀数据会造成高估;并且分析不同采样率融合数据估算的叶面积指数可知,点云密度越大,估算叶面积指数越大。由表3可以分析得到采样率为0.9和0.8时R2分别为0.89、0.93,ERMS分别为0.43、0.33,估算结果较接近落叶收集实测的叶面积指数值。当采样率从0.8再减小时,即点云密度更小时,估算精度逐渐下降,R2随之变小,ERMS变大。但以0.7采样率抽稀的融合数据估算叶面积指数的结果(R2为0.76,ERMS为0.63)依然比单独使用地基激光雷达数据估算叶面积指数的结果(R2为0.74,ERMS为0.63)和单独使用无人机载激光雷达数据估算叶面积指数的结果(R2为0.75,ERMS为0.66)好。采样率为0.6和0.5会造成明显低估。以5个不同采样率抽稀的融合数据估算叶面积指数的结果表明,当采样率为0.8时估算精度最高,而且就计算时间方面而言,采样之后的数据所需要的计算时间也更短。
图13 不同采样率数据估算有效叶面积指数
表3 不同采样率融合数据估算有效叶面积指数精度
单木分割提取树高、胸径、冠幅、冠面积如表4所示,由表可知单木分割提取的参数精度较高,可以用来进一步多元线性回归模型的建立。融合数据多元线性回归法估算真实叶面积指数结果精度为R2为0.86,ERMS为0.62;P<0.000 1,模型拟合结果显著。各个变量参数估计结果如表5。结果显示树高与叶面积指数估计相关性较小(t值较大),应用融合数据多元线性回归模型反演叶面积指数主要受胸径、冠幅、冠面积、冠体积影响。
表4 不同林分类型的样地实测与单木分割
表5 变量参数估计
两块样地无人机载激光雷达数据根据贝尔定律估算有效叶面积指数结果与LAI-2200实测有效叶面积指数对比如图14~15所示。结果表明,两块样地无人机载激光雷达数据估算有效叶面积指数与LAI-2200实测有效叶面积指数的相关系数均大于0.7,显著相关。针叶纯林叶面积指数的估算精度为R2为0.78,ERMS为0.55;针阔混交林叶面积指数的估算精度为R2为0.75,ERMS为0.66。就估算精度而言,两种林分类型单独使用无人机载激光雷达数据估算叶面积指数的能力与单独使用地基激光雷达数据估算叶面积指数的能力相差不大;就林分类型而言,不论单独使用地基激光雷达点云数据还是单独使用无人机载激光雷达数据估算叶面积指数,对针叶纯林叶面积指数的估算精度均高于对针阔混交林叶面积指数的估算精度;但就计算量与计算时间而言,无人机载激光雷达数据激光点数较少,计算需要的时间较短,对计算机的配置要求较低。
图14 针叶纯林UAV数据-LAI-2200估算有效叶面积指数比较
图15 针阔混交林UAV数据-LAI-2200估算有效叶面积指数比较
单独使用地基激光雷达数据估算有效叶面积指数和单独使用无人机载激光雷达数据估算有效叶面积指数时,针阔混交林的估算结果都没有针叶纯林的估算结果准确。原因是针阔混交林中阔叶对激光雷达点拦截更多,导致估算精度低。表6为两块样地不同方法估算叶面积指数对应的结果精度表。
表6 3种数据对应的方法估算叶面积指数结果精度
叶面积指数与森林生物量密切相关,是反映植被冠层结构最基本的参数之一,控制着植被生物、物理参数的变化[31],因此准确估算叶面积指数十分重要。本文对结合使用地基激光雷达点云数据与无人机载激光雷达数据估算叶面积指数方法进行了研究,对比了以不同采样率对融合数据进行抽稀再计算叶面积指数的结果差异,重点解决了使用单一类型的激光雷达数据造成信息缺失的问题。结合上述结果,分析得到以下结论:
(1)融合数据采用体素法对有效叶面积指数进行估算,估算精度较高(R2为0.84,ERMS为0.54),应用融合激光雷达点云数据进行单木分割提取单木参数(树高、胸径、冠幅、冠面积、冠体积)构建多元线性回归模型反演真实叶面积指数的结果精度也很高(R2为0.86,ERMS为0.62),两种方法结果精度相差不大,但均高于单独使用地基激光雷达数据和无人机载激光雷达数据估算时的估算精度。无人机载激光雷达数据体素法估算有效叶面积指数精度最低(R2均为0.76,ERMS分别为0.83、0.97)。融合数据可以更加充分的展现样地三维信息,弥补了单独使用地基激光雷达数据和单独使用无人机载激光雷达数据时信息获取不完全的缺陷。这一结果验证了在林分冠层密度较高的情况下,融合数据可以更加准确估算叶面积指数的论证。缺点是融合数据点云数量十分庞大,对计算机配置要求比较高,计算需要的时间较长。
(2)以不同采样率抽稀的激光雷达数据估算叶面积指数的结果可知,不是点云密度越高估算叶面积指数结果越准确,本文中针阔混交林融合数据的采样率为0.8时估算结果精度最高(R2为0.93,ERMS为0.33)。点云抽稀后需要的计算时间也大大缩短。未抽稀的融合数据估算叶面积指数的结果精度低于采样率为0.8的融合数据估算叶面积指数的结果精度主要原因有以下两个:一是融合数据中不仅包含地基激光雷达的噪点也包含无人机载激光雷达的噪点,去噪时公式(4)计算出的最大范围值会增加,导致去噪时一部分噪点保留,影响估算结果;二是两者数据融合后点云密度很大,张建鹏等[32]利用单木地基激光雷达数据探讨了点云密度对单木叶面积指数反演结果的影响,结果表明,点云密度越大,估测叶面积指数越大。未抽稀的融合数据激光点冗余,会造成高估,导致估算结果精度降低。
(3)不同林分类型方面,各单一类型激光雷达数据对针阔混交林叶面积指数的估算精度均没有对针叶纯林叶面积指数的估算精度高。原因是阔叶的叶面积与针叶相比较大,遮挡性更强,对激光的拦截更多,导致获取的林分三维信息更加不完全,于是在单独使用一种类型激光雷达数据对叶面积指数进行估算时,针叶纯林的估算结果精度较针阔混交林的估算结果精度更高。
张颖等[33]利用地基激光雷达对枝下高的高度进行提取,当相对高度在0.6~0.8时,枝条平均提取精度为0.775,原因是地基激光雷达数据很难获取树梢部分的点云数据,林分整体信息获取不完整。本文基于融合数据,应用体素法和多元线性回归法估算叶面积指数,估算精度均达到0.8以上,融合数据以0.8的采样率进行点云抽稀估算叶面积指数的精度最高(R2为0.93,ERMS为0.33),可以准确估算叶面积指数,说明融合数据可以很好地解决单一类型的激光雷达数据造成信息缺失的问题。但本文抽稀采用的是整体点云以相同采样率进行抽稀的方法,如果将林分分层按不同的采样率进行抽稀,估算结果应该更高。