林业科学  2018, Vol. 54 Issue (10): 98-107   PDF    
DOI: 10.11707/j.1001-7488.20181012
0

文章信息

李欢, 李明泽, 范文义, 王斌.
Li Huan, Li Mingze, Fan Wenyi, Wang Bin
基于机载激光雷达的林隙结构参数提取
Canopy Gap Structure Parameters Extraction Based on Light Detection and Ranging(LiDAR)
林业科学, 2018, 54(10): 98-107.
Scientia Silvae Sinicae, 2018, 54(10): 98-107.
DOI: 10.11707/j.1001-7488.20181012

文章历史

收稿日期:2017-08-07
修回日期:2017-12-14

作者相关文章

李欢
李明泽
范文义
王斌

基于机载激光雷达的林隙结构参数提取
李欢, 李明泽, 范文义, 王斌     
东北林业大学林学院 哈尔滨 150040
摘要:【目的】运用激光雷达技术识别提取林隙面积、边界木高和形状指数,为林隙结构参数调查提供技术支持。【方法】以机载激光雷达数据为数据源,以东北林业大学帽儿山实验林场内设置的5块林隙调查样地中的54个林隙为研究对象,通过对比普通克里金插值、反距离权重插值、样条插值、自然临近插值和局部多项式插值确定最优插值方法得到冠层高度模型,运用交替贯序滤波法对冠层高度模型进行林隙识别提取。【结果】通过对比5种插值方法的RMSE发现,普通克里金法适合插值DEM,反距离权重法适合插值DSM;采用交替贯序滤波法识别提取林隙时,林隙识别率为92.6%,运用配对检验法对提取的林隙面积、边界木高与野外调查数据进行检验,基于激光雷达技术提取的林隙面积和边界木高与野外实测值呈较强线性关系,R2分别为0.983和0.737,其中提取的林隙面积与实测值平均相对误差为15.78%,林隙边界木高与实测值平均相对误差为11.94%。【结论】基于机载激光雷达数据采用交替贯序滤波能有效识别林隙且能准确提取其结构参数;坡度、坡位及冠层结构复杂程度都会影响冠层高度模型的插值精度;林隙面积提取受林隙边界形状和林层结构影响,林隙边界木高随着样地地形、坡度增加误差也随之增大。机载激光雷达技术具有获取地形和树木三维结构信息的优势,可为林隙识别及其结构参数提取提供新的遥感方法。
关键词:机载激光雷达    林隙结构参数    数学形态学    CHM    
Canopy Gap Structure Parameters Extraction Based on Light Detection and Ranging(LiDAR)
Li Huan, Li Mingze , Fan Wenyi, Wang Bin    
College of Forestry, Northeast Forestry University Harbin 150040
Abstract: 【Objective】In order to provide technical support for the investigation of gap structure parameters, the LiDAR technology was used in this study to identify the gap structure parameters(including gap size, canopy height and gap shape index).【Method】Taking the airborne LiDAR as the data source, the 54 gaps in the five gap investigate plots were set up in the Mao'ershan experimental forest farm of Northeast Forestry University as the study area, The canopy height model was obtained by comparing the five interpolation method (ordinary kriging, inverse distance weight, splines, natural neighbor, local polynomial)to determine the optimal interpolation method. The canopy height model was used to identify and extract the gaps by alternative sequential filter.【Result】By comparing the RMSE of the five interpolation method, it is found that the ordinary kriging method is suitable for interpolating DEM, and the inverse distance weight method is suitable for interpolating DSM. The gap recognition rate was 92.6% when the gap was extracted by alternating sequential filter. Using the paired test method to check the gap size, canopy height with field survey data, linear correlation was observed between the values of LiDAR estimation and field investigation, and the R2 values of gap size and canopy height case were 0.983 and 0.737, respectively. Compared with field survey data, the average relative error of extracted gap size was 15.78% and canopy height was 11.94%.【Conclusion】Using alternating sequential filter based on airborne lidar data can effectively identify the gap and accurately extract its structural parameters. Interpolation accuracy of the canopy height model is affected by the slope, slope position and the complexity of the canopy structure; gap size is affected by boundary shape and forest structure, canopy height as the terrain gradient increases, the error increases. The airborne LiDAR technology has the advantage of acquiring three-dimensional structural information of terrain and trees, this study provides a new remote sensing method for the identification of gaps and the extraction of gap structure parameters.
Key words: light detection and ranging(LiDAR)    gap structural parameters    mathematical morphology    CHM    

Watt(1947)提出林隙的概念以来,与之相关的研究逐渐成为林学和生态学研究的热点和重点之一。林隙是一种经常发生的、源于自然和人为的树冠干扰,即主林层树木因为某种生物内因(如树木的年龄和根系特征)、气候外力(如大风、雪、水灾、泥石流等)或者二者综合作用而发生的树木翻倒、折干、虫害等事件,从而导致森林中相对连续的林冠层面出现间断,进而影响到林隙的环境因子变化、小尺度林分格局、群落更新演替过程的现象(费世民等,2009)。随着对林隙研究的不断深入,其概念也逐渐得到扩充。Runkle(19811982)在对美国东部原始老龄林进行研究的基础上,将林隙分为2类,即冠空隙与扩展林隙:冠空隙指直接处于林冠层空隙下的土地面积,也称为实际林隙(Worral et al., 1988);扩展林隙指由冠空隙周围树木树干围成的土地面积,包括冠空隙和边缘冠幅到树木树干基部围成的土地面积。当前这2种林隙概念都已被广泛应用,但是冠空隙的概念应用更多(Tyrrel et al., 1994),本研究中所指林隙为冠空隙。林隙面积、形状指数和边界木高是决定林隙环境异质性的3个林隙结构参数,本研究即针对这3个参数提取进行探讨。

传统测量林隙的方法多是在样方水平上进行,通过获取的实地数据得出调查结果,再推至整个景观或更大范围。传统测量方法具有一定局限性:1)样方范围很大时,数据收集耗时耗力;2)将实地数据推至整个景观,该过程忽略了空间异质性,不能很好解决尺度推绎问题。与传统林隙测量方法相比,利用遥感技术识别提取林隙及其结构参数(林隙面积、形状指数和边界木高)具有高效性和准确性。中等空间分辨率遥感数据在林隙识别中对较小林隙(30 m以下)的识别效果不好(Asner et al., 2004),高空间分辨率遥感数据的出现弥补了这一不足(Tanaka et al., 1997)。航空影像作为目前空间分辨率最高的遥感数据,在林隙识别方面具有广阔应用前景,但光学遥感仅能识别林隙而无法获取其结构信息;随着林业研究中对获取森林结构信息更精细的要求,激光雷达(LiDAR)技术被广泛应用。机载激光雷达是新兴的主动式遥感技术之一,能在多时空尺度上获取森林生态系统高精度的植被结构信息和三维地形特征(刘峰等,2015)。Zhang(2008)基于机载激光雷达技术,运用形态学滤波对美国大沼泽地国家公园红树林进行林隙识别,并将识别出的林隙与航空影像进行了对比;Gaulton等(2010)直接对点云数据进行操作识别林隙,并与插值得到树冠高度模型识别出的林隙进行了比较。以往研究只是通过识别出的林隙与航空影像或冠层高度模型中的林隙进行对比得出较高识别率,并未有样地尺度上林隙识别及结构参数的提取。运用机载激光雷达与航空影像能有效识别林隙并提取其结构参数,刘峰等(2015)以湖南省雪峰山南麓武冈林场为研究对象,利用机载激光雷达数据进行样地尺度识别并提取林隙结构参数,但在进行林隙提取时采用的是通过识别森林冠层求补集间接得出林隙的方法,并没有直接对林隙进行识别提取。

目前,基于遥感技术提取林隙大多仅对林隙识别进行分析探讨,并未涉及林隙结构参数的提取。鉴于此,本研究以黑龙江省帽儿山实验林场设置的5块林隙调查样地中3种森林类型(针叶林、阔叶林和针阔混交林)共54个林隙为研究对象,以激光雷达数据为数据源,运用形态学滤波法直接识别提取林隙结构参数,旨在:1)研究影响数字高程模型(digital elevation modle,DEM)和数字地表模型(digital surface model,DSM)插值精度的因素并选取最优插值方法;2)研究采用高度阈值及形态学滤波法能否准确识别出林隙;3)研究不同森林类型、不同坡位对林隙结构参数提取的影响。

1 研究区概况与数据来源 1.1 研究区概况

研究区位于黑龙江省尚志市帽儿山实验林场,45°15′—45°29′N,127°28′—127°43′ E,属长白山系支脉张广才岭西北部小岭余脉,林场总面积26 620 hm2,境内平均海拔428 m,最高海拔817 m,平均坡度14.2°,属低山丘陵缓坡地形,温带季风气候,夏季温暖多雨,冬季寒冷干燥,年降水量723.8 mm,土壤为典型暗棕壤;植被属长白山植物区系,原始地带性植被为阔叶红松(Pinus koraiensis)林,现为东北东部山区较典型的天然次生林区,次生林类型多样且具有代表性,群落类型有硬阔叶林、软阔叶林、针叶林和针阔混交林,平均林龄50年,主要树种有水曲柳(Fraxinus mandshurica)、胡桃楸(Juglans mandshurica)、紫椴(Tilia amurenisis)、榆(Ulmus pumila)、色木槭(Acer mono)、蒙古栎(Querus mongolica)、白桦(Betula phatyphylla)、山杨(Populus davidiana)、落叶松(Larix gmelinii)、樟子松(Pinus sylvestris var. mongolica)、红松等。

1.2 数据来源 1.2.1 遥感数据本研究所用机械激光雷达数据和正射影像数据获取时间为2015年9月13、14日,设计飞行航线31条,航线总长度798 km,飞行航高1.2 km,总覆盖范围

360 km2。机载激光雷达数据是由远程机载激光扫描仪LMS-Q680i获取的全波形数据,单条带点密度约2 pts·m-2,考虑多回波和重复覆盖,最大点云密度超过10 pts·m-2。LiDAR点云数据为LAS1.2格式,投影方式为横轴墨卡托(UTM)投影,参考椭球为WGS84,激光雷达数据包含每个点的坐标值、高度值、回波次数、航带和强度等信息。用于检验林隙的正射影像数据分辨率为0.2 m,投影方式与机载激光雷达数据一致。本数据由中国林业科学研究院LiCHy(LiDAR, CCD and Hyperspectral)机载遥感系统采集提供。

1.2.2 野外林隙调查数据

基于机载激光雷达扫描范围选取合适地区进行野外数据调查,调查时间为2015年9月,根据不同林分类型和坡位共设置10条宽20 m、长100~300 m不等的调查样地(表 1)。在调查样地中进行林隙调查,利用GPS记录并测量出林隙中心位置坐标和海拔;采用等角椭圆扇形八分法沿顺时针方向等间距测量出林隙中心到边缘的距离(朱教君等,2007),从而计算出林隙面积;采用激光测高器测量林隙边界木高。共调查林隙54个(表 2)。

表 1 林隙样地特征 Tab.1 Characteristics of the gap sampling plots
表 2 林隙调查数据 Tab.2 Gap survey data
2 研究方法 2.1 激光雷达数据预处理

激光雷达通过测量信号传播时间确定扫描仪与激光点的相对距离(梁欣廉等,2005),激光点云分类是预处理过程的重要步骤,如果点云未正确分类,将对之后的操作产生误差。常见的点云滤波方法有线性预测法、倾角滤波法和形态学滤波法等,要考虑区域地形与数据处理效率,地势陡峭、平缓及植被疏密都是要考虑的因素。本研究运用TerraScan软件的自适应TIN模型滤波(郑斌等,2010),对原始激光雷达点云数据进行分类得到数字高程模型(DEM),该算法可分离出激光雷达点云数据中的自然地形。在分类前需要剔除激光雷达点云数据中的噪声点,包括云、折射等造成的高程异常点。自适应TIN模型滤波通过反复迭代建立地面三角网来分离出地面点,首先从点云中选择一定密度的低点建立原始地表模型,因此在其运行之前先要剔除错误低点;再通过设置适当的迭代角度和迭代距离开始反复向上寻找点进行建模,由此将点云数据分为地面点云和非地面点云。

2.2 冠层高度模型的建立

冠层高度模型(canopy height model,CHM)也称树冠高度模型,在森林范围中,CHM是数字地表模型(DSM)与数字高程模型(DEM)的差值构成的。通过自适应TIN模型滤波将点云数据分为地面点云和非地面点云,地面点云用来生成DEM,机载激光雷达获取的首次回波点用来生成DSM。DEM和DSM建模时需要利用离散点云数据中的高程进行空间插值,实现这一目的的插值方法很多,本研究统一采用栅格建模。

栅格分辨率和所选取的插值方法都会影响CHM模型精度。对点云数据进行插值时设置的栅格分辨率与包含的点云个数呈线性关系,如果分辨率过大,容易造成栅格化过程中的信息丢失;反之,则会过度内插。本研究栅格分辨率按照点云数据中第1次回波点云密度来选取(刘峰等,2011),计算公式如下:

$ d = \sqrt {1/\lambda } 。$ (1)

式中:d为分辨率;λ为第1次回波点云密度。

地物表面模型插值是根据已知点高程对未知点高程进行估算的方法。根据点云数据密度和地形不同,不同插值方法得到的效果也各不相同。本研究通过普通克里金插值(ordinary Kriging,OK)、反距离权重插值(inverse distance weighted,IDW)、样条插值(splines,Spline)、自然临近插值(natural neighbor,NN)和局部多项式插值(local polynomial,LP)比较表 1中5块样地类型选取出适用的插值方法,采用交叉验证法比较不同插值方法的适用性,利用均方根误差(RMSE)评判不同插值方法的优劣(顾春雷等,2011):

$ {\rm{RMSE}} = \sqrt {\sum\limits_{i = 1}^n {{{({Z^{{\rm{inte}}}}-{Z^{{\rm{real}}}})}^2}/n} } 。$ (2)

式中:Zinte为插值得到的高程;Zreal为实际高程;n为用于交叉验证时的预测样点个数。

2.3 林隙识别 2.3.1 形态学滤波

数学形态学方法已在遥感图像处理多个方面取得了良好效果,如图像增强、图像分割、边缘检测和图像噪声处理等,其基本思想是利用构建的结构元素作为探针来获取分析图像目标形状特征和结构信息,腐蚀和膨胀是2个最基本的操作。假设CHM中的像素值为z[ik, jk],ij代表CHM像素的行和列,定义腐蚀和膨胀操作如下:

$ {\varepsilon _w}\left[{i, j} \right] = {\varepsilon _w}\left({z\left[{i, j} \right]} \right) = {\rm{mi}}{{\rm{n}}_{({i_k}, {j_k})}} \in w\left({z\left[{{i_k}, {j_k}} \right]} \right); $ (3)
$ {\delta _w}\left[{i, j} \right] = {\delta _w}\left({z\left[{i, j} \right]} \right) = {\rm{ma}}{{\rm{x}}_{({i_k}, {j_k})}} \in w\left({z\left[{{i_k}, {j_k}} \right]} \right)。$ (4)

膨胀可以填充图像中的小孔及图像边缘处的小凹陷部分,腐蚀可以消除小且无意义的物体,其他形态学操作都是基于这2个操作的组合形成的。形态学开、闭操作就是利用腐蚀和膨胀的先后顺序进行图像处理,开、闭操作如下:

$ {\gamma _w}\left[{i, j} \right] = {\gamma _w}\left({z\left[{i, j} \right]} \right) = {\delta _w}\left({{\varepsilon _w}\left({z\left[{i, j} \right]} \right)} \right); $ (5)
$ {\varphi _w}\left[{i, j} \right] = {\varphi _w}\left({z\left[{i, j} \right]} \right) = {\varepsilon _w}\left({{\delta _w}\left({z\left[{i, j} \right]} \right)} \right)。$ (6)

开操作是先进行腐蚀运算再进行膨胀运算,而闭操作正好相反。一般图像的滤波处理可通过一系列开-闭或闭-开交替的形态学滤波器实现,但是林隙大小不一,有些林隙较大,采用大结构元素的单一闭-开或开-闭滤波器会使图像中的目标产生偏差。本研究选取交替贯序滤波(alternative sequential filter,ASF)方法(Zhang, 2008),进行开和闭运算的迭代,从小结构元素开始逐渐增大结构元素的尺寸直到给定尺寸为止。对CHM进行交替贯序滤波根据待处理图像中林隙大小设置结构元素增大到的最大尺寸,尺寸过大将滤去较小林隙。

2.3.2 林隙识别验证

林隙本质上是森林冠层的某种状态,本研究从CHM的角度进行提取分析。试验主要分为3步:1)像元阈值设置,将获得的CHM进行二值图像转换,像元值小于等于5 m为0,5 m以上为1;2)林隙识别提取验证,对转换后的二值图像进行数学形态学滤波,滤去噪声并提取林隙;对提取出的林隙图像进行栅格转矢量处理,得到林隙多边形矢量图层(gaps polygon layer,GPL),采用0.2 m分辨率的遥感正射影像在GPL上结合林隙识别条件进行目视验证并勾绘出林隙面积,目视判读条件即林隙多边形面积在4~ 1 000 m2且多边形剖面与相邻树冠高度差在5 m以上,勾绘出的林隙多边形作为与LiDAR提取识别结果作对比的参考林隙;3)林隙结构特征参数提取,通过GPL得到林隙面积和周长;林隙边界木高通过单木点云的高程分位数反演,选择大于2/3林分高度的单木为林隙边界木。

根据Öhman等(2005)对林隙形状指数的定义,利用ArcGIS软件在得到的林隙矢量图层上计算出林隙形状指数:

$ {S_i} = {P_i}/2\;\sqrt {{\rm{ \mathsf{ π} }}{A_i}} 。$ (7)

式中:AiPi为林隙面积和周长。林隙形状指数越高,林隙内部受到的边缘效应越大。

2.4 精度验证方法

精度验证采用配对t检验和非参数检验中的符号秩检验。将每对数据求差值,对得到的数据进行正态性检验,服从正态分布的数据采用配对t检验;反之,采用符号秩检验。配对t检验实际与单样本t检验相似,即差值的样本均数d所代表的未知总体均数ud与已知总体均数u0=0的比较,检验统计量公式如下:

$ t = \frac{{\bar d-{u_d}}}{{{S_{\bar d}}}} = \frac{{\bar d-0}}{{{S_d}/\sqrt n }} = \frac{{\bar d-0}}{{{S_d}/\sqrt n }} = \frac{{\bar d}}{{{S_d}/\sqrt n }}, v = n - 1。$ (8)

式中:d为每对数据的差值;d为差值的样本均数;Sd为差值的标准差;Sd为差值样本均属的标准误;n为对子数;v为自由度。

符号秩检验是非参数检验中对位置参数进行推断常用的方法,其不但利用了观测值和零假设的中心位置之差的符号进行检验,而且利用了差值大小所包含的信息(吴丹宇等,2008)。

3 结果与分析 3.1 不同插值方法的比较

根据式(1)计算出不同样地插值所需分辨率,A、B、C、D、E 5块样地的分辨率分别为0.604、0.637、0.629、0.570和0.606 m。为统一分辨率,本研究均取分辨率为0.6 m。

利用每块样地总样本点数的90%进行插值,剩余的10%进行验证,根据验证精度选取合适的插值方法进行插值;分别比较普通克里金插值(OK)、反距离权重插值(IDW)、样条插值(Spline)、自然临近插值(NN)、局部多项式插值(LP)在各样地相同分辨率时DEM和DSM的误差。由表 3可知,地面点云在插值生成DEM时利用普通克里金插值的RMSE普遍小于用其他方法,而首次回波点云在插值生成DSM时利用反距离权重插值的RMSE普遍小于其他方法,说明在对比的5种插值方法中,普通克里金法适合于插值DEM,而反距离权重法适合插值DSM;总体而言,DEM的5种插值方法有效性由好到差依次为OK>NN>SP>IDW>LP,DSM的5种插值方法有效性由好到差依次为IDW>NN>OK>LP>SP,与顾春雷等(2011)Guo等(2010)的研究结果一致,因此本研究分别选用OK和IDW对DEM和DSM进行插值。

表 3 5种插值方法的有效性交叉验证结果 Tab.3 Cross-validation results of 5 interpolation methods

运用普通克里金法插值DEM时,山谷针阔混交林(A样地)的RMSE最小,为0.09 m,陡坡阔叶林(C样地)的RMSE最大,为0.20 m;坡度对插值精度有一定影响,坡度越小、坡位越平缓插值效果越好。而运用反距离权重法插值DSM时,山脊针阔混交林(B样地)的RMSE最小,为3.46 m,陡坡针叶林(E样地)的RMSE最大,为4.87 m,与DEM插值结果有些不同,最大值都为陡坡样地,最小值都为山脊样地。总体上DSM的RMSE比DEM的RMSE大,且DEM的RMSE随着样地地形坡度增加而增加,而DSM随着样地地形坡度、树种混交程度增加而增加。

3.2 林隙识别结果

对LiDAR提取出的林隙多边形矢量图层进行可视化检查可知,LiDAR对林隙的总识别率为92.6%(表 4),其中B、C、E样地的林隙识别率最高,A、D次之。山脊、山坡阔叶林和针叶林样地的林隙面积较大且分布较分散,所以识别率高。山谷样地和陡坡针阔混交林(A、D样地)中林层多为复层结构,较小林隙未能识别,造成4处漏判。

表 4 基于LiDAR的林隙识别 Tab.4 Gaps recognition using LiDAR

图 1A~F给出了运用交替贯序滤波法识别出各个样地的部分林隙和参考林隙结果,可以看出,基于机械激光雷达数据采用高度阈值及形态学滤波的交替贯序滤波方法能较好识别出各个样地中的林隙,较为规则的林隙识别匹配结果较好,而不太规则的林隙识别匹配结果较差。如山脊针阔混交林样地(图 1B)中1块面积较大的林隙形状非常不规则且周围林层结构复杂,未能很好地将林隙周围边界木的树冠区分开,以至于部分树冠被识别成空隙;而山坡阔叶林样地(图 1D)中大面积形状较规则的林隙识别匹配结果较好。

图 1 5块样地的林隙提取结果 Figure 1 Examples gap results of five gap investigate plots A.山谷针阔混交林A;B. 山脊针阔混交林B;C. 山坡阔叶林C;D.山坡阔叶林D;E.山坡针阔混交林E;F.山坡针叶林F。 A.Valley mixed forest A; B. Ridge mixed forest B; C. Hillside broad leaf forest C; D. Hillside broad leaf forest D; E. Hillside mixed forest E; F. Hillside needle leaf forest F.
3.3 林隙结构特征提取精度验证

表 5可知,样地A、D、E配对面积差值正态分布显著性均大于0.05,说明利用配对t检验具有统计学意义,样地B、C配对面积差值正态分布显著性小于0.05,不能采用配对t检验,所以利用符号秩检验比较。配对检验的零假设是面积差值的总体均值为零,除样地A外,其他样地P均大于0.05,样地总体P也大于0.05,说明LiDAR与野外观测的林隙面积没有显著差异,而样地A的P为0.044,说明样地A的林隙面积与LiDAR提取的林隙面积有一定差异。分别以LiDAR和野外观测的林隙面积为因变量和自变量进行线性回归分析,得到的相关系数(R2)为0.983,说明二者具有强线性相关性,回归方程斜率大于1,说明LiDAR估测的林隙面积通常比野外观测的大。

表 5 LiDAR与野外观测林隙面积的配对检验 Tab.5 LiDAR-field survey gap area canopy height paired samples test and regression

样地A利用LiDAR提取的林隙面积与野外观测面积有差异,但差异不大,可能是因为样地A较其他样地树种混交程度、林层结构复杂,林隙形状不规则。5块样地的平均相对误差在13%~20%之间,阔叶林样地C与针叶林样地E相比,针叶林样地的平均相对误差更小,这是因为阔叶林树冠形状比针叶林树冠更加不规则,所以在提取面积时平均相对误差较大;同为针阔混交林的样地A、B、D相比,样地D的平均相对误差最大,为20.74%在之前林隙识别,时样地D就有3个林隙未识别说明样地D林层结构和树种混交程度较复杂。

根据相关研究结果,本研究选择林隙边界木点云百分位高度H90反演边界木高(刘峰等,2014)。5块样地的边界木高差值正态分布显著性都大于0.05,所以采用配对t检验法检验,P均大于0.05,说明LiDAR估测的5块样地的林隙边界木高与野外观测高度没有显著差异;将所有边界木进行检验,二者配对的高度差值正态分布显著性小于0.05,所以采用符号秩检验法检验,林隙边界木P为0.105(>0.05),接受零假设,即二者测量的林隙边界木高没有显著差异。总体而言,LiDAR估测反演的值比野外观测的高,回归模型的R2为0.737,二者具有一定相关性(表 5)。根据平均相对误差可以看出,LiDAR估测林隙边界木高时与地形地势有关,样地A、B、D都为针阔混交林,其中山谷A样地相对误差最小为8.35%,山脊B样地相对误差最大为18.21%,中坡D样地介于二者之间,坡度越陡,对提取边界木高影响越大;对比同为陡坡的样地C、E可知,林分类型为阔叶的C样地小于针叶的E样地。

基于LiDAR数据在5块样地共识别出50个林隙,根据提取出的林隙面积和林隙周长可求出林隙形状,5块样地林隙形状指数均值为1.16,最大值为1.94,最小值为1.01,集中分布在1.06~1.24之间(表 6),即样地内多为较小形状指数、边缘效应不太显著的林隙。

表 6 林隙形状指数基本特征 Tab.6 Characteristic of gap shape index
3.4 林隙区域专题信息

基于LiDAR数据在17.36 hm2的帽儿山实验林场老山施业区部分林班中共识别出219个林隙,通过目视判读发现有22个林隙未能识别出,林隙识别率为90.87%,林隙密度为12.62个·hm-2,提取出的林隙面积在4.21~1 240.63 m2之间,均值为31.41 m2。将目视判读勾绘出的林隙面积作为真实值与提取的林隙面积进行配对检验(表 7),由于差值不服从正态分布所以采用非参数检验,P为0.893(>0.05),说明基于LiDAR数据提取此区域的林隙面积与目视判读勾绘出的林隙面积没有显著差异,相关系数R2为0.98,具有较强线性关系,平均相对误差为25.92%。

表 7 区域林隙配对检验 Tab.7 Area of gap paired samples test and regression

从林隙识别结果(图 2)来看,研究区内林隙形状较规则、面积较大的林隙能被很好识别且匹配程度较高,而形状不规则且林层结构复杂的林隙识别出后未能很好勾绘出边界形状,面积匹配程度较低。图中有22个林隙未能识别出,其原因有2个:1)林隙面积较小,形状不规则且大多处于林层结构较复杂的谷地;2)林隙形成木复杂而未能识别出。由此看出,本研究方法在区域尺度上不受地形地势的影响,识别效果较好,但林层结构、林隙形状和形成木类型对林隙识别提取有一定影响。

图 2 试验区林隙识别结果 Figure 2 The results of gaps in experimental area
4 讨论 4.1 影响插值精度的因素

在插值地面点云和首次回波点云时,本研究根据有效性分别选取普通克里金和反距离权重对DEM和DSM进行插值,应用于区域尺度林隙识别中效果较好,这说明普通克里金和反距离权重2种插值方法也适用于大尺度林隙识别。地形、坡度、坡位等对点云插值效果有一定影响,运用普通克里金法插值DEM时,坡度和坡位对精度的影响最为明显,针阔混交林A样地的RMSE最小,其坡位平缓且坡度最小,阔叶林C样地的RMSE最大,其坡位较陡且坡度较高;而运用反距离权重插值DSM时,针阔混交林B样地的RMSE最小,其坡度比山谷地高。造成这一现象的原因是2种方法所用插值数据不同,DSM是首次回波点云而DEM为地面点云,首次回波点云数据量远多于地面点云,林地冠层结构越复杂,插值精度越低。由于DEM插值时所用分类后的地面点数据,所以地形、坡度对其插值精度影响较大,而DSM插值时除了地形、坡度对其造成影响外,林地内冠层结构的复杂程度也是影响因素之一。

4.2 影响林隙识别的因素

林隙识别率和匹配程度是林隙结构参数正确提取与否的关键。本研究共调查54块林隙,其中4块未识别(样地A中1块,样地D中3块),林隙识别率为92.6%;样地A和D都是针阔混交林,林层结构复杂,面积较小的林隙不容易被识别。依据刘峰等(2015)方法对提取出的林隙多边形进行验证,林隙多边形与相邻树冠高度差在5 m以上,此外本研究野外调查的林隙形成木都为折干或砍伐后的伐桩,所以将转换的阈值设置为5 m,通过LiDAR数据提取出的林隙与在高分辨率遥感图像勾绘的参考林隙(图 1)进行对比,发现利用高度阈值与形态学滤波方法识别提取林隙不受地形地势的影响,适合于森林中任何地形条件。这是因为在进行交替贯序滤波时已将图像设置了高度阈值,转换成了二值图像,所以高度信息、地形地势不会影响林隙的识别提取;而影响林隙提取的因素主要是林隙的边界形状和林层结构,这2个因素对林隙识别影响不是很大,主要影响林隙面积提取的正确性,林隙形状越不规则、森林结构越复杂,林隙的提取难度越大。

4.3 影响林隙结构参数提取的因素

基于LiDAR点云数据采用交替贯序滤波方法提取出林隙结构参数并验证,此方法提取林隙面积时几乎不受地形、坡度影响,但是受林层结构、树种混交程度以及林隙形成木的影响较大,本研究样地的林隙形成木类型都为折干或砍伐后的伐桩,所以采用高度阈值可以有效筛选出林隙,而林隙边界木高的提取受地形、坡度的影响较大,LiDAR估测反演值比野外观测值高与刘峰等(2015)结果相反,可能的原因有2点:1) LiDAR点云密度不同,刘峰等(2015)点云平均密度为4~6 pts·m-2,而本研究考虑多回波和重复覆盖,最大点云密度超过10 pts·m-2;2)野外测量树高时,由于地形因素和人工操作可能存在一定误差。阔叶林样地平均相对误差小于针叶林样地与刘峰等(2015)研究结果一致,因为针叶树和阔叶树冠形的垂直结构变化差异,针叶树树冠自上而下形状变化较快,可能导致LiDAR点云采样时未能碰到树冠顶点。交替贯序滤波方法在提取林隙面积时需要对其迭代次数和结构元素大小进行设置,结构元素尺寸过大将会滤去面积较小的林隙,所以采用交替贯序滤波时要根据研究区域林隙大小来调整参数,做到因地制宜。

5 结论

本研究以东北林业大学帽儿山实验林场内设置的5块林隙调查样地中的54个林隙为研究对象,探讨了基于LiDAR数据插值出冠层高度模型,采用形态学交替贯序滤波识别并提取林隙结构参数的过程。生成冠层高度模型时,普通克里金法适合插值DEM,反距离权重法适合插值DSM;总体上DSM的RMSE比DEM的RMSE大,且DEM的RMSE随着样地地形、坡度增加而增加,而DSM随着地形、坡度、树种混交程度增加而增加。采用交替贯序滤波法对调查的林隙进行识别,其识别率为92.6%,效果较好,且提取的林隙面积和边界木高与野外测量并无显著差异。

在生成冠层高度模型时,坡度、坡位及冠层结构复杂程度都会影响插值精度。使用高度阈值转换为二值图像进行交替贯序滤波识别并提取林隙结构参数时,林隙面积提取不受地形因素的影响,但是林隙边界形状和林层结构对其有重要影响,林隙边界木高随着地形、坡度增加误差也随之增大。本研究识别提取林隙结构参数采用高度阈值和交替贯序滤波法,应进一步研究此方法是否适用于识别提取形成木类型更加复杂的林隙结构参数。机载LiDAR数据具有在多时空尺度上获取森林生态系统高精度的植被结构信息和三维地形特征的能力,是高效描述林隙特征的有效工具,本研究不仅为林隙遥感识别提供了参考,而且为野外林隙结构参数测量给予了技术支持。

参考文献(References)
费世民, 何亚平, 何飞, 等. 2009. 关于森林林窗中几个问题的综述和展望. 四川林业科技, 30(4): 28-37.
(Fei S M, He Y P, He F, et al. 2009. Review and prospects on some questions of forest gaps. Journal of Sichuan Forestry Science & Technology, 30(4): 28-37. DOI:10.3969/j.issn.1003-5508.2009.04.005 [in Chinese])
顾春雷, 杨漾, 朱志春. 2011. 几种建立DEM模型插值方法精度的交叉验证. 测绘与空间地理信息, 34(5): 100-102.
(Gu C L, Yang Y, Zhu Z C. 2011. Accuracy cross-validation of several interpolation methods of DEM. Geomatics & Spatial Ineformation Technology, 34(5): 100-102. [in Chinese])
刘峰, 谭畅, 王红, 等. 2015. 基于机载激光雷达的中亚热带常绿阔叶林林窗特征. 应用生态学报, 26(12): 3611-3618.
(Liu F, Tan C, Wang H, et al. 2015. Characterization of mid-subtropical evergreen broad-leaved forest gap based on light detection and ranging(LiDAR). Chinese Journal of Applied Ecology, 26(12): 3611-3618. [in Chinese])
刘峰, 杨志高, 龚健雅. 2011. 基于标记控制分水岭算法的树冠高程模型分割. 中国农学通报, 27(19): 49-54.
(Liu F, Yang Z G, Gong J Y. 2011. Canopy height model segmentation based on marker controlled watershed algorithm. Chinese Agricultural Science Bulletin, 27(19): 49-54. [in Chinese])
刘峰, 谭畅, 雷丕锋. 2014. 中亚热带森林单木地上生物量的机载激光雷达估测. 应用生态学报, 25(11): 3229-3236.
(Liu F, Tan C, Lei P F. 2014. Estimating individual tree aboveground biomass of the mid-subtropical forest using airborne LiDAR technology. Chinese Journal of Applied Ecology, 25(11): 3229-3236. [in Chinese])
梁欣廉, 张继贤, 李海涛, 等. 2005. 激光雷达数据特点. 遥感信息, (3): 71-76.
(Lang X L, Zhang J H, Li H T, et al. 2005. LiDAR data characteristics. Remote Sensing Information, (3): 71-76. DOI:10.3969/j.issn.1000-3177.2005.03.019 [in Chinese])
吴丹宇, 肖新平. 2008. Wilcoxon符号秩检验的P-值的数值解法. 中国青年信息与管理学者大会.
(Wu D Y, Xiao X P. 2008. Numerical solution of P-value of Wilcoxon's symbol rank test. China Youth Information and Management Scholars Conference. [in Chinese])
朱教君, 刘世荣. 2007. 森林干扰生态研究. 北京: 中国林业出版社.
(Zhu J J, Liu S S. 2007. Studies on disturbance ecology of forest. Beijing: Chinese Forestry Publishing House. [in Chinese])
郑斌, 周继存. 2010.基于TerraScan的机载雷达激光数据分类.地理信息与物联网论坛暨江苏省测绘学会2010年学术年会论文集.
(Zheng B, Zhou J Z. 2010. Airborne LiDAR data classification based on TerraScan. Geographic Information and Internet of Things Forum and Jiangsu Institute of Surveying and Mapping 2010 annual meeting of the Proceedings. [in Chinese])
Asner G P, Keller M, Rodrigo Pereira J, et al. 2004. Canopy damage and recovery after selective logging in Amazonia: field and satellite studies. Ecological Applications, 14(SP4): S280-S298. DOI:10.1890/01-6019
Gaulton R, Malthus T J. 2010. LiDAR mapping of canopy gaps in continuous cover forests: a comparison of canopy height model and point cloud based techniques. International Journal of Remote Sensing, 31(5): 1193-1211. DOI:10.1080/01431160903380565
Guo Q, Li W, Yu H. 2010. Effects of topographic variability and LiDAR sampling density on several DEM interpolation methods. Photogrammetric Engineering & Remote Sensing, 76(6): 701-712.
ÖhmanK, LämåsT. 2005. Reducing forest fragmentation in long-term forest planning by using the shape index. Forest Ecology and Management, 212(1): 346-357.
Runkle J R. 1981. Gap regeneration in some old-growth forests of the eastern United States. Ecology, 62: 1041-1051. DOI:10.2307/1937003
Runkle J R. 1982. Gap patterns of disturbance in some old-growth mesic forests of eastern North America. Ecology, 63: 1533-1546. DOI:10.2307/1938878
Tyrrel L E, Crow T R. 1994. Structural characteristics of old-growth hemlock-hardwood forests in relation to age. Ecology, 75: 370-386. DOI:10.2307/1939541
Tanaka H, Nakashizuka T. 1997. Fifteen year of canopy dynamics analyzed by aerial photographs in a temperate deciduous forest. Journal of Ecology, 58: 835-841.
Watt A S. 1947. Pattern and process in the plant community. Journal of ecology, 35(1/2): 1-22. DOI:10.2307/2256497
Worral J J, Harrington T C. 1988. Etiology of canopy gaps in spruce-fir forests of Crawford Notch, New Hampshire. Can J For Res, 18: 1463-1469. DOI:10.1139/x88-225
Zhang K. 2008. Identification of gaps in mangrove forests with airborne LiDAR. Remote Sensing of Environment, 112(5): 2309-2325. DOI:10.1016/j.rse.2007.10.003