基于高光谱影像的舟曲曲瓦沟树种识别

周利鹏　马金辉

摘要 正确识别舟曲地区沟内森林树种是进行植被保护的基础，也为准确计算泥石流发生的降雨阈值提供参考。利用HYPERION高光谱影像，采用基于纯净像元指数（PPI）的端元提取方法，提取了6类端元，用波谱特征拟合方法（SFF）、波谱角分类（SAM）方法和二进制编码方法（BE）识别出核桃、矮灌木、橡树、刺柏、冷杉和灌木蒿。同时利用基于最小噪声分离（MNF）的最小能量约束法（CEM）分离出各自的分布范围。

关键词 高光谱；树种识别；Hyperion；波谱分析；混合像元分解

中图分类号 S126 文献标识码 A 文章编号 0517-6611（2014）16-05298-04

森林作为陆地生态系统的主体，在治水保土中发挥着重要作用。森林既能减少林地土壤或森林流域的水分入收量，也能削减或削弱降雨动能。不同的森林覆盖类型对降雨水分的截留量和对雨滴动能的削弱效果是不相同的[1]。甘肃白龙江泥石流地质灾害频发，泥石流高度发育，对当地人民的生命和财产安全造成极大威胁，因此，泥石流监测预警工作十分重要，降雨作为泥石流激发的决定性因素也被用于泥石流预测工作中[2]。因此正确识别流域内森林物种不仅是当地植被保护的基础和依据[3]，也是计算泥石流爆发降雨阈值的重要因素。传统的野外调查方法费时费力、成本较高，并且无法涉及到地形中难于到达的区域。采用大比例航片进行判读也存在费时费力、成本昂贵的问题。随着遥感技术的发展，遥感方法被应用于树种识别，但受限于遥感数据的低光谱分辨率和森林树种光谱特征的相似性，只能粗略地将其划分为植被、非植被，或者划分为针叶、阔叶[4-6]。高光谱影像具有较高的光谱分辨率，通常有数十到数百个波段，它能够反映地物间细微的光谱差异，在森林树种识别中有着广阔的应用前景[7]。Darvishsefat等利用HyMap高光谱数据和线性分离模型进行森林树种类型的识别研究，并与航空高光谱数据结果进行对比，证明了其在落叶松类复杂混交林的识别效果[8]。Aspinall利用波段组合、Logistics回归、建立光谱信息模型等方法进行森林主要树种类型的识别，均取得了与地面数据相吻合的结果[9]。宫鹏等利用高分辨率光谱仪（PSD1000）实地测定数据对美国加州乃至北美西海岸的6种主要树种进行识别，结果表明高光谱数据在树种识别上具有较强能力[10]。刘秀英等利用地物光谱仪采集数据，采用波段组合的方法，成功识别出杉木、雪松、小叶樟树和桂花这几个树种[11]。陈尔学等利用Hyperion高光谱数据，并以SPOT-5数据和地面观测数据为辅助，比较了几种高光谱统计模式识别方法，结果表明，对高光谱数据进行降维，并采用二阶统计量估计方法，进而结合空间上下文信息和光谱信息的分类算法（ECHO）可以有效地提高森林树种的识别精度[12]。王志辉等利用波段指数、分段主成分分析、光谱微分（一阶、二阶微分）和包络线去除等5种方法从叶片和图像角度验证了高光谱树种识别的可行性[13-14]。

Hyperion影像具有较高的光谱分辨率，能够反映地物目标的波谱特征。然而其空间分辨率只有30 m，因此在影像内混合像元普遍存在，同一像元的光谱特征往往受到几种地物的共同影响，因而加大了树种识别难度。为此，笔者利用HYPERION高光谱影像对舟曲西部曲瓦沟内的森林树种进行识别，探讨其在该地区森林树种识别上的可行性。

1 研究区概况及数据预处理

1.1 研究区概况 研究区域位于甘肃白龙江流域中上游，地处甘南藏族自治州舟曲县西部，与本州迭部县相接（图1）。该区域地处南秦岭山地，构造上属于西秦岭构造带的延伸部分，新构造运动活跃。气候属于温带半湿润气候，年降雨量800 mm左右，主要集中在5～10月份[15]，平均海拔2 886 m。

1.2 数据预处理 Hyperion是搭载在美国NASA地球观测卫星EO.1（Earth Observing.1）上的高光谱成像光谱仪。EO.1是NASA“新千年计划”的一部分，2000年11月21日在范登堡空军基地发射升空。EO.1卫星轨道与Landsat 7 基本相同，为太阳同步轨道，轨道高度为705 km，轨道倾角为98.2°，轨道周期为98.9 min，比Landsat 7 晚1 min过赤道。EO.1上搭载了3种传感器，除Hyperion外，还有高级陆地成像仪ALI、大气校正仪AC[16]。Hyperion是第1个星载民用的成像光谱仪，应用领域十分广泛。Hyperion传感器以推扫的方式获取可见光、近红外（VNIR，400～1 000 nm）和短波红外（SWIR，900～2 500 nm）光谱数据。Hyperion产品分两级：L0和L1，L0是原始数据，仅用来生成L1产品[17]。Hyperion L1产品只有一种数据格式HDF，波段存储格式为BIL。L1产品又分为L1A、L1B、L1R。2002年以前处理的数据为L1A，L1B产品由美国TRW处理而成，L1R产品由美国USGS处理生成。L1A与L1B、L1R最大的不同在于L1A产品没有纠正VNIR波段与SWIR波段之间的空间错位问题，而L1B和L1R产品中的VNIR波段与SWIR波段之间的空间位置是经过纠正的，用户无需再进行匹配[17-18]。该研究所用数据为L1R产品，该产品不能直接使用，必须进行数据预处理，流程如图2所示。

数据预处理流程 L1R产品有242个波段，波段范围356～2 577 nm，其中只有198个波段经过辐射定标处理，为8～57的可见光近红外波段（VNIR），77～224的短波红外波段（SWIR），由于VNIR 56～57与SWIR 77～78重叠，实际只有196个独立波段。在这些波段中，121～130、165～184和222～224波段受水气影响较大，极少包含地面信息，需要剔除。剩下的波段为8～57、79～120、131～164、185～221，共163个波段[17-19]。由于在Hyperion L1产品中存在坏线和条纹，坏线的DN值为零或者极小值，条纹的值不为零，一般表现为明显高于或低于周围列。对163个波段的影像逐一检查，并记录存在坏线和条纹的波段和列号。将坏线普遍存在的第1列和第256列删除，对其余坏线采用其相邻列的平均值代替达到修复的效果[18]，并采用“全局去条纹”的方法去除Hyperion影像上的条纹 [18，20]。Smile效应是指由于推扫式扫描传感器光谱扭曲产生的像元的波长从中心位置向两侧偏移[21]。Smile效应降低了影像的分类精度，必须加以纠正。该研究采用Goodenough等提出的Column Mean Adjusted in MNF Space方法[22]，分别对VNIR波段和SWIR波段进行MNF变换，取变换结果的前15个波段（后面的波段基本为噪声），经MNF逆变换，再将生成的VNIR和SWIR图像合在一起，将Smile效應纠正，并消除光谱曲线的锯齿现象[18-19]。该研究采用ENVI软件的FLAASH大气校正模块，由于Hyperion影像在生成时像元值采用了扩大因子，VNIR波段为40，SWIR波段为80，影像的单位为w/（m2·um·sr），而FLAASH大气校正模块需要输入的量纲为μm/（cm2·nm·sr），两者存在一个10倍的因子系数，所以在进行大气校正时应输入一个VNIR波段为400、SWIR波段为800的比例因子[17]。完成大气校正后图像的辐射值是乘以10 000的，采用Band Math工具将辐射值除以10 000转换为反射率。以校正好的TM影像为基准，采用人机交互选择地面控制点的方式对生成的影像进行几何校正。

2 树种识别

2.1 最小噪声分离（Minimum noise fraction，MNF） MNF是将一幅多波段图像的主要信息集中在前几个波段中，达到降维、分离数据中噪声的目的，以减少后处理中的计算量的方法。MNF是一种线性变换，本质是含有两次叠置分析的主成分分析。将包含163个波段的Hyperion影像全部进行MNF分析，可以看到，影像信息主要集中在前20個波段，后面的波段基本是噪声。

2.2 纯净像元指数（Pixel Purity Index ，PPI） ENVI的PPI工具通过N维散点图将所有像元迭代映射为一个随机单位向量，每次映射的极值单元被记录下来，将每个像元被标记为极值的总次数生成一幅“像元纯度图像”，像元的DN值表示像元被标记为极值的次数，迭代的次数越多，就越容易发现纯净像元，如图3所示。此次研究将迭代次数设置为10 000次，由于输入数据为MNF结果，所以阈值系数选择2和3即可，阈值系数越大，能找到更多的极值像元，但是包含的不纯净像元也越多。

图3 PPI图像2.3 N维散点图 利用n.D Visualizer工具将PPI值较大的像元点投影到MNF变换后的主成分空间，组成20维的空间散点图，设置适当的速度，通过交互式旋转，将处于点云拐点或末端的像元确定为端元。当数据维数高于4或5时，交互式识别和定义类别变得比较困难，可对类别进行Collapse变换降低维数。当端元识别完成以后，提取端元的波谱曲线并保存，图4为提取的端元波谱。

图4 端元波谱曲线2.4 波谱分析 对ENVI标准波谱库中的USGS植被波谱库进行重采样，使之与所提取的端元波谱相匹配。采用波谱特征拟合方法（SFF）、波谱角分类（SAM）方法和二进制编码方法（BE），并对3种匹配方法的权重进行不同的设置，将获取的端元波谱与重采样后的USGS植被波谱库的物质波谱进行光谱匹配分析。根据匹配结果评分确定两种波谱的相似度，进而确定端元的物质类别。很多树种在一个波长范围内是相似的，在其他范围内可能具有较大差异，因此应该使用含有诊断性特征的波长范围进行分析。表1为提取的第1类端元波谱同USGS植被波谱库进行光谱分析后的得分，核桃树与未知端元波谱相似度最高，得分为1.859，因此将之判定为核桃树。

2.5 最小能量约束法（Constrained Energy Minimization，CEM） CEM的原理是使用有限脉冲响应线性滤波器（Finite Impulse Response，FIR）和约束条件来最小化平均输出能量，以抑制影像中的噪声和非目标端元光谱信号，即抑制背景光谱，定义目标约束条件来分离目标光谱。CEM使用MNF变换结果作为输入文件，计算过程中选择协方差矩阵计算。每个端元波谱对比每个像元得到一幅灰度图像，值越大表示越接近目标。选择图像值最大的部分即可得到核桃树的分布范围（图5）。

3 结果与分析

运用上述方法，对研究区域Hyperion影像进行数据处理与分析研究，尝试进行树种识别与填图工作，共提取6类植物端元光谱，经过光谱分析，将6类端元识别出来，分别为核桃树、矮灌木、橡树、刺柏、冷杉和灌木蒿。利用矿物填图的方法对识别出的6类植物进行分离，得到6类植物的分布图（图6）。

为了验证识别结果的准确性，在研究区域进行了实地调查。在感兴趣的区域分别采集了24个点，其中21个为实采点，3个点为因地形原因无法达到而根据实际位置进行估计得到的点。利用ENVI的混淆矩阵对分类结果进行验证，得到分类精度统计结果（表2）。

由表2可知，第1类端元核桃树的识别结果很准确，实地采集点都落在分类范围内，然而实际核桃林分布面积要比分类结果小很多。由于影像获得时间较早（2004年）和近年来核桃林面积的缩小，分类精度显著高于实际。

第2类端元矮灌木的种类识别精度为80%，其分布范围与实际较为吻合，在低海拔山坡下部也有同类矮灌木分布，然而面积较小。

第3类端元识别为橡树，在实地调查中发现，在相应区

图6 6类端元分离填图结果域内种植有青冈树，青冈树为栎属，别称橡树、橡子树、栎树，这表明其被识别为橡树的结果是具有一定准确性的。得到的青冈树分布精度为60%，青冈树在离沟道较近的山坡坝上以及村后的缓坡上也有分布，但在识别范围中未被识别出来，其面积较小以及与巨大杨树混交是其未被识别的主要原因。

第4类端元识别为刺柏，分布在海拔较高的地带，由于地形原因未能进行实地采点验证，通过山中村民了解到在海拔较高地带确实分布大片柏树林，其实际分布范围不详。

第5类端元识别为冷杉，在研究区域内冷杉以成片树林或者单独一棵或数棵的形式分布，在分类结果所示冷杉识别精度为71.4%，海拔较低的山坡上成片冷杉林以及零星分布的冷杉均未被识别，这与其分布面积较小或者零星存在有关。

第6类端元识别为灌木蒿，在实际调查中该类植被主要分布在废弃或休耕的田地中，或者在因崩塌或滑坡产生的缓坡上，识别精度为80%，一些面积较小或者植被覆盖度较低的崩塌蒿草区未被识别。在沟谷的两侧山坡上长有大量的灌木和蒿草，其种类繁多，分布交错复杂，在光谱特征上没有反映出诊断性特征，所以未能得出有意义的结果。另外，在沟谷内和较低的缓坡上存在柿子树、柳树以及杨树，但其分布范围非常小所以未能将其识别。

4 结论

该研究利用Hyperion高光谱影像，采用基于MNF与PPI的端元提取技术，在研究区域内共提取出6类端元，并运用波谱分析技术，识别出核桃、矮灌木、橡树（青冈树）、刺柏、冷杉和灌木蒿，然后利用CEM进行混合像元分解，将6种树木和灌草分离出来，并得到其分布范围。在种类识别上与实际情况较为吻合，在分布范围上平均分类精度可以达到78.3%。由于森林树种的光谱具有极大的相似性，环境、地形因素以及不同树种混交的相互效应使树种的光谱特征产生变异[23]，而且在光谱分析中使用的是美国植物标准波谱数据库，它与我国树种存在的差异性使得光谱识别变得困难，因此建立和完善我国地区的树种标准波谱数据库在高光谱树种识别工作中尤为重要。由于Hyperion高光谱影像空间分辨率不高（30 m），在野生森林中多存在树木混交现象，所以在端元提取以及混合像元分解上只利用了较少的光谱信息，尚有未被提取和识别的植被，因此对Hyperion高光谱数据的挖掘和分析仍有较大的空间。

参考文献

[1] 刘向东，吴钦孝，赵鸿雁.森林植被垂直截留作用与水土保持[J].水土保持研究，1994，1（3）：8-13.

[2] 陈景武我国暴雨泥石流预报研究[C]//第四届全国泥石流学术会讨论论文集.兰州：甘肃文化出版社，1994.

[3] 刘旭升，张晓丽.森林植被遥感分类研究进展与对策[J].林业资源管理，2004（1）：61-64.

[4] 谭炳香.高光谱遥感森林类型识别及其郁闭度定量估测研究[D].北京：中国林业科学研究院，2006.

[5] 谭炳香.高光谱遥感森林应用研究探讨[J].世界林业研究，2003，16（2）：33-37.

[6] 谭炳香.高光谱遥感森林信息提取研究进展[J].林业科学研究，2008，21（S1）：105-111.

[7] MARTIN M E，NEWMAN S D，ABER J D，et al.Determining forest species composition using high spectral resolution remote sensing data[J].Remote Sensing of Environment，1998，65（3）：249-254.

[8] DARVISHSEFAT A A，KELLENBERGER T W，ITTEN K I.Application of hyperspectral data for forest stand mapping[C]//Symposium on Geospatial Theory.Otawa，Processing and Applitions，2002.

[9] ASPINALL R J.Use of logistic regression forvalidation of maps of the spatial distribution ofvegetation species derived from high spatial resolution hyperspectral remotely sensed data[J].Ecological Modelling，2002，157：301-312.

[10] 浦瑞良，宫鹏.高光谱遥感及其应用[M].北京：高等教育出版社，2000.

[11] 刘秀英，林辉，熊建利，等.森林树种高光谱波段的选择[J].遥感信息，2005，80（4）：41-45.

[12] 陈尔学，李增元，谭炳香等.高光谱数据森林类型统计模式识别方法比较评价[J].林业科学，2007，43（1）：84-90.

[13] 王志辉，丁丽霞.基于叶片高光谱特性分析的树种识别[J].光谱学与光谱分析，2010，30（7）：1825-1829.

[14] 王志辉.高光谱遥感在森林树种识别中的应用[D].杭州：浙江农林大学环境与资源学院，2011.

[15] 魏新功，王振国，包红霞.降水原因造成的舟曲县地质灾害分析[J].甘肃科技，2008，24（21）：84-88.

[16] PAMELA BARRY.EO-1/ Hyperion Science Data Users Guide[R].TRW Space，Defense & InformationSystems，One Space Park Redondo Beach，CA 90278，2001.

[17] Richard Beck，EO-1 User Guide v.2.3[EB/OL].http：//eo1.usgs.gov & http：//eo1.gsfs.nasa.gov，2003.

[18] 谭炳香，李增元，陈尔学，等.EO-1 Hyperion高光谱数据的预处理[J]. 遥感信息，2005（6）：36-41.

[19] 梁继，王建.Hyperion高光谱影像的分析与处理[J].冰川冻土，2009，31（2）：247-252.

[20] 周清，祝民强.EO-1Hyperion高光谱数据FLAASH模块大气纠正及评价[J].测绘与空间地理信息，2011，34（6）：149-151.

[21] KRUSE F A，BOARDMAN J W，HUNTINGTON J F.Comparison of airborne hyperspectral data and EO-1 hyperion for mineral mapping[J].IEEE Transactions on Geoscience and Remote Sencing，2003，41（6）：1388-1400.

[22] GOODENOUGH D G，DYK A，NIEMANN K O，et al.Processing hyperion and ALI for forest classification[J].IEEE Transactions on Geoscience and Remote Sencing，2003，41（6）：1321-1331.

[23]曾慶伟，武红敢.基于高光谱遥感技术的森林树种识别研究进展[J].林业资源管理，2009（5）：109-114.安徽农业科学，Journal of Anhui Agri. Sci.2014，42（16）：5324-5325责任编辑 胡剑胜 责任校对 况玲玲