利用小光斑激光雷达估测林分参数和生物量1)
2012-06-13刘东起范文义李明泽
刘东起 范文义 李明泽
(东北林业大学,哈尔滨,150040)
森林资源地面调查的传统方法周期比较长,数据的时效性较差[1],遥感技术的出现和迅速发展,为准确、快速地获取大范围森林资源调查数据提供了一种潜力巨大的技术手段。利用激光作为遥感设备可追溯到30多年以前,但机载激光雷达测量技术直到最近十几年才取得重大进展,并研制出精确可靠的雷达测量传感器。国外已经就激光雷达在林业上的应用进行了大量成功的尝试[2-4],并估测了树高、胸高断面积及林分蓄积量等林分特征[5-8]。在林分尺度上,激光雷达已成功用于大尺度上的林分平均高、平均胸径、冠层垂直结构、冠层体积、郁闭度、地上生物量、蓄积量、胸高断面积、叶面积指数等参数的估计[9-11]以及林区变化的检测[12],相关算法经过数十年的发展也较为成熟。如Sr-Onge等[13]用ALTM获取的小光斑RADAR数据估测单株木树高,采用Laplacian of Gaussian(LoG)算法并结合空间过滤和边缘提取方法提取单株树冠,并用得到的树冠多边形中RADAR回波最大值与实际树高建立模型。一些国家也正在考虑用RADAR技术进行森林资源调查[14-17]。近几年,国内也逐渐开展了激光雷达在林业上的应用研究,如国家“948”项目“高精度激光雷达树高测量技术引进”,国家“863”项目“机载激光雷达数据处理软件平台”等[18]。覃驭楚等[19]利用Axelsson提出的渐进式不规则三角网加密法,分别采用反距离权重方法、样条插值方法以及普通克里格方法对小光斑激光雷达的多次回波数据进行滤波,对冠层与林下高程分别进行了重建并做了精度检验,研究结果证明:RADAR数据获取的树高可以达到实地调查的精度,其误差绝对平均值均在1.4 m以内。李奇等[20]依据雷达等式,通过对波形数据的高斯分解来模拟激光回波,计算树高、林下地形、冠层体积、地表反射率、植被反射率以及森林郁闭度等森林参数,这种算法对森林结构参数获取的自动化程度较高,基本不需要人工干预,但其没有实测参数进行对比,无法量化精度。唐菲菲等[21]采用融合有标记约束的分水岭分割和流域跟踪分割两种图像分割方法,以从机载激光雷达数据中分离单株木,提取单株木树高和树冠大小信息为目的,提出了一种新的单株木识别思路,并以美国某地区实地采集激光雷达数据为例对该方法进行了验证,试验结果表明:该方法通过增加边缘检测范围的约束条件,能够有效避免过分割现象,并通过使用约束条件,减少了在其检测范围内的目标数量,从而避免不必要的检测干扰,较传统方法能快速准确地识别单株木。激光雷达的自身优势也使其受到越来越广泛的应用,但关于小光斑激光雷达的研究还比较少,且精度也较低,同时又局限于单一遥感数据的研究。笔者则利用激光雷达数据的特点与TerraScan的强大功能对小兴安岭凉水地区的航拍数据进行研究,估测了树高等林分参数以及不同林分类型的生物量,为后续的结合多遥感数据研究做准备。
1 研究区域概况
研究区选在黑龙江省伊春市带岭区凉水国家级自然保护区。地理位置东经 128°47'8″~128°57'19″,北纬47°06'49″~47°16'10″。保护区面积12 133 hm2,森林总蓄积量170万m3,森林覆被率98%。保护区植物区系属泛北极植被区、中国—日本森林植物亚区东北地区、长白植物亚区小兴安岭南部区。本区的地带性植被是以红松(Plnus koralensls)为主的温带针阔叶混交林,属典型阔叶—红松林分布亚区。该区地处欧亚大陆东缘,小兴安岭南坡达里带岭支脉的东坡,海拔高度280~707 m,为典型的低山丘陵地貌,最高山峰为保护区北部的岭来东山,海拔为707 m,由北向南逐渐降低至该区的东南端,海拔仅有280.0 m,山地坡度为10°~15°。境内地带性土壤为山地暗棕壤,即由地带性上坡暗棕色森林土(简称暗棕壤坡)的亚类与非地带性土壤如草甸土、沼泽土及其亚类所组成。局部地段可出现20°~40°的陡坡。区内主要河流为凉水河,由北向南贯流全境。该区具有明显的温带大陆性季风气候特征,年平均气温-0.3℃ ,年平均降水量676.0 mm,年平均相对湿度78%,无霜期100~120 d。气候特点:春天来的迟缓,降水少,多大风;夏季短促,气温较高,降雨集中;秋季降温急剧,常出现早霜;冬季漫长,严寒干燥而多风雪。研究区地理位置及飞行区域见图1。
图1 研究区地理位置及飞行区域
2 研究方法
2.1 雷达数据获取
此次飞行使用的激光雷达系统为 LiteMapper5600,激光扫描仪为RieglLMS—Q560,CCD相机型号为DigiCAM—H/22。该系统集激光测距、全球定位系统(GPS)和惯性导航系统(INS)于一体,包括一个单束窄带激光器和一个接收系统。飞行地速180 km/h,总覆盖面积约为 208.8 km2,数据获取当天天气晴朗无(少)云。飞行航高1 000 m,设计飞行航线26条,激光脉冲长度3.5 ns,波长1 550 nm,最大频率100 kHz,脉冲重复频率为50 kHz,采样间隔1 ns,激光离散角0.5 mrad,扫描角±30°,地表定位精度水平0.1 m,垂直精度0.15 m,获得的平均点云密度>2 point/m2,部分小范围飞行区达到>5 point/m2(航高300 m)。此外,为保证POS系统的定位精度,飞行前在测区附近的已知大地控制点上建立了GPS基准站。基准站上架设高精度的GPS信号接收机,该接收机与机载POS设备同步记录,进行动态相位差分定位。
2.2 地面数据获取
2009年,与飞行同步进行了部分地面数据的获取,并于2010年8月对飞行区域进行了地面数据的补充调查。文中用到的每木检尺数据是2010年8月获取的。地面调查时沿着山坡走势布设了51块样地,视样地类型分别设置边长为20 m×20 m、20 m×30 m、30 m×30 m的方形样地进行每木检尺:胸径尺测量单木胸径;用DGPS定位样地中心点坐标,罗盘确定单木角度,皮尺丈量单木距离中心点长度并凭借角度、距离与中心点位置坐标推算单木位置;皮尺测量东西及南北两个主方向冠幅;激光测距/测高仪(TruPulse360)测量树高,该仪器测距精度300 m以内误差为-0.3 m,角度测量精度为±0.25,保证了高度测量的准确性。此次地面调查共测1 804棵样木,去除测量值明显错误的285株林木和缺少单木位置的25株林木,最后剩余1 494株,树种主要包括:红松、云杉(Picea asperataMast.)、冷杉(Abies)、落叶松(LarixMill.)、白桦(Betula platyphyllaSuk.)等。其基木统计量均值如表1。
表1 单木基本统计量
3 结果与分析
3.1 机载雷达数据处理
采用瑞典Peter Ax-lesson提出的不规则三角网法在 Terrscan V4.005(Terrasolid,Helsinki,Finland)中通过迭代、孤立点、阈值等算法对点云产品进行地面点和非地面点的分离[22-23]。该方法先选择种子点构建一个初始的地面三角网,剩余的点为非地面点,然后逐渐从非地面点中选择满足条件的点向初始三角网中添加构成新的网,迭代计算直到所有的点被分为地面点和非地面点。再对非地面点数据子集操作,采用孤立点、阈值等算法进行非地面点的分类。自动滤波完成后有些点的分类是错的,需要手动进行进一步分类。但需尽量减少手动操作,以节省时间及成本[22]。
对已分类点云数据中的地面点通过不规则三角网(TIN)插值运算生成数字高程模型DEM(Digital Elevation Model)见图2,首回波点插值生成数字表面模型DSM(Digital Surface Model)见图3,DSM与DEM进行栅格差值运算得到高程归一化后的数字冠层高度模型CHM(Canopy Height Model)见图4。经过上述操作即可实现对点云数据的分类,提取数字高程模型DEM、数字表面模型DSM与树冠高模型CHM。在此基础上可估测单木树高、冠幅、胸径等参数,亦可估测林分平均高、林分平均胸径以及蓄积量、生物量等参数。
图2 提取的部分DEM
图3 提取的部分DSM
图4 部分CHM的三维渲染图
3.2 树高和冠幅的估测与分析
3.2.1 树高的估测
采用外业调查确定单木位置,在单木位置周围1.5 m(考虑到郁闭条件下的单木定位精度)范围内搜索CHM中对应的树冠顶点[24],取范围内最大值。估测单木树高、冠幅等参数。
剔除奇异点后的实测树高与估测树高散点图如图5所示。
图5 实测树高与估测树高散点图
对实测树高与估测树高的回归分析的结果显示,平均百分比误差(Mean Percent Error)为-7.119 601%;平均绝对百分比误差(Mean Absolute Percent Error)为9.041 874%;精度(Precision)为 97.642 19,表明RADAR数据很好地估测了样地内的树高。
3.2.2 冠幅的估测
实测树高与实测冠幅的散点图如图6所示。
图6 实测树高与实测冠幅拟合图
从图6中可以看出,实测树高与实测冠幅之间的关系很稳定,且存在线性关系。为此本文建立了树高和冠幅的回归模型并推导出参数,然后在此基础上对树冠进行了估测。回归模型为:
式中:CD表示冠幅;H为树高;k为斜率;b为截距。通过实测树高与冠幅的线性回归分析确定方程的参数,k=0.107,b=3.492 6。将估测的树高代入(1)式即可以得到冠幅的估测值。实测冠幅与估测冠幅的回归分析结果为:平均百分比误差1.251 971%;平均绝对百分比误差9.568 91%;精度93.175 76。回归模型的检验精度达93%以上,表明该方法可以很好地估测样地内单木的冠幅。
3.3 胸径与生物量的估测与分析
3.3.1 胸径的估测
对胸径的估测采用树高曲线方程
进行推导,得胸径
式中:h为树高;a0、a1为方程参数,依据外业实测数据建模计算得a0=4.098,a1=2.878;d为直径。将估测树高作为已知,带入参数公式(3)计算估测胸径,剔除部分奇异点后的实测胸径与估测胸径的散点图如图7所示,实测胸径与估测胸径的回归分析结果为:平均百分比误差-11.835 55%;平均绝对百分比误差 10.756 82%;精度 92.393 15。
从图7可以看出,剔除异常点后实测胸径与估测胸径具有明显的线性关系;且理论精度为92%,能够达到应用要求。
图7 实测胸径与估测胸径的散点图
3.3.2 生物量的估测
参照Ross Nelson等[25]人的方法,采用统计分析的方式通过建立回归模型估测生物量,由于外业调查选取的样地林分类型均能很好代表其所在小班的林型特点,而所选取的标准木亦均为干型匀称的林木且均按照林分类型选取,故可作为标准木测定材积,为此文中采用平均标准木法测算不同林型的蓄积,并最终援引范文义等[26-29]对小兴安岭评价时所使用的生物量转换因子法利用RADAR数据对林分生物量和蓄积量等进行了估测。所需公式为:
式中:vi为样木材积;f1.3为胸高形数,按不同树高依圆锥体查表而定;g1.3为胸高断面积;Vj为不同林型各自的平均蓄积量(m3/hm2);Gj为不同林型各自的总断面积;n为样木株树;Sj为不同林分类型各自的总面积;B为单位面积生物量(t/hm2);aj和bj为转换参数,依不同林分类型确定[29]。不同林分类型估测的平均单位蓄积量与生物量见表2。
表2 不同林分类型估测的平均单位蓄积量、生物量及精度
4 讨论
树高是一个十分重要的林分参数,同时也是激光雷达所能够估测的最直接的林分参数,树高的估测精度直接影响胸径、蓄积量以及生物量的林分参数的估测精度。文中结合雷达数据自身优点的同时充分利用了TerraScan的强大功能,估测的树高精度达到97%以上。对于冠幅的估测,由于没有采用传统的与大比例尺航空相片相结合的办法进行树冠边缘分割,而是直接采用外业调查确定单木位置的方式通过回归模型进行估测,精度在93%以上,能够满足实际应用的需求。
胸径的估测是在树高估测的基础上,运用树高曲线方程推导而得。由于树木在实际生长的各阶段呈现不同的规律,其树高与冠幅不会像估测值那样理想,因为冠幅并不是一直随树高的生长而生长,然而,由于采用了大量的实测数据为依托能够覆盖研究区域不同林型、不同树种在不同阶段所呈现出的树高—冠幅之间的关系,所以估测的结果也比较理想,冠幅估测的检验精度为92%,亦可以满足应用需求。
[1]庞勇,于信芳,李增元,等.星载激光雷达波形长度提取与林业应用潜力分析[J]林业科学,2006,42(7):137-140.
[2]赵丽琼,张晓丽,孙红梅.激光雷达数据在森林参数获取中的应用[J].世界林业研究,2010,23(2):61-64.
[3]Marr D,Hildrith E.Theory of edge detection[J].Proceedings of the Royal Society of London.Series B,Biological Sciences,1980,207:187-217.
[4]马利群,李爱农.激光雷达在森林垂直结构参数估算中的应用[J].世界林业研究,2011,24(1):41-44.
[5]Reese H,Mats N,Granqvist P T,等.运用卫星数据和来自国家森林清查的野外数据进行全国森林变量估计[J].AMBIO-人类环境杂志,2003,32(8):539-544.
[6]Ross N,Krabill W,Maclean G.Detemining forest canopy characteristics using airborne laser data[J].Remote Sensing of Environment,1984,15(3):201-212.
[7]Naesset E,Bjerknes K O.Estimating tree heights and mumber of stems in young forest stands using airborne laser scanner data[J].Rmote Sensing of Envirorunent,2001,78(3):328-340.
[8]Means J E,Acker S A,Fitt B J,et al.Predicting forest stand characteristics with airborne scanning LIDAB[J].Photogrammetric Engineering and Remote Sensing,2000,66(11):1367-1371.
[9]Mallet C,Bretar F.Full-waveform topgraphic lidar:State-of-theart[J].Journal of Photogrammetry and Remote Sensing,2009,64(1):1-16.
[10]Lefsky M A,Cohen W B,Acker S A,et al.Lidar remote sensing of the canopy structure and biophysical properties of Douglasfir western hemlock forests[J].Remote Sensing of Environment,1999,70(3):339-361.
[11]Dubayah R O,Drake J B.Lidar remote sensing for forestry[J].Journal of Forestry,2000,98(6):44-46.
[12]Duong H V,Lindenbergh R,Pfeifer N,et al.Single and two epoch analysis of ICESat full-waveform data over forested areas[J].International Journal of Remote Sensing,2008,29(5):1453 -1473.
[13]庞勇,李增元,陈尔学,等.激光雷达技术及其在林业上的应用[J].林业科学,2005,41(3):129-136.
[14]马茂江,张文,万国礼.德国森林资源调查与监测对我国的启示[J].四川林勘设计,2008(2):9-13.
[15]Clark N A,Wynne R H,Schmoldt D L.A review of past research on dendrometers[J].Forest Science,2000,46(4):570-576.
[16]殷鸣放,谭希彬.国际模式森林的研究进展[J].辽宁林业科技,2001(4):24-27.
[17]刘平.美国森林植被模拟系统(FVS)在北京地区人工林上的应用研究[D].北京:北京林业大学,2007.
[18]徐光彩.机载LIADR波形数据处理及分类研究[D].南京:南京林业大学,2010.
[19]覃驭楚,吴运超,牛铮,等.基于小光斑激光雷达数据的稀疏林地冠层高程重建[J].自然资源学报,2008,23(3):507-510.
[20]李奇,马洪超,邬建伟,等.机载小光斑LIDAR的森林参数评估[J].林业资源管理,2008(1):74-81.
[21]唐菲菲,刘星,张亚利,等.基于机载激光雷达数据识别单株木的新方法[J].遥感技术与应用,2011,26(2):196-202.
[22]张煜,窦延娟,张晓东.机载激光雷达数据采集及数据处理[J].长江科学院院报,2010,27(1):13-16.
[23]黄金浪.基于TerraScan的LiDAR数据处理[J].测绘通报,2007(10):13-17.
[24]刘清旺.机载激光雷达森林参数估测方法研究[D].北京:中国林业科学研究院,2009.
[25]Nelson R,Short A,Valenti M.Measuring biomass and carbon in Delaware using an airborne profiling lidar[J].Scandinavian Journal of Forest Research,2004(19)6:500-511.
[26]李峰,刘桂英,王力刚.黑龙江省森林碳汇价值评价及碳汇潜力分析[J].防护林科技,2011(1):87-88.
[27]王维枫,雷渊才,王雪峰,等.森林生物量模型综述[J].西北林学院学报,2008,23(2):58-63.
[28]刘志华,常禹,陈宏伟.基于遥感、地理信息系统和人工神经网络的呼中林区森林蓄积量估测[J].应用生态学报,2008,19(9):1891-1896.
[29]杨金明,范文义.小兴安岭主要树种生物量的理论模型[J].东北林业大学学报,2011,39(3):46-48.