Sentinel-2A卫星数据处理分析及在干旱河谷提取中的应用

杨 斌， 李 丹， 高桂胜， 陈 财， 王 磊

(西南科技大学环境与资源学院，绵阳 621010)

0 引言

哨兵-2A(Sentinel-2A)卫星作为欧洲空间局(European Space Agency，ESA)“哥白尼”计划系列中的一颗重要的光学遥感卫星，于2015年6月23日发射升空，卫星运行在高度为786 km、倾角为98.5°的太阳同步轨道上。Sentinel-2系列计划发射2颗卫星，构建双星对地观测平台，主要用于全球生态环境的监测[1]，并且该系列卫星与Landsat系统在对地观测领域互为补充，具有较高空间分辨率、免费共享使用等特点，因而对该卫星的深入研究与应用具有重要意义。

Sentinel-2A卫星具有13个光谱波段，涵盖可见光、近红外和短波红外波段，拥有10 m，20 m和60 m3个不同尺度的空间分辨率，重访周期为10 d，Sentinel -2B卫星发射升空后，Sentinel-2这2颗卫星每隔5 d就可对全球N84°～S56°范围内的所有地物覆盖一遍[2]。Sentinel-2A卫星自发射升空后经过近半年的运行和调试，从2015年12月3日起，正式向全球用户提供免费下载服务[3]，国外研究人员前期则曾模拟Sentinel-2数据波段进行了小麦土豆的叶面积指数提取[4]及植被叶绿素含量的评估[5-6]，而后又利用Sentinel-2A卫星数据进行了温室效应变化监测[7]、红边光谱指数适应性分析[8]和云检测分析[9]，而国内研究人员对于Sentinel-2A卫星数据的应用研究还不多见。为进一步了解认识Sentinel-2A卫星影像数据的基本特点，通过对该影像进行数据处理与分析，选取研究区域尝试运用该数据相关指标进行地学环境领域综合应用，为了解该遥感卫星影像数据的应用潜力提供较为全面的基础信息。

1 Sentinel-2A卫星概述

1.1 卫星数据基本特征

Sentinel-2A卫星属于全新的高空间分辨率多光谱成像卫星，主要用于包括陆地植被[11]、土壤以及水资源、内河水道和沿海区在内的全球陆地观测以及紧急救援服务[12-14]。选取目前应用较为广泛且参数相近的陆地资源观测遥感卫星Landsat8，分别从波段设置、波长范围、幅宽大小、空间分辨率、时间分辨率和辐射分辨率等基本特征参数对这2种卫星数据进行对比，如表1所示。

表1 2种卫星数据主要参数对比Tab.1 Comparison of spectral bands data between two satellites

从表1中可以得出Sentinel-2A数据除在空间分辨率、时间分辨率及幅宽方面具有一定的优势之外，在红光边缘波段、近红外波段均比Landsat8数据具有更细致、更全面的波谱分辨率[15-16]，其中，Sentinel-2A波段中的B1是气溶胶波段，B5，B6，B7和B8a是植被红边波段，B9是水汽波段，B10是卷云波段，B11和B12适合于识别提取雪、云、冰信息(图1)，由此可见，利用Sentinel-2A数据在全球植被生态特征监测与分析方面具有更强的应用潜力。

图1 Sentinel-2A卫星数据波段分布Fig.1 Band distribution of Sentinel-2A satellite

1.2 数据组织形式及产品格式

Sentinel-2A卫星数据的分发处理模式与Sentinel系列卫星相一致，分别通过地面控制中心、地面数据处理中心及数据分发中心3个环节进行。按照ESA发布的Sentinel-2A数据产品等级依次划分为L0，L1A，L1B，L1C和L2A这5种不同层次的数据格式[17]，其中L0级别为卫星获取的初始压缩数据，并提供生产L1A级别数据所需的所有信息，L1A级别是在L0级别基础上，对所有波段进行解压，并对短波红外波段进行重新排列处理； L1B级别数据是在L1A级别数据的基础上，经过辐射校正，并基于有效的地面控制点构建物理几何精校正模型，并附带到该级别产品数据中； L1C级别数据则是在L1B级别数据的基础上，通过重采样、正射校正和大气表观反射率计算处理等操作，生产云、陆地/水体掩模数据。其中，L0，L1A和L1B级别的数据不对用户发布，L1C级别数据则通过ESA SciHub网站对用户免费公开发布。

ESA对外共享发布的L1C级别数据都是按照100 km2作为最小瓦片单元进行存储，该数据经过UTM/WGS84投影正射校正处理，且最基本的文件单元包括影像数据文件(IMG_ DATA )、影像数据辅助文件(AUX_DATA)、影像质量指标文件(QI_DATA)以及元数据信息(Metadata)，对L1C级别基本数据的说明分析如表2所示。

表2 L1C级别数据说明Tab.2 Data illustrate of L1C basic unit

根据ESA发布的Sentinel-2A用户手册指出，L1C产品的获取、处理、归档及发布均由Sentinel-2卫星数据地面数据分发中心提供，而 L2A级别数据的生成处理与应用分析则需用户自行利用L1C级别数据和Sentinel-2 Toolbox(SNAP软件工具箱)处理完成。

1.3 L2A级别数据生成

通过对Sentinel-2A用户手册分析得出，L1C级别数据虽然经过了正射校正、大气表观反射率处理等操作，但在光谱特征、地形变换等方面仍然具有一定的缺陷，而L2A级别数据处理则基于查找表法(LUT)利用已经计算出的不同复杂大气条件下的气溶胶光学厚度、反射率和相函数等指标，进行大气下垫面发射率校正处理，最大程度考虑了大气对遥感数据的影响(即大气精校正)。由于 L2A级别数据是在L1C级别数据基础上用户自行生成，因此需对L2A级别数据生成进行探讨。

L2A级别数据生成是基于L1C级别数据调用SNAP软件Sentinel-2 Toolbox中的sen2cor处理模块，进行大气校正和几何精校正，生成精度更高的大气下垫面反射率影像数据，还附加生成了气溶胶厚度分布图、水汽分布图和场景分类图等环境质量分析数据[18]。Sentinel-2 Toolbox主要针对于Sentinel-2系列光学影像数据开发的可视化遥感处理工具箱软件，其处理操作可通过桌面可视化和命令行运行界面2种方式调用，并提供了基于Java和 Python语言的开发接口与环境，便于用户自行开发定制。根据ESA网站信息及发布的Sentinel-2卫星数据用户指南[19]，Sentinel-2 Toolbox集成在SNAP软件中，而在处理生成L2A级别数据之前，需先安装SNAP软件再运用Python语言命令行安装sen2cor处理模块，然后在SNAP软件中调用sen2cor进行L2A级别数据的生成。

生成的L2A级别数据除完成大气下垫面反射率校正(bottom-of-atmosphere，BOA)、卷云计算(cirrus)和地形校正(terrain)等操作外，还将波段保留为10个(去除B1，B9和B10)，且空间分辨率可以根据需求选择10 m，20 m和60 m的任何一种类型，数据的精度和应用价值均得到了显著提升[19]。此外，SNAP软件还可利用L2A级别数据生成叶面积指数(leaf area index，LAI)、吸收光合有效辐射比例(fraction of absorbed photosynthetically active radiation，FAPAR)、叶绿素含量(chlorophyll content in the leaf，CAB)、冠层含水量(canopy water content，CWC)和植被覆盖度(fraction of vegetation cover，FVC)等生物物理参数，可为大面积进行植被生态环境的潜力研究与评价提供重要数据。

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

2.1 研究区概况及数据获取

本文选取四川省茂县与黑水县交界区域为研究区。该区经纬度坐标范围为E103°17′～103°44′，N31°47′～31°58′，长41.61 km，宽19.89 km，总面积约848 km2。该区地处岷江上游黑水河支流区域，具有山地地形和青藏高原季风气候特征，干旱河谷比较发育。

遥感数据采用2016年5月5日Sentinel-2A L1C级别数据，投影/坐标系为UTM/WGS84，整景数据的含云量为14.84%，但研究区范围内云量较少，符合研究需求。

研究区范围及其影像如图2所示。

(a) Sentinel-2A B4(R)，B3(G)，B2(B)真彩色合成影像

(a) 研究区DEM

图2研究区Sentinel-2A数据范围及其遥感影像

Fig.2DistributionofSentinel-2Adataintheresearcharea

2.2 数据预处理

利用裁剪后的L1C级别数据，将其导入SNAP软件中，调用SNAP软件中sen2cor处理模块，选取输出数据精度为10 m，分别对研究区数据进行大气下垫面发射率校正、卷云和地形校正处理，生成研究区L2A级别数据及气溶胶光学厚度、水汽分布和场景分类信息，如图3所示。对处理后的L2A影像数据各波段基本特征值统计结果如表3所示。

(a) Sentinel-2A B8(R)，B4(G)，B3(B)假彩色合成影像 (b) 气溶胶光学厚度

(c) 水汽分布 (d) 场景分类

从表3可以看出，B2，B3，B4和B8的标准差较大，在第一主成分分量中所占比例也较高，说明其信息量较为丰富，数据较为稳定，从而可得出红边波段在第一主成分分量中所占比例较为稳定，对于生态环境监测具有较好的应用价值。

2.3 数据产品分析

分别从数据生成机理、流程与应用效果等方面对L2A级别数据生成的气溶胶光学厚度、水汽分布和场景分类数据类型与影像数据的光谱特征进行重点分析[20]。

1)气溶胶光学厚度图层是利用浓密植被算法(dense dark vegetation)，使用短波红外波段(B12)、红光波段(B4)和蓝光波段(B2)生成，该图层能有效反演大气透明度。图3(b)中可反映出研究区内高山无云区的气溶胶光学厚度较低。

2)水汽分布图层是利用大气预校正微分吸收算法(atmospheric pre-corrected differential absorption，APDA)，使用Sentinel-2A数据中的窄近红外波段(B8a)和水汽波段(B9)生成。其中，将B8a设置为大气窗口区参考通道，将B9 设置为水汽吸收过程中的测量通道，吸收深度的计算是在假设测量通道的表面反射率与参考通道相同的情况下获得。图3(c)中可反映出研究区河谷区内水汽含量较大。

3)场景分类图层是利用阈值测试法(threshold test)和光谱大气表观反射率特性，检测出研究区云、雪、阴影和地物类型。其中算法可获得4种不同类型的云类型(低概率云层、中概率云层、高概率云层和薄卷云)，并可生成云层阴影、植被、裸地、水体、雪和黑暗区域(盲区)等类型。图3(d)中显示出了研究区内10种不同的场景类型。

4)获得的L2A级别影像数据波段将原有的13个波段变为10个波段。其光谱特征也经过了大气下垫面发射率校正、卷云和地形校正处理，从生成的20 m空间分辨率L2A级别数据[21]中分别获取研究区内同一地点的水体、植被和裸地3种不同地物光谱特征曲线。图4中反映出L1C和L2A级别数据在光谱特征中的差异性。根据经验得出，L2A级别数据的光谱特性与真实地物光谱特性更为接近。

(a) L1C级别数据(b) L2A级别数据

图4典型地物光谱特征分析

Fig.4Spectralcharacteristicsanalysisoftypicalobject

3 生物量提取分析

基于L2A级别数据，利用SNAP软件中生物物理量处理器模块[22](biophysical processor)生成LAI，FAPAR，CAB，CWC和FVC这5种因子，并利用波段运算方法提取出研究区NDVI(图5)。该模块是将事先处理好的区域植被特征与Sentinel-2数据冠层表观(top of canopy)反射率相关联构建一个综合数据库，并通过训练定标好的神经网络去评估遥感图像中每个冠层的特性，而模型中每个神经网络元的构成则涉及到L2A级别数据中的B3，B4，B5，B6，B7，B8a，B11，B12和天顶角余弦值、太阳高度角余弦值、相对方位角余弦值等11个输入参数和5个具有正切S形曲线(sigmoid)传递函数的神经元辅助参数。此模型算法最大程度精细计算出各植被特征与冠层表观发射率之间的关系，可为生态环境定量遥感分析提供强有力的数据支持。

(a) LAI (b) FAPAR

(c) CAB (d) CWC

(e) FVC (f) NDVI

图5研究区植被生物量因子数据分布

Fig.5Distributionofvegetable-biophysicalindexdataintheresearcharea

得出的5种生物物理量因子中，LAI是众多植被-大气相互作用模型(特别涉及到碳循环和水循环模型)中的一个关键参数[23]，图5(a)反映出河谷地带的LAI显著低于海拔较高的山坡及山顶地区； FAPAR是描述植被结构以及与之相关的物质与能量交换过程的基本生理变量，也是遥感估算陆地生态系统植被第一性生产力(net primary productivity,NPP)的重要参数，该指标也能有效表征植被生长状态和演化过程，图5(b)中也反映出云雪与植被的显著反差，还能有效显示出河谷地带植被长势较差； CAB能够真实有效地指示出地表植被的生长和营养状况(受氮素胁迫程度)、光合作用能力、所处的生长发育阶段及是否受到自然环境破坏或病虫害影响等状态信息，图5(c)中CAB对植被特征区分精细度较高，并反映出黑水河上游区域(左下方)植被长势较好； CWC能够有效评估植被的生长状况，对干旱监测和生态环境改善等具有现实意义，但图5(d)中CWC对于区域植被特征显示不明显； 而FVC直接表征地表植被覆盖的多少，能大致反映出地表资源与生态环境的好坏程度，图5(e)中FVC能显著区分出裸露地区与植被覆盖度较高地区的差异性。综上所述，这5种因子对于区域植被情况特征分布均具有一定的指示作用，但具体的应用效果还有待进一步探讨。

4 干旱河谷提取与分析

由于该研究区位于青藏高原东缘地带的岷江上游流域，长年受西南横断山区季风气候的影响，气候温暖、干燥度与辐射度强、雨热同步等特点，且河谷区内植被发育情况较高海拔山坡差，形成了山地典型垂直分异植被特性。从遥感假彩色图像上可明显看出干旱河谷区内影像呈现出暗灰色调，其边界在图3(a)中的显示较为明显。

通过Sentinel-2A提取出的5种生物物理量因子和NDVI值能有效反映研究区内植被发育长势情况，对干旱河谷边界的定量分析与提取具有较好的指示作用[24]。选取基于知识经验的决策树分类法，将上述5种生物物理量因子和NDVI参数导入到ENVI软件中进行图像特征分析，如表4所示。

表4 各植被生物量图像特征分析Tab.4 Analysis of image feature based on every vegetable-biophysical data

结合各指标影像特征反映得出标准差较高且信息量较大的因子为CAB。并通过查阅该区域干旱河谷前期研究资料和野外实际调研[25]，得出研究区范围内干旱河谷分布与海拔高度分布具有极其紧密的联系，该区干旱河谷海拔主要集中分布在600～2 400 m区间的河谷范围内[25]。利用研究区30 m空间分辨率的DEM(从美国USGS网站下载)和CAB因子，分别对比假彩色影像选择该因子数据的合理阈值(本文阈值选为32.368 9)，在ENVI软件中调用决策树分类法，按图6所示流程定量提取出研究区的干旱河谷分布边界与面积，其结果如图7所示。

图6 决策树分类流程Fig.6 Flowchart of decision tree classification

图7 干旱河谷提取分布Fig.7 Distribution of dry river valleys

通过上述分析与提取，将决策树分类结果导入到ArcGIS 软件中，运用ArcGIS空间分析功能，计算出该研究区内干旱河谷的面积为114.147 7 km2，占研究区总面积的13.46%。对干旱河谷空间分布进行的分析表明，干旱河谷区主要集中分布在黑水河河谷坡度小于25°的范围内。对照该区域的地质图，可以看出干旱河谷区内灰岩、千枚岩和石英岩较为发育。结合包维楷等[26]、杨兆平等[27]前期研究成果与作者以往的研究积累及野外实地考察，选取6个野外观测样点(图7)，进行实地观测验证，提取精度达到90%以上。研究结果表明，干旱河谷提取结果能较好表征干旱河谷分布情况，为岷江上游干旱河谷大面积国土资源普查和生态环境变迁分析研究提供了方法支持。

5 结论

Sentinel-2A光学遥感卫星作为对地观测陆地资源卫星系列大家庭的一员，以其涵盖面积大、时间和空间分辨率高、数据免费共享等优势受到广大遥感应用工作者的广泛关注。本次遥感应用的主要结论如下：

1)从Sentinel-2A遥感数据的基本属性、光谱特征、产品特征和应用特征4个方面进行系统分析，利用SNAP软件sen2cor模块进行数据处理，处理结果可以较好满足后续应用需要。

2)以四川省岷江上游流域黑水河区域为研究区进行了遥感生物量因子提取分析，并利用分析后的生物量因子和DEM数据，采用决策树分类法模型进行了干旱河谷区域定量提取。结果表明，Sentinel-2A系列卫星数据光谱特征丰富、应用前景广阔、数据可扩展性较强，可为植被生态环境定量分析和地学综合应用提供更多的遥感数据资源。

志谢： 感谢ESA提供的Sentinel-2A数据资源及学习文档，感谢http: //forum.step. esa.int/c/s2tbx/sen2cor论坛网站上国内外同行给予的帮助！