基于ICESat-GLAS数据和回波仿真原理识别森林类型

蔡龙涛，邢艳秋，黄佳鹏，崔 阳，秦 磊，马建明，赵霄洋

（1.东北林业大学 工程技术学院，黑龙江 哈尔滨 150040；2.杭州极遥科技有限公司，浙江 杭州 311200）

森林作为陆地生态系统重要的组成部分，存储了陆地上近45%的碳，并贡献了约50%的陆地净初级生产[1]。森林生物量的估测研究，对陆地生态系统碳储量及碳循环研究提供了重要的理论支撑。森林生物量估测时受森林水平结构参数和森林垂直结构参数影响较大。而森林类型作为森林水平结构参数的重要组成部分，对森林类型的精确识别在森林生物量定量化估测方面具有重要的理论和现实意义。传统野外实地调查耗时久，费用高，难以在短期内实现大区域森林类型识别。因此，为实现大区域森林类型快速、高精度识别，就需要引入一种全新的森林类型识别技术与方法。

随着对地观测技术的快速发展，遥感技术特别是卫星遥感技术逐渐成为获取地面信息的重要手段和途径，且因激光测高在获取地表三维信息上的快速和高精度特点，使得卫星遥感技术逐步成为新型地理空间信息获取手段[2]。2003年1月冰、云和陆地高程卫星/地球科学激光测高系统（the Ice,Cloud,and land Elevation satellite/Geoscience Laser Altimeter System,ICESat/GLAS）发射升空，2009年10月结束数据采集，期间获取了大量对地观测数据，在树高估测[3-7]和生物量反演[9-11]方面，都取得了不错的应用效果。但在森林类型识别方面，受光斑直径影响，GLAS 波形数据作为大尺度测量级数据[12]，回波波形中携带了大量噪声信息，当携带噪声信息的波形数据应用于森林类型识别时，难以实现森林类型的有效识别。

为提高森林类型识别精度，国内外学者在森林类型识别方面做了大量研究。Ranson 等[13]对地形坡度在5°以下的林地进行了研究，研究过程中把林地分为4 种类型：阔叶林、常绿针叶林、落叶针叶林和混交林，以此分析不同森林类型与GLAS 回波波形特征参数之间的对应关系。Li 等[14]引入了支持向量机（Support vector machine，SVM）分类法，对针叶林、阔叶林和混交林3 种森林类型进行识别，总体识别精度为53.66%。而刘美爽等[15]剔除了混交林这一森林类型，仅对针叶林和阔叶林两组森林类型进行识别研究，总体识别精度达到85.24%；Zhang 等[16]依据地形坡度把林地分为两类，仅对地形坡度在14.2°以下的GLAS 回波波形进行数据处理，识别森林类型，其中，针叶林、阔叶林和混交林森林类型总体识别精度为75.56%，针叶林和阔叶林森林类型总体识别精度为90.57%。综上所述，目前提高森林类型识别精度的方法多为引入新型分类法、减小森林类型个数和剔除高坡度地形法，然而这些方法并未从回波波形获取原理方面分析GLAS 发射波波形与森林类型作用方式，研究不同森林类型与GLAS 波形特征参数之间的对应关系，进而提取与森林类型相关的波形特征参数，用于森林类型识别。

依据以上研究的不足之处，本研究结合回波仿真原理分析GLAS 发射波波形和森林植被作用方式，对GLAS 回波波形进行了理论分析，并据此提取与森林类型相关的波形特征参数R-cafit1-47和K1-47，然后与其他波形特征参数建立参数组合。再利用建立的波形特征参数组合分别对针叶林和阔叶林，以及针叶林、阔叶林和混交林进行森林类型识别，并提高了森林类型识别精度。此外，本研究基于GLAS 波形数据提出的波形特征参数为森林类型识别提供了参数支持；结合回波仿真原理对GLAS 回波波形获取方式进行理论分析，对后续波形研究具有指导意义。

1 研究区概况和研究方法

1.1 研究区域概况

研究区位于吉林省汪清林业局管辖区，该区域属于长白山系中低山区（129°56′～131°04′E，43°05′～43°40′N），地处寒温带，总面积3.04×105hm2，南北长约60 km，东西长约85 km，地面高程变化范围为360～1 477 m，地形坡度变化范围为0～45°（图1）。

图1 研究区位置及野外样地分布（圆形代表采样点）Fig.1 Location of study area and the distribution of field sampling plots(Cycles represents field sampling plots)

林区内森林覆盖率达到95.95%，深山区林分以针叶林、阔叶林和混交林为主，带状分布于海拔500～1 100 m 之间。针叶树主要有红松Pinus koraiensis、云杉Peiceaspp.和臭松Symplocarpus foetidus、阔叶松Larix gmelini，阔叶树多为椴树Tiliaspp.、蒙古栎Quercus momglicaFisch.、枫桦Betula costaTrautv.、色木Acer monoMaxim 和白桦Betula platyphylla等。

1.2 研究数据

1.2.1 GLAS 波形数据

ICESat-1/GLAS 是第一个极地轨道大光斑激光雷达卫星，该卫星共提供15 个数据产品：GLA01、GLA02、…、GLA15。其中，GLA01 数据产品记录了GLAS 回波波形数据；GLA14 数据产品记录了GLAS 波形数据对应的地面光斑地理位置和高程数据。GLAS 波形数据可从美国国家冰雪数据中心（http://nsidc.org/data/ice-sat/）下载，美国国家冰雪数据中心拥有ICESat 卫星2003年至2009年采集的有效波形数据。

本研究所用波形数据为2006年10月、2009年10月和2010年10月的GLA01 和GLA14 波形数据产品，通过经纬度查找样地所对应的GLAS回波波形。

1.2.2 野外调查数据

所用数据分别于2006年9月、2007年9月和2010年9月在长白山汪清林业局研究区采集获取，共采集286 个GLAS 光斑点。野外采样时利用分层随机采样法，对具有GLAS 激光光斑覆盖的样地进行数据采集。在实地勘测过程中，以针叶林、阔叶林和针阔混交林为研究对象，使用GPS 对已选定的激光光斑采样点进行定位，把GLAS 激光光斑中心点作为地面调查样地的圆心，建立水平投影面积为500 m2的圆形样地，记录样地内植被分布情况、植被类型和植被覆盖度。结合我国森林资源调查主要技术规定将森林类型划分为针叶林、阔叶林和混交林，把针叶林蓄积量比例占总蓄积量65%以上的林分定义为针叶林，阔叶林蓄积量比例占总蓄积量的65%以上的林分定义为阔叶林，任何一个树种蓄积量占总蓄积量不到65%的林分定义为混交林。

1.3 结合回波仿真原理选取波形特征参数

不同森林类型林分对应的GLAS 回波波形受GLAS 发射波（指发射波波形和LPA 分布）及光斑内林分冠层结构影响，呈现出不同的波形起伏。如阔叶树叶片平展且阔大，冠形多为椭球体、上半球体和下半球体形；针叶树叶片如针，呈扁平尖细状，冠形多为圆锥体形。由阔叶树和针叶树组成的针叶林、阔叶林和混交林与发射波作用后在GLAS 回波波形上产生的波形起伏有所不同。本小节结合回波仿真原理[17]及不同森林类型冠层结构的差异性分析GLAS 发射波波形与不同森林类型的作用方式，探究不同森林类型对GLAS 回波波形的影响，并从中选取与森林类型相关的波形特征参数，用于森林类型识别。

回波仿真原理实际为发射波函数与地表响应函数的卷积，即GLAS 回波波形为GLAS 发射波波形函数与地表响应函数的卷积。本研究依据卷积原理把GLAS 回波仿真模型函数细分为3 个函数，如式（1）所示。3 个函数分别为有效回波波形前47（发射波波形长度48-1）帧回波波形能量值函数：随波形帧数增加，对应的回波能量值方程自变量个数逐渐增加；有效回波波形中第47～N帧回波波形能量值函数：随帧数的增加，所对应的回波能量值方程自变量个数不变；有效回波波形中第N～N+47 帧回波波形能量值函数：随帧数的增加，所对应的回波能量值方程自变量个数逐渐减小。其中，有效回波波形指携带地物信息的回波波形。

式（1）中，wave 为有效回波波形；N为以0.15 m为间隔对光斑内地物信息进行分层，地物信息最高点到最低点之间的层数；m为有效回波波形前47 帧中第m个回波波形能量值对应的帧数；n为有效回波波形中第n个回波波形能量值对应的帧数，其中，N＞n＞47；p为有效回波波形中第p个回波波形能量值对应的帧数，其中，p＞N；f(i)为GLAS 测高系统发射波波形中第i帧发射波波形能量值；按照网格划分原理，对光斑内地物以0.15 m 为间隔逐层分割，共分为N+47 层；对每层地物覆盖面面积进行统计，第i层地物覆盖面面积为area(i)。

另外，不同森林类型冠层水平分布有所不同，如阔叶林林分树叶多集中于冠层上半部分，而针叶林林分树叶多集中于冠层下半部分。当GLAS发射波波形到达林分冠层上半部分（有效回波波形前47 帧）之后，受回波波形能量方程影响，使得阔叶林林分冠层波峰斜率小于针叶林；当GLAS发射波波形到达林分冠层下半部分之后，受冠层反射和透射作用影响，发射波波形需穿过阔叶林林分上半部分冠层后方可获取阔叶林下半部分冠层反射回波，导致阔叶林林分冠层对应波形斜率值小于针叶林。可利用不同森林类型冠层水平分布和垂直分布不同引起的波形差异进而提取波形特征参数，识别森林类型。

本研究选用的研究区为天然林区，林区内存在幼龄林，而幼龄林林木树高值较低，使得光斑内存在最大树高低于7.05 m 的情况，反映在GLAS 回波波形上为波形起波点到地面回波波形起波点之间的波形长度小于47 帧。针对这种情况，本研究利用波形特征参数Rfit1代替R-cafit1-47，波形特征参数K1代替K1-47，作为森林类型参数组合因子，用于识别森林类型。其中，Rfit1为GLAS 有效回波波形中第一个高斯拟合波形能量值之和与冠层回波波形能量值之和的比值；R-cafit1-47为GLAS 有效回波波形中前47 帧回波波形能量值之和与冠层回波波形能量值之和的比值；K1为GLAS 回波波形起波点与第一个高斯拟合波形波峰点之间斜率值；K1-47为GLAS有效回波波形中前47帧波形斜率值。

1.4 波形特征参数的提取

GLAS 回波波形中含有噪声数据，波形特征参数提取之前需进行去噪处理，以前学者多采用小波变换[18-20]或高斯算法[21-24]进行去噪。考虑到森林冠层垂直分布差异，本研究选用高斯拟合方程对GLAS 回波波形进行去噪处理，然后基于拟合后的波形进行波形特征参数提取。

1.4.1 地面回波波形的提取

GLAS 光斑内地形坡度较大时能够导致回波波形出现波形展宽[25-26]现象，影响冠层回波能量值在总回波能量值中的比重，增大森林类型分类难度。因此，在提取波形特征参数之前，需剔除地面回波波形。在地形相对平坦的地区，Allouis 等[27]把高斯拟合波形中最后一个高斯波形作为地面回波波形。在地形坡度小于14.04°的光斑中，Zhang 等[16]把有效回波波形中右半部分回波波形能量值最大处定义为地面，所处的拟合波形定义为地面回波波形。

本研究在进行森林类型分类时，并未剔除高坡度GLAS 回波波形。另外，研究时把有效回波波形中右半部分回波波形能量值最大处定义为地面回波波形，所处的拟合波形定义为第一个地面回波波形。又因光斑内地形呈无规律起伏，反映在GLAS 回波波形上为地面回波波峰之后可能存在第二个或第三个地面回波波形，本研究对所有地面回波波形进行叠加，组成新的地面回波波形，并把该地面回波波形作为GLAS 光斑内实际地面回波波形。

1.4.2 参数的提取

GLAS 回波波形高斯拟合之后，依据拟合波形提取波形特征参数：R-cafit1-47和K1-47，公式分别如式（2）和式（5）所示：

式（2）～（5）中，R-cafit1-47为有效波形中前47帧波形能量值之和；Efit为拟合波形能量值之和，Efit(i)有效回波波形中，第i帧回波波形能量值；EG为有效回波波形中地面回波波形能量值之和；ground_sig 为拟合波形中地面回波波形起始点处帧数值；K1-47为有效回波波形中前47 帧波形斜率值；wave 为有效回波波形。

当GLAS 回波波形起波点到地面回波波形起波点之间波形长度小于47 帧时，本研究利用Rfit1代替R-cafit1-47，K1代替K1-47，作为波形特征参数值，用于森林类型识别，Efit1和K1计算公式分别如式（7）和式（8）所示。

式（6）～（9）中，Efi1t为GLAS 回波波形中第一个高斯拟合波形能量值之和；Efit为GLAS 回波波形能量值之和；m为第一个高斯拟合波形长度值；K1为第一个高斯拟合波形前倾角斜率值；A1为第一个高斯拟合波形振幅值；σ1为第一个高斯拟合波形标准差；Efit(i)为GLAS 有效回波波形中，第i帧回波波形能量值。

另外，本研究还提取了Zhang 等（2011）[16]提出的波形特征参数AGS、SGS 和MSGS，用于森林类型识别。其中，AGS、SGS 和MSGS 计算方法分别如式（10）、（11）和（13）所示。

式（10）～（12）中，n为GLAS 冠层回波波形高斯拟合个数；Ai为冠层回波波形中第i阶高斯拟合波形振幅值；σ1为第i阶高斯拟合波形标准差。

式（13）～（14）中，Si为GLAS 冠层回波波形中第i阶高斯拟合波形斜率值；Ei为GLAS 冠层回波波形中第i阶高斯拟合波形能量值；ET为GLAS冠层回波波形能量值之和；n为GLAS 冠层回波波形高斯拟合阶数。

1.5 波形特征参数组合

结合回波波形仿真原理，对GLAS 回波波形获取方式进行理论分析，提取波形特征参数：R-cafit1-47和K1-47。然后与以前学者提出的波形特征参数AGS、SGS、MSGS 联合，建立波形特征参数组合。并引入支持向量机分类算法对不同森林类型进行识别研究。其中，波形特征参数组合如表1所示。

表1 支持向量机分类中需要测试的波形特征参数组合Table 1 Waveform feature parameter combination to be tested in support vector machine classification

2 结果与分析

2006年、2007年和2010年3年野外实测调查，共收集了286 组样地数据。其中阔叶林样地137 组，针叶林样地50 组，混交林样为60 组，其他类型样地如水地、裸地和草地等共计39 组。本研究对针叶林、阔叶林和混交林进行森林类型识别时，从阔叶林随机选取96 组样地数据作为训练集数据，41 组数据作为测试集数据；从针叶林随机选取35组样地数据作为训练集数据，15 组数据作为测试集数据；从混交林随机选取42 组样地作为训练集数据，18 组数据作为测试集数据。

2.1 针叶林和阔叶林森林类型识别结果与分析

把结合回波仿真原理提取的波形特征参数：R-cafit1-47和K1-47，与其他森林类型相关波形特征参数：AGS、SGS 和MSGS 进行联合，建立波形特征参数组合，利用支持向量机分类算法对针叶林和阔叶林进行识别。其中，表2为不同波形特征参数组合条件下针叶林和阔叶林森林类型识别结果；表3为波形特征参数组合R-cafit1-47、K1-47对针叶林和阔叶林森林类型识别结果；表4为波形特征参数组合R-cafit1-47、K1-47、AGS 和MSGS对针叶林和阔叶林森林类型识别结果。

表2 不同参数组合针叶林和阔叶林森林类型分类结果Table 2 Classification results of forest types in broad-leaved forests and coniferous forests with different parameters

表3 基于R-cafit1-47、K1-47 的森林类型分类结果Table 3 Classification results of forest types based on R-cafit1-47and K1-47

表4 基于R-cafit1-47、K1-47、AGS 和MSGS 的森林类型分类结果Table 4 Classification results of forest types based on R-cafit1-47,K1-47,AGS and MSGS

从表2可以看出，在对针叶林和阔叶林进行森林类型识别时，结合回波仿真原理提取的波形特征参数组合R-cafit1-47、K1-47森林类型识别精度为92.86%，略高于波形特征参数AGS、SGS 组合森林类型识别精度85.71%，以及AGS、MSGS 组合森林类型识别精度85.71%。分析其原因为提取的波形特征参数R-cafit1-47、K1-47在林分的林层结构中对应于林分冠层顶部位置，而针叶林和阔叶林在林分冠层顶部位置结构特征差异较为明显，不同的结构特征反映在GLAS 回波波形上为波形特征呈现较大差异，利用针叶林和阔叶林不同的结构特征在GLAS 回波波形上的差异性提取的波形特征参数，有利于森林类型识别；而波形特征参数AGS、SGS 组合以及AGS、MSGS 组合在林分的林层结构中对应于整个冠层结构，如GLAS回波波形中波峰个数与林分中林层个数具有较大相关性，各个波峰能量值也与其对应林层覆盖面面积存在较高相关性，虽然针叶林和阔叶林在林层结构上存在差异，然而在林木密度较大情况下，不同林木枝叶间的交叉分布使得针叶林和阔叶林林层结构差异不明显，反映在GLAS 回波波形上为针叶林和阔叶林波形特征存在较大相似性，利用波形特征存在较大相似性的回波波形进行波形特征参数提取并进行森林类型识别，其总体分类精度会相对较低。

表3详细展示了波形特征参数组合R-cafit1-47、K1-47对针叶林和阔叶林森林类型识别结果，其中，阔叶林森林类型识别精度为95.12%，略高于针叶林森林类型识别精度86.67%。分析其原因为相对于针叶树，阔叶树冠层覆盖面较大，反映在GLAS回波波形上其波形特征更为明显，利用GLAS 回波波形提取波形特征参数，更有利于阔叶林的识别；另外，按照森林类型定义标准，针叶林中同样存在阔叶树，当针叶林中存在多颗阔叶树时，针叶林对应的GLAS 回波波形特征参数与阔叶林更为相似，容易出现针叶林错分为阔叶林的现象，使得针叶林林分识别精度较低。

为提高针叶林和阔叶林森林类型识别精度，本研究把波形特征参数R-cafit1-47、K1-47与波形特征参数AGS、SGS 和MSGS 进行联合，重新建立波形特征参数组合R-cafit1-47、K1-47、AGS、MSGS与R-cafit1-47、K1-47、AGS、SGS，然后进行森林类型识别。其中，两种参数组合森林类型识别精度都 为94.64%，而R-cafit1-47、K1-47、AGS、MSGS参数组合森林类型识别Kappa 系数值为0.871 6，略高于R-cafit1-47、K1-47、AGS、SGS 参数组合。分析其原因为相对于波形特征参数SGS，波形特征参数MSGS加入了回波波形中不同波峰能量值，增大了波形特征参数与森林类型相关性，当波形特征参数MSGS 用于森林类型识别时，森林类型识别结果与真实结果更为相近。

2.2 阔叶林、针叶林和混交林森林类型识别结果与分析

本研究还利用波形特征参数组合对针叶林、阔叶林和混交林进行了森林类型识别。其中，表5为不同波形特征参数组合条件下针叶林、阔叶林和混交林森林类型识别结果；表6为波形特征参数组合R-cafit1-47和K1-47针叶林、阔叶林和混交林森林类型识别结果；表7为波形特征参数组合R-cafit1-47、K1-47、AGS 和MSGS 针叶林、阔叶林和混交林森林类型识别结果。

表5 不同参数组合针叶林、阔叶林和混交林森林类型分类结果Table 5 Classification results of forest types in broad-leaved forest mixed forest and coniferous forests with different parameters

从表5可以看出，结合回波仿真原理提取的波形特征参数R-cafit1-47、K1-47与其他波形特征参数AGS、SGS 和MSGS 建立的参数组合用于森林类型识别时，识别精度较高，其中R-cafit1-47、K1-47、AGS 和MSGS 参数组合森林类型识别精度最高，其值为89.19%。而波形特征参数R-cafit1-47、K1-47组合森林类型识别精度低于波形特征参数AGS、MSGS 和AGS、SGS 组合，分析其原因为波形特征参数R-cafit1-47、K1-47的值为林分中冠层最高处位置对应的冠形特征参数值，即林分中优势树种冠形特征参数值，而混交林林分中优势树种不唯一，可能是阔叶树，也可能是针叶树，这种树种的不唯一性使得波形特征参数R-cafit1-47、K1-47组合进行森林类型识别时，容易把混交林误分为阔叶林或针叶林，即表6中混交林森林类型识别结果，使得森林类型识别精度较低。

表6 基于R-cafit1-47和K1-47 的森林类型分类结果Table 6 Classification results of forest types based on R-cafit1-47and K1-47

对比表6和表7可以看出，利用波形特征参数R-cafit1-47、K1-47组合进行森林类型识别时混交林森林类型识别精度较低，而利用波形特征参数R-cafit1-47、K1-47、AGS、MSGS 组合进行森林类型识别时混交林森林类型识别精度较高。其原因是波形特征参数AGS、MSGS 在林分的林层结构中对应于整个冠层结构，在波形特征参数组合R-cafit1-47、K1-47中加入波形特征参数AGS、MSGS，能够弥补形特征参数组合R-cafit1-47、K1-47在整体冠层方面分析的不足，提高混交林森林类型识别精度。此外，波形特征参数R-cafit1-47、K1-47组合对针叶林进行识别时其识别精度低于波形特征参数R-cafit1-47、K1-47、AGS、MSGS 组合。其原因是波形特征参数AGS、MSGS，能够弥补波形特征参数R-cafit1-47、K1-47组合在冠层方面分析的不足，然而在林木密度较大情况下枝叶的交叉分布使得3 种森林类型之间林层结构较为相似，即不同森林类型波形特征参数AGS、MSGS 的值较为相近，这种情况增大了针叶林森林类型识别难度，导致针叶林森林类型识别精度有所下降。

表7 基于R-cafit1-47、K1-47、AGS 和MSGS 的森林类型分类结果Table 7 Classification results of forest types based on R-cafit1-47,K1-47,AGSand MSGS

3 结论与讨论

以吉林省汪清林业局管辖区为研究区，利用ICESat-GLAS 波形数据，结合回波仿真原理，提取了波形特征参数R-cafit1-47、K1-47，并与其他同森林类型相关的波形特征参数AGS、SGS、MSGS进行联合，建立波形特征参数组合，通过支持向量机分类法对任意坡度下针叶林和阔叶林，以及针叶林、阔叶林和混交林进行森林类型识别。研究结论如下：

1）结合回波仿真原理提取的波形特征参数组合R-cafit1-47、K1-47对针叶林和阔叶林进行森林类型识别，森林类型识别精度优于波形特征参数组合AGS、SGS 和AGS、MSGS；说明波形特征参数组合R-cafit1-47、K1-47在针叶林和阔叶林森林类型识别方面，有较高优势。

2）利用波形特征参数组合R-cafit1-47、K1-47对针叶林、阔叶林和混交林进行森林类型识别，森林类型识别精度低于波形特征参数组合AGS、SGS 和AGS、MSGS。说明波形特征参数组合R-cafit1-47、K1-47在针叶林、阔叶林和混交林森林类型识别方面能力相对较弱。

3）联合波形特征参数R-cafit1-47、K1-47与其他波形特征参数AGS、SGS 和AGS、MSGS，建立新的波形特征参数组合，相对于波形特征参数组合R-cafit1-47、K1-47，AGS、SGS 和AGS、MSGS，较大地提高了针叶林和阔叶林，以及针叶林、阔叶林和混交林森林类型识别精度。

结合回波仿真原理与森林冠形特征对GLAS 回波波形进行理论分析，提取与森林类型相关的波形特征参数，并与Zhang 等[16]提出的波形特征参数建立波形特征参数组合，从理论上弥补了本研究波形特征参数在森林类型识别上的不足，而波形特征参数组合森林类型识别精度值验证了本研究建立的波形特征参数组合在森林类型识别上的优势。利用建立的波形特征参数组合对针叶林和阔叶林进行森林类型识别，较本研究提出的波形特征参数分类精度提高了1.78%，较Zhang 等[16]提出的波形特征参数分类精度提高了8.91%；波形特征参数组合在针叶林、阔叶林和混交林森林类型识别方面较本研究提出的波形特征参数分类精度提高了12.16%，较Zhang 等[16]提出的波形特征参数分类精度提高了8.11%，有效提高了森林类型分类精度。而加入混交林之后，虽然森林类型总体分类精度有所提高，但是精度值仍然较低。为进一步提高森林类型识别精度，后续将采用模糊模式识别分类法对森林类型进行识别研究，以期提高森林类型分类精度。