基于机载大光斑激光雷达的森林冠层高度估测

孙忠秋,吴发云,高显连,高金萍,胡 杨

(1.国家林业和草原局调查规划设计院,北京 100714；2.北京林业大学,北京 100091)

森林冠层高度是反映森林垂直结构的重要参数,对定量化估测森林蓄积量、生物量及动态变化的研究具有重要的意义[1-3]。传统森林调查中,林分平均高一直通过打样地的方法进行人工测量,不但费时费力,而且难以对大面积森林区域的树高进行采集[4]。被动光学遥感影像的获取,在一定程度上提高了森林水平信息估测的能力,却限制了对森林的垂直结构信息进行提取,从而无法满足森林资源探测的需要[5-7]。微波遥感具有一定的树高估算能力,但其信号易在高郁闭林区达到饱和,从而在冠层高度估测方面有欠缺[8]。激光雷达作为一种主动遥感技术,其在获取森林垂直结构信息方面有其它技术无可比拟的优势[9-10]。根据传感器发射光斑直径的大小,激光雷达又可分为大光斑激光雷达(8～70m)和小光斑激光雷达(<1m)2种数据样式[11]。近30年来,国外对大、小光斑激光雷达在林业探测中的应用进行了大量的研究,成功地获得了森林平均高、单木平均高、郁闭度、冠层垂直结构、生物量、蓄积量等一系列森林监测参数。在国内,中国林业科学院组织设计的机载小光斑激光雷达系统(LiCHy)获得了大量的森林样地数据,精确地探测了研究区的林木高度、郁闭度等森林参数[12]。然而,受飞行成本和数据量存储等因素的限制,小光斑激光雷达很难对大区域的森林资源进行探测[13]。

星载激光雷达作为全球唯一一个由卫星搭载的大光斑激光雷达系统,其可通过发射到达地面直径为70m的光斑对陆地地物进行探测。在森林探测领域,已有很多学者在全球很多区域对其接收的回波波形进行分析,并成功地应用于森林冠层高度估测[2,14-20],表现出优良的森林探测性能。目前,国家林业和草原局卫星林业应用中心组织研发设计的机载林业探测大光斑激光雷达系统在湖南张家界地区进行了大面积的飞行试验,该系统通过连续发射直径30m的光斑对地面信息进行了数据采集。本文首次对机载林业探测大光斑激光雷达系统收集的波形数据进行了噪声估计、信号始末判断、森林冠层高度计算,利用样地数据对波形数据估测的森林冠层高度进行精度验证、探讨了机载林业探测大光斑波形数据估测森林冠层高度的潜力。

1 研究区概况

研究区位于湖南省西北部桃源县境内,地理位置为28° 55′N,111° 29′E,东部地势平坦、海拔50m上下,西部地势起伏较大,海拔200～1 000m,南部与北部地势均由西向东倾斜,海拔一般在50～200m之间。该区属中亚热带季风湿润气候,年平均降水量为1 447.9 mm,年平均相对湿度为82%,四季冷热分明,干湿两季明显。该区树种类型丰富,分布广泛。树种主要有栎类、马尾松、杉木和樟类。截至2014年,桃源县森林面积达到 227 740.1hm2,蓄积量达到 11 598 043hm2,其中乔木林面积为 185 993.2hm2、蓄积量为11 164 995m3,分别占桃源县森林总量的81.67%和96.27%[21]。

2 研究方法

2.1 数据来源

1) 机载林业探测大光斑激光雷达数据。机载林业探测大光斑激光雷达系统是目前我国唯一的机载大光斑激光雷达设备。该设备于2017年12月下旬,在张家界荷花机场附近区域开展了首次森林资源探测试验,该次飞行试验共包含2个架次,航线依据树种类型主要覆盖了桃源县、花垣县、保靖县、吉首市、永顺县、沅陵县和古丈县,航线总长超过980km,航高超过3 000m,航时总计11h,获取了研究区大量有效的大光斑数据(图1)。为了便于分析,本文只选取了桃源县境内光斑号为2002～2011的连续10个大光斑波形数据。

2) 野外调查数据。地面调查样地于2017年底进行,预先划分各主要树种的树高级和郁闭度级,对湖南省7个主要树种在大光斑飞行区域投影下的样地进行了地面调查,共调查了588个投影半径为15m的地面样地。调查过程中,利用全站仪与差分GPS相结合的方式精确定位调查样地样圆中心和样木位置。利用超声波测高仪测量样地每株样木的树高和枝下高。利用激光测距仪、胸径尺测量每株样木的冠幅和胸径,并记录样地的基本信息[22-23]。

图1 研究区大光斑波形位置

2.2 数据处理和分析流程

利用机载林业探测大光斑激光记录的回波数据、实测样地数据对光斑尺度森林冠层高度进行估测、验证,具体流程如图2所示。

图2 机载林业探测大光斑数据处理流程

2.2.1波形数据处理

机载林业探测大光斑激光雷达设备完整记录回波信号,信号采样间隔1ns。为了减少数据量,提高处理与计算的速度,在找到回波能量最大值后,在最大值前后分别取200 ns和300 ns数据作为回波信号产品保存下来。估测森林冠层高度,首先要对回波信号进行处理,其包括噪声估计,阈值计算,波形起始点判断,地面点估计等几个步骤。

1) 噪声估计。在大光斑设备接收到地面或者树冠回波之前,最开始接收到的是信号的背景噪声。本文经过对200个回波波形分析其前后噪声,发现其前后背景噪声值存在一定的差异。因此,本文对大光斑波形信号开始前的噪声和结束后的噪声分别进行计算。

采用波形信号前50帧数据(即有效波形范围内的第1～50个记录,有效波形共500帧数据)进行均值、方差的估计。即:

Noisebegin=f(w1,w2,w3,…,w50)

MEAN_Noisebegin=mean(Noisebegin)

VAR_Noisebegin=var(Noisebegin)

采用波形信号后50帧数据(即有效波形范围内的第451～500个记录,有效波形共500帧数据)进行均值、方差的估计。即:

Noiseend=f(w451,w452,w453,…,w500)

MEAN_Noiseend=mean(Noiseend)

VAR_Noiseend=var(Noiseend)

式中：wi为第i个回波信号的电压值,f为进行方差估计的函数。MEAN_Noisebegin,VAR_Noisebegin,MEAN_Noiseend,VAR_Noiseend分别为起始噪声信号和结束噪声信号的均值及标准差。

2) 信号阈值计算。考虑到机载林业探测大光斑激光雷达接收的信号尚属于首次处理分析,信号的波动等情况与GLAS接收的信号有一定的相似性,因此本文采用庞勇等[23]介绍的方法,计算信号的阈值。采用双阈值的策略。

thresholdbeg=MEAN_Noisebegin+ 2VAR_Noisebegin

thresholdend=MEAN_Noiseend+ 2VAR_Noiseend

式中：thresholdbeg用于提取信号起始位置,thresholdend用于提取地面位置。

3) 信号起始位置判断。波形信号的起始位置代表研究区内最大树高的位置。定义信号开始位置为波形信号连续3帧均大于该信号起始阈值的位置。即:

Sigbegin=(Wi>thresholdbeg)&(Wi+1>thresholdbeg)&(Wi+2>thresholdbeg)

4) 地面波峰判断。考虑到研究区森林覆盖率较高,地面回波的能量较小,所以激光发射到地面反射回波后信号趋于结束。因此利用上述方法进行地面波峰位置的判断。定义地面位置为波形信号连续3帧均大于该信号结束阈值的位置。

Sigground=(Wj-2>thresholdbeg)&(Wj-1>thresholdbeg)&(Wj>thresholdbeg)

5) 森林冠层位置判断。森林冠层位置即为激光能量发射到森林冠层表面返回时的信号位置。在森林密集区域,森林冠层处的回波信号能量最强,因此森林冠层的位置是回波波形第1个峰值位置。

Sigcanopy=Max(Wi)

6) 森林冠层高度估测。根据判断出的森林冠层信号位置和地面信号位置(图3),即可计算森林冠层回波位置与地面回波信号位置之间的高度H。即:

图3 森林冠层高度估测

H=(Number(ground)-Number(canopy))×

Resolution

式中:H为森林冠层高度,Number(ground)为信号结束位置Sigend光斑记录的位置,Number(canopy)为森林冠层位置Sigbegin光斑记录的位置。Resolution为激光回波记录的分辨率。在本文中,激光信号按照1ns的时间间隔记录,则Resolution为15cm。

2.2.2样地数据处理

地面调查样地采用WGS84坐标系,UTM 49N投影,选择距离大光斑波形投影最近的样地,对样地内调查的样木数据进行如下计算,为了选出最优的大光斑估测值,共计算7种高值。

1) 平均树高,样地内所有单木高度的平均值

式中:Hmean为平均树高；hi为单木树高；N为样木总数。

2) 胸高断面积加权平均高

式中:Hbasal为胸高断面积加权平均高；D2为胸高断面积；hi为单木树高；N为样木总数。

3) 树冠面积加权平均高

4) 优势木平均高D

在优势树种(组)中选择3株胸径最大的林木测定其树高,取算术平均值作为小班的优势木平均高。

式中:Hdominant-D为优势木平均高D；hDmax为优势树种胸径最大树高；hDmax-1为优势树种胸径次大树高；hDmax-2为优势树种胸径第三大树高。

5) 优势木平均高H

在优势树种(组)中选择3株树高最高的林木测定其树高,取算术平均值作为小班的优势木平均高。

式中:Hdominant-H为优势木平均高H,hHmax为优势树种最大树高,hHmax-1为优势树种次大树高,hHmax-2为优势树种第三树高。

6) 优势树种算术平均高

式中：Hdominant-mean为优势树种算数平均高,hi-dominant为优势树种单木高,n为优势树种株数。

7) 优势树种平均木平均高

测每棵优势树种的胸径,计算平均胸径,找到胸径接近林分平均值的3棵,然后测树高,计算平均值。

式中:Hdominant-Dmean为优势树种平均木平均高；hDmean,hDmean-1,hDmean+1为胸径最接近平均胸径对应树高。

3 结果与分析

根据大光斑和样地的位置,对平原区域采集的大光斑波形数据进行森林冠层高度估测。大光斑的光斑直径为30m,每秒重频40次,光斑的地面覆盖重叠度达到86.97%,因此，取连续相邻的10个大光斑波形数据估测该组光斑对应样地的森林冠层高度,并与地面样地数据进行对比验证。光斑号2002～2006 原始波形如图4所示。

图4 全波形原始信号

对原始波形进行噪声估计、信号开始位置判断、结束位置判断,得到回波地面位置和冠层位置,从而估测出森林冠层高度,如表1所示。

由表1可知,10组大光斑回波数据对样地冠层高度均值为10.32m,RMSE=1.31m,可见大光斑回波波形估测森林冠层高度较稳定。

表1 大光斑波形分解结果统计

对样地内样木检尺数据进行处理,并与森林各参数进行对比,结果如表2所示。

表2 森林各参数估测精度

通过直方图5可知,机载林业探测大光斑回波波形对森林各参数的估测能力依次为胸高断面积加权平均高、优势树种优势木平均高、优势树种算术平均高、算术平均高、树冠面积加权平均高、优势木平均高D和优势木平均高H。其中,对森林参数胸高断面积加权平均高估测能力最好,绝对误差0.47 m,相对误差4.36%,RMSE为1.4 m。

图5 森林各参数估测精度直方图

4 结论与讨论

1) 本文首次对国家林业和草原局卫星林业应用中心设计研发的机载林业探测大光斑激光雷达系统采集的大光斑波形数据处理方法进行了探讨,对大光斑覆盖下的森林冠层高度进行了估测,并利用地面样地调查数据对大光斑波形估测的森林冠层高度进行了验证。结果表明,大光斑激光雷达对平原地区森林冠层高度估测能力较强,通过对10组大光斑波形数据分析,大光斑对7种森林冠层高度估测均有不同程度的估测能力,其中以胸高断面积加权平均高、优势树种平均木平均高估测效果最好,相对误差分别为4.36%和8.29%,RMSE为1.40 m和1.55 m。对优势木平均高H、优势木平均高D估测能力最差,相对误差为19.81%和22.00%,RMSE为2.99 m和3.34 m。

2) 由于是首次对机载林业探测大光斑波形数据进行分析,本文只选取了平原地区的大光斑波形数据对森林冠层高度进行估测,结果只能说明大光斑波形数据对平原地区的森林有较好的估测能力,但波形还受地形等因素的影响[25],所以,考虑不同地形条件估测森林冠层高度将是下一步研究分析的方向。