陈宗铸 ,杨 琦,雷金睿,陈小花,李苑菱
(海南省林业科学研究所,海南 海口 571100)
传统森林空间结构地面测定方法仅能得到一些点上数据,很难及时获取区域范围森林参数的空间分布信息,而快速发展的遥感技术能够满足此应用需求[1-3]。目前,遥感技术已被广泛应用于林业资源管理、动态监测与分析、灾害监测、预测、评估等领域,遥感图像光谱信息具有良好的综合性和现势性,能够反映森林空间分布的变化特征。但是,大多数遥感传感器只能提供水平方向的信息,垂直方向的信息(例如林木垂直方向的信息)很难从遥感图像中获得,而LiDAR技术的出现引起了自然资源管理和研究方面的重大变化。激光雷达(Light detection and ranging,LiDAR) 在森林垂直结构测量方面具有无可比拟的优势,将遥感图像、激光雷达数据与地面测量得到的样点数据相结合,可以更精准地获取大范围森林类型及三维空间结构的分布特征,通过多时相遥感数据序列分析森林结构的动态变化规律[4]。激光扫描仪发出的激光可以部分穿过植被的空隙到达地面,系统接收植被和地面反射激光,获得森林垂直结构参数。然而,LiDAR技术自身的相关特性限制了获取的垂直信息精度。有研究表明,在LiDAR数据中,中密度至高密度冠层林的林分高度通常被低估,由于激光脉冲的采样间距,激光脉冲刚好打到树冠顶点的概率相对较小[5-6],该因素是影响LiDAR数据估计参数精度的重要因素。
国内外对基于遥感的森林参数估计均进行了大量的研究,吴楠、李增元、吕国屏、胥喆等多数学者对国内外林业遥感应用的研究概况与进展进行充分的归纳与总结[7-10],目前雷达、激光雷达等主动遥感方式是反演林分结构较为有效的方法[11],国内许多学者对森林蓄积量的遥感估测开展了大量的研究工作[11-23],王雪军等研究了基于遥感大样地抽样调查森林面积监测,实践了遥感样地判读与地面样地抽样调查相结合的大样地调查方法的应用,发现影像分辨率对判读精度影响较大[24];韦雪花研究了轻小型航空遥感森林几何参数提取,结果表明,轻小型机载LiDAR数据和高分航空影像提取的样地尺度的林分平均高、平均冠幅以及单木尺度的树高和冠幅等参数,精度能达到72%~88%[25];刘清旺等结合了激光雷达和ALHIS对湖北典型亚热带森林开展了航空遥感试验,结果分析表明,激光雷达对林分平均高的估测精度达到90.67%,ALHIS能够获取教高分辨率的植被遥感特征数据[26];林怡、陈健、刘鲁霞等利用LiDAR数据进行单木识别与胸径、树高等数据提取[27-30],但提取数据量较大,其算法时间成本较高,因此仍有待提高。总体来说,基于LiDAR的参数估计效果要高于基于遥感影像的参数估计,因为林业管理者对森林垂直结构(如林冠高度)的估算是最感兴趣的[31],而LiDAR与遥感影像的结合可以提高参数的估计精度[32]。目前基于LiDAR数据对森林林冠高度的估测研究已成为一个新的研究热点,国内越来越多的学者相继开展了相关研究工作[33-37]。
国外多名学者对平均树高与冠层模型高度之间的预期差异开展了详细研究[38-41]。众多研究表明,使用几何概率来估计给定平均冠幅的激光命中的平均垂直位置,然后将该位置与研究区域林分的实际树高进行比较,并通过将计算的差值与激光估计树高相加来计算平均树高,提高了估计值与实测值的相关性[39-41]。但这些研究多集中于特定的区域或特定的树种,推广至其他区域或树种需要重新对算法进行验证[39,42]。
LiDAR技术在大面积森林作业研究中有着非常突出的成效,然而国内目前开展基于LiDAR技术在热带地区森林空间结构方面的研究却很少,主要是由于热带地区的森林冠层郁闭度较高,其林下植被复杂多样,这使得LiDAR在勘测过程中可能会无法完全穿透树枝树叶到达地面而提前反射从而造成估测结果产生偏差。本研究在总结前人大量已有研究成果的基础上,以海南省三亚市热带区域森林为研究对象,并基于LiDAR点云数据的动态搜索窗口算法,提出了一种改进的冠层估计方法,并提出了LiDAR数据估计的平均高度与林分平均高度之间的关系,总结出LiDAR数据估计的平均树高与实地调查平均树高之间的关系,通过这种关系,可以从LiDAR数据中快速估计其他区域热带森林的平均高度。
三亚市地处海南岛最南端,地理坐标为 18°09′34″~ 18°37′27″N、108°56′30″~109°48′28″E,东邻陵水,北依保亭,西毗乐东,南临南海,陆地总面积1 919.58 km2;其北靠高山,南临大海,地势自北向南逐渐倾斜,形成狭长状的多角形,境内海岸线长258.65 km,有大小港湾19个,大小岛屿40个;市区三面环山,北有抱坡岭,东有大会岭、虎豹岭和高岭,南有南边岭,形成环抱之势,山脉连绵起伏、层次分明;其地处低纬度,属热带海洋性季风气候区,年均气温25.7 ℃,全年日照时间2 534 h,年均降水量1 347.5 mm。境内有中、小河流多达12条,主要分为中、东、西部三个水系。中部水系以三亚河为主,包括大茅水;东部水系以藤桥河为主,包括藤桥西河和藤桥东河;西部水系以宁远河为主,为三亚第一,海南第四大河。三亚市地表水年总径流量11.5亿m3,地下水动储量1.42亿m3,并具有丰富的地热资源。
三亚市动、植物资源较为丰富,有关调查资料显示,境内有植物(含栽培品种)共1 838种,其中属国家重点保护的野生植物有坡垒、无翼坡垒、红榄李、油楠、青梅、白桫椤等;境内分布的陆生野生动物多达300余种,其中兽类50种,鸟类300余种,属国家重点保护的野生动物有云豹、孔雀雉、海南山鹧鸪、白鹇、穿山甲、猕猴、水鹿等。
实验区域选择在三亚市福万水库北侧山地树林,范围为50 hm2,共实地采集高程点48个,单株树木69株,同时采集了树木的胸径、树高,以及实地照片。
采集LiDAR数据的设备为德国TopoSys公司的Harrier 68i,飞机采用运-5。激光扫描仪的型号为Riegl LMS-Q680,最大脉冲频率400 kHz,扫描角度45°~60°,该仪器的GPS、IMU、激光器、相机高度集成、固化在设备箱内。该技术是以激光为测量媒介的新型遥感技术,LiDAR数据点密度为1.51点/m2。
较早提出激光LiDAR数据完整处理过程的是Hoffman和Jain[43],后来的数据处理研究也大多按照这种过程进行。总体而言,激光扫描系统获取的LiDAR数据处理一般按照如图1所示流程进行,包括数据采集、数据预处理、点云分割(或称滤波分类)、地物分类、对象建模5个步骤。基于这些步骤,对于森林应用,本研究的数据处理过程包括数据采集、数据预处理、DEM和DSM提取、CHM生成、平均树高等参数估计、精度验证及结果改进等步骤,如图1所示。
图1 实验区域实地核查数据(红色为高程点、绿色为单株树木)Fig.1 Field verification data in experimental area (Red is elevation point, green is single tree)
数据收集和预处理不在本研究讨论。本研究重点研究DEM和DSM提取,CHM生成和树高估计。
鉴于林地地形的复杂性,本研究采用改进的移动窗口差分算法。该算法首先对点云数据进行网格化,获得每个网格内的最高点、最低点,作为初始DSM、DEM参考面,然后根据初始参考面,通过分别设定不同的阈值,对剩余未分类的点云进行细分,以将DSM、DEM进行细化。最后将DEM、DSM点云集合存入文本文件输出。算法的具体思想和步骤如下:
2.1.1 拟合DSM与DEM初始参考面
首先在三维激光点云文件中自动找出、坐标的最小值和最大值,即(Xmin,Ymin)和(Xmax,Ymax)。选择格网间隔大小为2 m,然后分别计算出坐标的最大值和最小值的差值并取整。以左上角坐标为起点,建立格网间隔为2 m的矩阵,然后将激光点云的点转化到格网中,将每个格网中高程最小的(格网中无高程的不计入)点标记为地面点,作为DEM初始参考面的第1批参考点;每个格网中高程最大的点标记为植被最高处的点,作为DSM初始参考面的第1批参考点。
图2 LiDAR数据应用处理流程Fig.2 LiDAR Data application process
水平和垂直方向分别依次移动格网起始点位置(例如1/n格网间隔),重新建立格网间距2 m的网格矩阵,重复以上步骤,获得DEM和DSM初始参考面的第2批参考点;水平和垂直方向分别移动n-1次后,获得DEM和DSM的离散值初始参考面。
离散DEM和DSM初始参考面提取完成后,利用IDW插值方法,插值格网大小选择1/n格网间隔,得到连续的DEM和DSM初始参考面模型E(i,j),T(i,j)。
图3 点云转化到格网示意图像Fig.3 Point cloud convert to grid diagram
2.1.2 DSM点与DEM点精确分类
利用DEM和DSM初始参考模型,对初始的激光点云进行精确分类。
对于第一步中未被标记为DEM和DSM的点集合Z中的每个点,分别计算该点与E(i,j),T(i,j)初始参考面的差值:
其中,dZEn|dZTn分别为该点与DEM初始参考面、DSM初始参考面的高程差值,Zn为该点高程,Es(i,j)和Ts(i,j)为点Z对应DEM和DSM参考面网格(i,j)处3×3邻域内高程平均值。
计算公式如下:
设置地面点与植被点判断阈值,根据如下规则判断该点的分类:
2.1.3 插值生成DSM、DEM模型
集合Z中的所有点判断完毕后,地面点、植被点分别与对应类型的初始参考面标记点合并,形成DSM、DEM点集合,最后通过反距离加权插值(IDW)形成研究区域最终的DSM、DEM模型。
冠层高度模型(Canopy height model,CHM)是一个代表植被高度的连续数字数据集,是一个表达植被距离地面高度的表面模型,能够提供树木冠层的水平和垂直分布情况[44]。根据CHM的定义,它一般是由DSM与DEM进行差值运算生成。森林调查中很多的植被参数都可以从CHM中直接或间接获得,例如树高、冠幅、郁闭度、蓄积量和生物量等,但是获得的这些参数通常会比实际值低。
本研究利用GIS中的栅格空间分析原理,对DSM与DEM进行差值运算生成CHM,通过去除了地面点后,对CHM进行统计得出平均树高。基于平均树高估计冠幅、郁闭度、蓄积量和生物量等指标作为下一步的研究内容。
基于LiDAR的平均树高估算分两个步骤,首先估算研究区域内林木数量及单株位置分布,然后根据CHM获得每株单树对应位置的高度值作为该树的树高,最后计算所有树高的平均值,得到整个区域的平均树高[47-48]。
2.3.1 单株树木位置与数量的估算
以CHM为基础,在一定半径范围内查找最高值处的位置点,作为一株单树的估计位置,该点的高度作为该树的树高。
2.3.2 研究区域的平均树高估算
由于树枝、树叶、地面植被等的影响,激光点很难直接打到地面,在到达地面前就被提取反射回去,造成LiDAR地面高度有一定的误差。本研究在估计平均树高时,设定树木最低值HL,统计大于等于最低值HL的树高的平均值作为平均树高估计值,以减少上述因素对林木高度统计的影响。
实验区域的LiDAR点云数据经过系统的处理,本次选择阈值大小为0.2 m。将最后所得到的激光点进行内插,格网间距为2 m,即可得到较为精确的DEM、DSM:经过DEM提取、DSM提取、IDW插值、CHM生成、平均树高计算等步骤,得到了实验区域的DEM、DSM、CHM、单株树木位置估计图等。如下图所示:
图4 原始点云数据(Z值渲染)Fig.4 Raw point cloud data (Z-value rendering)
图 5 DSM渲染图像及图例Fig.5 DSM rendering and illustration
根据CHM提取平均树高,与核实采样点进行对比分析。通过设置不同的最低树高,平均树高的估计结果如表1所示:
图6 DEM渲染图像及图例Fig.6 DEM Rendering and illustration
图7 CHM渲染图像及图例Fig.7 CHM Rendering and illustration
图8 单株树木位置与数量的估算Fig.8 Estimation of location and quantity of single tree
表1 LiDAR点云数据提取平均树高Table 1 LiDAR point cloud data extraction average tree height (m)
图9 不同最低树高设置下LiDAR平均树高计算Fig.9 Calculation of LiDAR average tree height under different minimum tree height setting
通过表1可以得出以下结论:LiDAR点云获得的树高比实际值低。由于LiDAR数据树顶的错失及点云的低密度,树高被低估0~2 m。按以上统计结果,通过最低树高的值不断增加,实测值与估计值的差值逐渐缩小,最低树高设置为6 m时,差值为0 m。此表分析结果可作为经验值,用于后期LiDAR数据树高改正。
造成LiDAR数据估计的树高比实际低的可能原因分析如下:
(1)树顶采集点缺失。由于激光采样密度的原因,树顶可能不会被激光打到,形成采集数据树顶缺失,造成估计树高偏低。
(2)激光无法完全穿透树枝树叶到达地面。虽然激光有一定的穿透性,但是对于在热带森林地区树枝茂密、林下低矮植被繁盛的情况下,较难到达实际地面,大部分被提前反射,同样造成估计树高偏低。
基于以上因素,后期推广利用LiDAR数据估计平均树高时,应进行相关实地外业核查,获得相关改正数,使估算结果更加精确。
本研究提出了一种大面积森林区域机载激光扫描仪冠层表面重建和森林参数预测方法,改进了激光雷达数据处理的动态搜索窗口算法,并用C#、ArcGIS等软件进行了编码实现。研究选择三亚市福万水库北侧山地50 hm2的热带森林作为试验区,并将估算结果与现场实测数据进行了比较验证。结果表明,经过算法的改进,大面积森林的DEM、DSM模型均可以从LiDAR数据中快速提取,且能够较好地重建森林冠高模型;并从CHM中提取估计森林平均树高,发现LiDAR点云获得的树高比实际值低,实际树高被低估0~2 m,但最低树高的值在不断增加时,其实测值与估计值的差值也在逐渐缩小。
本研究利用试验区低密度LiDAR数据开展了大面积森林的数字地面模型(DEM)、数字表面模型(DSM)和冠高模型(CHM)的提取方法研究,发现激光雷达技术在估测大面积山区热带森林参数方面是可行的,未来可以将其集中应用在更多的参数如森林冠幅、冠层密度、生物量及蓄积量等指标的估测和改进方法中。机载激光雷达技术是目前提取信息技术途径的关键之一,我国在这一研究领域仍存在许多技术难题,本研究中在数据源使用和研究方法选择上仍有待进一步深入探索完善,还存在许多需要解决的问题。主要体现在以下几个方面:
(1)研究所采用激光点云的密度仅为1.51个/m2,密度较低。对比高密度(大于 3个/m2)点云数据在树冠模型等林木参数提取时效果稍差,说明点云密度对林木参数提取精度仍有一定的影响,然而如果一味地增加点云密度则数据获取成本也会大大增加,且点云数据的采样密度和数据处理过程中各种参数的设置均会对结果产生一定的影响[47]。因此,定量地分析点云密度与森林冠高模型、林分平均高估测精度之间的关系,对优化在森林资源调查应用中LiDAR系统参数的设置调配具有一定的指导作用。
(2)不同的森林类型、海拔、坡度、坡向等均与林木参数有密切的相关,在进行森林冠层表面重建和参数预测反演时,应将海拔、坡度、坡向等因素考虑到模型中。
(3)需要不断的改进与更新点云提取方法和LiDAR数据算法,使其具备更高的识别精度与效果,提高估测精度,且同时兼具低时间成本的算法,大大增加LiDAR数据替代部分人工样地调查的可能性。
[1]Masek JG, Hayes DJ, Hughes MJ, et al.The role of remote sensing in process-scaling studies of managed forest ecosystems[J].Forest Ecology and Management, 2015,355: 109-123.
[2]Ginzler C, Hobi M.Countrywidestereo-image matching for updating digital surface models in the framework of the Swiss national forest inventory[J].Remote Sensing, 2015,7(4): 4343-4370.
[3]Waser L T, Fischer C, Wang Z,et al.Wall-to-wall forest mapping based on digital surface models from image-based point clouds and a NFI forest definition[J].Forests, 2015,6(12): 4510-4528.
[4]刘清旺,李世明,李增元,等.无人机激光雷达与摄影测量林业应用研究进展[J].林业科学,2017,53(7):134-148.
[5]Naesset E.Estimating timber volume of forest stands using airborne laser scanner data[J].Remote Sensing of Environment,1997, 61(2): 246-253.
[6]Nilsson M.Estimation of tree heights and stand volume using an airborne lidar system[J].Remote Sensing of Environment, 1996,56(1): 1-7.
[7]吴 楠,李增元,廖声熙,等.国内外林业遥感应用研究概况与展望[J].世界林业研究, 2017,30(6):34-40.
[8]吕国屏,廖承锐,高媛赟,等.激光雷达技术在矿山生态环境监测中的应用[J].生态与农村环境学报,2017,33(7):577-585.
[9]胥 喆,舒清态,杨凯博,等.星载激光雷达在林业上的应用研究进展[J].福建林业科技,2017,44(1):141-148.
[10]赖旭东,李咏旭,陈佩奇,等.机载激光雷达技术现状及展望[J].地理空间信息,2017,15(8):1-4.
[11]李增元,刘清旺,庞 勇,等.激光雷达森林参数反演研究进展[J].遥感学报,2016,20(5):1138-1150.
[12]王东亮.低空遥感数据处理与非生长季草地生物量反演研究[D].北京:中国农业科学院,2016.
[13]邢艳秋,姚松涛,李梦颖,等.基于机载全波形LiDAR数据的森林地上生物量估测算法研究[J].森林工程,2017,33(4):21-26.
[14]姚兴成,曲恬甜,常文静,等.基于MODIS数据和植被特征估算草地生物量[J].中国生态农业学报,2017,25(4):530-541.
[15]庞 勇,蒙诗栎,李增元,等.机载高分辨率遥感影像的傅氏纹理因子估测温带森林地上生物量[J].林业科学,2017,53(3):94-104.
[16]邢元军,刘晓农,宋亚斌,等.国产高分辨率遥感影像融合方法比较与分析[J].中南林业科技大学学报,2016,36(10):83-88.
[17]林岳峰.光学与雷达遥感数据联合反演森林生物量方法研究[D].成都:电子科技大学,2016.
[18]刘海清.森林蓄积量遥感估测的应用研究[D].西安:西安科技大学, 2009.
[19]曹明兰,张力小,王 强.无人机遥感影像中行道树信息快速提取[J].中南林业科技大学学报,2016,36(10):89-93.
[20]郄广平.高分辨率遥感影像的森林结构参数反演[D].长沙:中南林业科技大学, 2011.
[21]谢进金.基于遥感抽样技术的平南县森林蓄积量估测[D].长沙:中南林业科技大学, 2011.
[22]曾 晶,张晓丽.高分一号遥感影像下崂山林场林分生物量反演估算研究[J].中南林业科技大学学报, 2016, 36(1):46-51.
[23]涂云燕.森林蓄积量遥感估测研究[D].北京:北京林业大学,2013.
[24]王雪军,马 炜,黄国胜,等.基于遥感大样地抽样调查的森林面积监测[J].北京林业大学学报, 2015, 37(11): 1-9.
[25]韦雪花.轻小型航空遥感森林几何参数提取研究[D].北京:北京林业大学, 2013.
[26]刘清旺,谭炳香,胡凯龙,等.机载激光雷达和高光谱组合系统的亚热带森林估测遥感试验[J].高技术通讯,2016,26(3):264-274.
[27]林 怡,季昊巍,叶 勤,等.基于LiDAR点云的单棵树木提取方法研究[J].计算机测量与控制,2017,25(6):142-147..
[28]陈 健.基于地基激光雷达的不同森林类型单木胸径与树高提取[D].合肥:安徽农业大学,2016.
[29]刘鲁霞,庞 勇,李增元,等.基于地基激光雷达的亚热带森林单木胸径与树高提取[J].林业科学,2016,52(2):26-37.
[30]陶江玥,刘丽娟,丁友丽,等.高光谱和LiDAR数据融合在树种识别上的应用[J].绿色科技,2017(8):212-214.
[31]刘 峰,谭 畅,王 红,等.基于LiDAR的亚热带次生林林窗对幼树更新影响分析[J].农业机械学报,2017,48(3):198-204.
[32]高 婷,李卫忠,赵鹏祥,等.基于ArboLiDAR的大野口林区森林参数估测[J].西北林学院学报,2017,32(4):172-177,223.
[33]冯静静,张晓丽,刘会玲,等.基于灰度梯度图像分割的单木树冠提取研究[J].北京林业大学学报,2017,39(3):16-23.
[34]吴金汝,王新闯,陆凤连,等.利用波形去噪算法的森林最大冠顶高估测研究[J].测绘科学,2017,42(7):22-28.
[35]孙 浩,刘晋浩,黄青青,等.基于二维激光扫描的立木胸径计算方法性能分析[J].农业机械学报,2017,48(8):186-191.
[36]张卫正,董寿银,王国飞,等.基于机载LiDAR数据的林木冠层投影面积与体积测量[J].农业机械学报,2016,47(1):304-309.
[37]段祝庚,肖化顺,袁伟湘,等.基于离散点云数据的森林冠层高度模型插值方法[J].林业科学,2016,52(9):86-94.
[38]Lim K, Treitz P, Wulder M,et al.LiDAR remote sensing of forest structure[J].Progress in Physical Geography, 2003,27(1):88-106.
[39]Magnussen S, Boudewyn P.Derivations of stand heights from airborne laser scanner data with canopy-based quantile estimators[J].Canadian Journal of Forest Research, 1998, 28(7):1016-1031.
[40]Schardt M, Ziegler M, Wimmer A,et al.Assessment of forest parameters by means of laser scanning[J].International Archives of Photogrammetry Remote Sensing and Spatial Information Sciences, 2002, 34(3/A): 302-309.
[41]Sileshi G W.A critical review of forest biomass estimation models, common mistakes and corrective measures[J].Forest Ecology and Management, 2014, 329: 237-254.
[42]穆喜云,张秋良,刘清旺,等.基于激光雷达的大兴安岭典型森林生物量制图技术研究[J].遥感技术与应用,2015(2): 220-225.
[43]Hoffman R, Jain A K.Segmentation and classification of range images[J].IEEE Transactions on Pattern Analysis and Machine Intelligence, 1987, (5): 608-620.
[44]Wulder M A, Bater C W, Coops N C,et al.The role of LiDAR in sustainable forest management[J].The Forestry Chronicle, 2008,84(6): 807-826.
[45]陈 浩.激光雷达测量技术及其应用[J].电子技术与软件工程,2017(3):125-125.
[46]滕继峣,秦 凯,汪云甲,等.基于激光雷达观测的大气边界层自动识别局部最优点算法[J].光谱学与光谱分析,2017,37(2): 361-367.
[47]梁 宏,刘长长.机载激光雷达点云数据处理研究[J].科技视界,2017(9):144.
[48]庄永健,冯仲科,李亚藏,等.激光雷达技术测树方法原理与应用[J].林业资源管理,2016(6):116-119.