李海翠 况润元 宋子豪
(江西理工大学土木与测绘工程学院,赣州 341000)
水色遥感是根据水体在可见光或近红外谱段的吸收与散射的光谱特性,利用传感器测得的离水辐亮度对水体的有关现象与参数进行解译的遥感技术,该项技术可以对各项水色参数做出评价并提供决策服务[1]。水色卫星传感器接收到的总辐射量 80%以上来自大气,而来自水体的辐射只有 3%~15%[2],大气辐射贡献是卫星遥感数据反演水体悬浮泥沙和叶绿素 a浓度等水色成分的限制因素之一[3],因此,消除大气辐射的干扰来获得有效的离水辐射信号,实现卫星影像的大气校正,是水色遥感信息提取的关键技术之一[4]。Morel和Prieur[5]根据光学性质将水体分为两类:Ⅰ类水体(多为开阔大洋)和Ⅱ类水体(多为近岸海域或内陆湖泊),Ⅰ类水体的光学特性主要是由水中叶绿素所决定的,Ⅱ类水体是由悬浮无机物、有色可溶性有机物和叶绿素共同决定的。目前针对Ⅰ类水体已有较为成熟的清洁水体的大气校正算法[2],而由于Ⅱ类水体组成复杂、光学特性多变,适用于Ⅰ类水体大气校正的假设条件已不成立,Ⅱ类水体通用的大气校正业务化算法还需要深入研究[6]。
与中、低分辨率的卫星传感器(如:中分辨率成像光谱仪)相比,多光谱成像仪(Multispectral Imager,MSI)得到的影像卫星可以重采样至10m的高空间分辨率,具有更好的空间信息,与其它中分辨率传感器相比,MSI具有更好的信噪比,并且有短波红外波段用于大气校正,该卫星是欧空局“哥白尼计划”中发射的第2组卫星,其搭载的MSI载荷延续并扩展了陆地卫星系列(Landsat,包括ThematicMapper-TM、Enhanced Thematic Mapper Plus-ETM+和Operational Land Imager-OLI三种传感器)卫星和SPOT系列卫星的对地观测任务[7]。目前国外专家学者已经将“哨兵二号”(Sentinel-2)卫星MSI影像应用到II类水体的大气校正研究中。Souza Martins等[8]基于Sentinel-2 MSI影像数据对亚马逊河漫滩湖泊的大气校正方法进行研究,分析评估6S(Second Simulation of the Satellite Signal in the Solar Spectrum)、ACOLITE(Atmospheric Correction for OLI lite model)和Sen2cor(Sentinel 2 correction)对该湖泊的适用性。Marcela和 Ana等[9]基于Sentinel-2 MSI数据对西班牙内陆水域进行大气校正算法评估,对ACOLITE、C2RCC(Case 2 Regional Coast Color)、iCOR(image correction for atmospheric effects)、Sen2cor和 Polymer(Polynomial based algorithm applied to MERIS)等五种大气校正方法的精度进行评估分析。但在国内,利用Sentinel-2 MSI等高分辨率影像对Ⅱ类水体的大气校正研究还比较少,国内有不少学者利用Sentinel-2 MSI影像数据对不同地物以及太湖进行研究,但都没有探究大气校正方法在鄱阳湖的适用性如何。
因此,本文以鄱阳湖为研究区域,以Sentinel-2 MSI卫星影像的大气上层表观反射率(TOA)产品为数据源,利用 Sen2cor、6S、FLAASH(Fast Line-of-sight Atmospheric Analysis of Hypercubes)、ACOLITE和QUAC(Quick Atmospheric Correction)[10]五种大气校正方法对鄱阳湖水体进行大气校正研究,并通过地面同步实测光谱数据进行对比分析,对五种大气校正方法后得到的反射率精度进行验证,综合评价各种算法的优缺点,得出适用于鄱阳湖的大气校正方法,以期为提高未来水质参数反演精度提供参考,为环境管理部门提供数据支持。
鄱阳湖是中国最大的淡水湖泊,位于江西省北部(东经115°49′~116°46′,北纬 28°11′~29°51′)(如图1 所示),长江中下游交接处的南岸,是长江流域一个重要的过水性、吞吐型、季节性的内陆湖泊[11]。
图1 研究区域与实测站点图Fig.1 Diagram of research area and measured site
Sentinel-2B于2017年3月7日发射,该卫星主要应用于陆地和沿海地区监测和应急服务的极轨光学任务,与之前的SPOT卫星以及Landsat卫星等类似光学传感器相比,其光谱范围和性能都有所提高。该卫星装备有一个多光谱仪器(MSI),具有从蓝色(Blue)到短波红外(Shortwave Infared,SWIR)的13个光谱波段(如表1所示),地面分辨率分别为10m、20m和60m,MSI的目标是测量13个光谱波段的反射率[12]。标准MSI影像是1C级产品(L1C),具有UTM/WGS84(World Geodetic System 1984世界大地测量系统)投影且进行了辐射和几何校正。本文选用 2018年 10月 29日的Sentinel-2B L1C的数据,这3景影像完全包含研究区域鄱阳湖,并且研究区域以及其周围地区上空为晴朗无云。关于Sentinel-2B MSI光谱响应函数(如图2所示)的信息可以从欧洲航天局(the European Space Agency,ESA Sentinel-2 MSI Technical Guides)获得。在大气校正前需对遥感影像进行重采样(利用“哨兵”系列卫星数据应用平台重采样至10m分辨率)、图像镶嵌、辐射定标以及图像裁剪等预处理工作。
图2 Sentinel-2B VNIR光谱响应函数Fig.2 Sentinel-2B VNIR spectral response function
表1 Sentinel-2B MSI传感器相关参数Tab.1 Sentinel-2B MSI sensor parameters
实测数据与影像的时间差会影响反射率对比,文献[8]的结果指出,在水和环境条件没有出现快速变化的情况下,使用±3天或至多±8天的测量数据进行比较是合理的。本文所使用的实测数据为2018年10月28日及29日所测得,因此Sentinel-2 MSI图像和现场测量数据之间的时间差距原则上不被考虑作为比较的因素。通过品质控制选取 13个站点(O1~O13)数据(见图1),将其作为大气校正结果的检验数据,利用美国ASD(Analytical Specral Devices)公司生产的便携式光谱仪测量水体遥感反射率。测量是在晴朗且风速低的天气条件下开展,每个站点至少测量十次以保证测量的精确度[13]。本文用遥感反射率来分析比较校正精度,水体遥感反射率Rrs定义为
式中Lw为离水辐亮度,Lw=Lsfc-rLsky,其中Lsfc为光谱仪测量的海面辐亮度,r是气水界面对天空光的反射率(取定值0.026),Lsky为天空方向的辐亮度;Es()λ为海面入射辐照度,为波长,Lp()λ为测得的标准板辐亮度,Pp()λ为标准板发射率,要求在10%~35%之间[14]。
式中yi代表Sentinel-2 MSI传感器第i波段的实测光谱反射率;λsi是波段i的起始波长;λei是波段i的终止波长;Oi(λ)是波段i在波长λ处的光谱响应函数。
在测量得到的光谱曲线中,选择400~900nm的数据计算遥感反射率,利用地面实测数据和Sentinel-2 MSI的光谱响应函数,按照式(2)计算模拟得到B2~B8A对应波段的等效遥感反射率。实测遥感反射率和模拟得到的等效遥感反射率如图3所示。从图3中可以看出,13个站点的光谱曲线形状基本一致,但遥感反射率数值大小因站点水体浑浊程度不一样而具有明显差异。
图3 实测遥感反射率和拟合后MSI波段反射率Fig.3 Measured remote sensing reflectance and MSI band reflectance after fitting
本文采用Sen2cor、6S、FLAASH、ACOLITE以及QUAC五种大气校正方法。由于大气校正过程中需要输入较多参数,为了保证校正结果的可比性,在参数设置上保持一致。根据获取的研究区Sentinel-2数据时相选取中纬度夏季(mid-latitude summer,MLS)和农村气溶胶,能见度数据采用默认数据为40km,地物高度采用鄱阳湖多年平均水位12.68m[15]。现简单介绍五种大气校正方法的原理如下:
Sen2cor是欧空局官方提供的Sentinel-2卫星影像的大气校正方法,用于生成MSI陆地产品(级别2A),该算法是一种半经验算法,集成了基于图像检索和LibRadtran(辐射传输计算软件包)[13]模型的查找表(Look-Up-Tables,LUTs),以消除MSI图像的大气影响[16]。6S模型[17]是在1986年Tamre等人开发的 5S模型的基础上进行了改进,该模型对大气散射率、大气吸收率、大气透过率的计算进行改进,增加交叉辐射项及双向反射分布函数(Bi-direction Reflectance Distribution Function,BRDF)估值算法,对辐射传输方程的求解算法进行了效率优化。FLAASH模型是集成在遥感图像处理平台(The Environment for Visualizing Images,ENVI)中,该模型是对MODTRAN4+(MODerate spectral resolution atmospheric TRANsmittance)[18]模型改进后提出的,可以有效的消除大部分大气反射地面物体的影响并获得准确的地面反射特性,适于对 0.4~3波长范围内影像的大气校正[19]。ACOLITE[20]算法用于研究 Landsat 8-OLI和Sentinel-2 MSI图像的大气校正,适用于海洋和内陆水的研究,此算法能有效去除分子和气溶胶组分在清澈和浑浊水域中的散射效应,算法中的瑞利散射是利用6SV(a version of second simulation of the satellite signal in the solar spectrum)模型的LUTs来校正,而气溶胶散射是根据清澈水体的近红外波段(842nm和865nm)和浑浊水体的SWIR波段(1 610nm和2 130nm)来进行估计的[21]。QUAC算法是基于暗目标的算法,直接根据场景中观测到的像元光谱确定大气校正参数,而不需要辅助信息,该算法集成在ENVI软件中的QUAC工具中,通过自动从图像上收集不同物质的波谱信息来获取经验值以完成影像的快速大气校正[22-23]。
为了比较Sen2cor、6S、FLAASH、ACOLITE和QUAC五种大气校正算法的效果。选取实测光谱拟合后的同步等效遥感反射率共计13个点来评价五种大气校正结果的精度,利用决定系数R2(用于评估校正值和实测值的符合程度),和平均相对误差MRE来验证校正精度。
式中yi为ASD实测的第i个波段遥感反射率;′是大气校正后的第i个波段遥感反射率;是实测遥感反射率平均值;n是样点总数目。(注:R2可以取0到1的值。R2越接近1,大气校正值与地面真实值之间的相关性越强,表明该方法越好。)
利用Sen2cor、6S、FLAASH、ACOLITE和QUAC方法对鄱阳湖Sentinel-2 MSI影像进行大气校正,以鄱阳湖 2018年 10月 28日和 29日的同步实测光谱数据为基础,对校正结果进行对比分析,如图4所示。
图4 校正后遥感反射率与实测遥感反射率散点图Fig.4 Scatter diagram of corrected remote sensing reflectance and measured remote sensing reflectance
图4(图4(a)~(e)横坐标对应的是yi,纵坐标对应的是yi′,图4(f)横坐标对应的是λ,纵坐标对应的是Rrs)是五种大气校正方法的校正结果与实测遥感反射率的对比图。由图4可知,6S、Sen2cor、ACOLITE方法在反射率为0~0.03范围之内,校正后遥感反射率整体高于实测遥感反射率,在 0.03~0.06范围之内则校正后遥感反射率整体低于实测遥感反射率;FLAASH算法在反射率在0~0.02范围内,校正后遥感反射率整体高于实测遥感反射率,在 0.02~0.06范围内,整体低于实测遥感反射率。QUAC算法校正后遥感反射率整体都高于实测遥感反射率。整体而言Sen2cor算法(图4(a))在这五种算法中校正精度最高,其平均相对误差为29.37%,决定系数R2为0.969 41,其次是FLAASH算法(图4(c)),平均相对误差是30.11%,决定系数R2为0.939 11,ACOLITE算法(图4(d))和6S算法(图4(b))校正结果相似,其平均相对误差分别为52.91%和48.14%,决定系数分别为0.930 3和0.929 8。QUAC算法(图4(e))大气校正结果表现效果最差,平均相对误差为53.03%,决定系数R2为0.908 92。为了更好的查看五种大气校正方法之间的对比,特随机选取某一采样点来进行说明,图4(f)是3号采样点五种大气校正方法与实测遥感反射率之间的对比图,可以得出QUAC算法整体比实测光谱曲线高,而其他四种算法均与实测光谱曲线较为接近。
对于单波段,从表2和图5中可以看出6S算法和ACOLITE算法都在560nm和665nm的性能是最佳的,Sen2cor算法的整体性能都较好,最佳性能在443nm和490nm;FLAASH算法各波段都处于中等偏上的性能,最佳性能在665nm和842nm;QUAC算法的整体校正效果不佳,但在492nm其性能是最佳的。在865nm波段,FLAASH算法和Sen2cor算法误差相对较小,但也都超过70%,其他三种算法的平均相对误差都高于100%,可以看出五种算法在865nm的校正结果并不理想。在B2~B5可见光波段,五种算法的校正精度都较高,从B6波段过渡到B7~B8A(近红外波段)其校正精度都较差,主要原因是:由于获取的近红外波段的光谱数据信噪比比较低,获取的水体信号中干扰比较严重[24]。
表2 各波段大气校正精度Tab.2 Atmospheric correction accuracy of each band
从图5中可知,在490~705nm波段范围内,除QUAC算法之外的四种算法的校正后遥感反射率整体都低于实测光谱遥感反射率。在560nm处,离散程度远高于490nm、665nm以及705nm,其次是490nm波段,665nm和705nm更靠近1:1线。在740~865nm波段范围内,五种校正算法的反射率都高于实测光谱的遥感反射率,波长越长其离1:1线就越离散,其校正精度就越差。
图5 各波段散点图Fig.5 Scatter diagram of each band
上述结果的误差不仅包含大气校正方法校正误差,还包含实测数据的测量误差、测量期间水体平静状态以及定标结果本身误差三方面的误差。实测数据的测量误差经过计算,其平均相对误差范围在1%~2%,方差在0.001~0.002;Mobley[25]在1999年的研究成果中表明风速与水气界面反射率r呈正比,水面越平静则r越小,测量当天鄱阳湖水面处于较为平静状态(风速在0~5m/s), 取值为0.026是较为精确的;可见光及近红外在轨绝对辐射定标精度在 5%左右[26]。对于大气校正算法本身的特点来看,Sen2cor算法从图像本身获取大气参数,不需要人为输入复杂的参数,所以在这五种大气校正方法中表现最好,但需要占用大量的计算机内存[27]。FLAASH模型输入的参数较6S模型简单,较容易获取校正所需的参数,其校正精度比6S模型表现好[28];6S模型在陆地大气校正中精度较高,但对于水体的校正精度还有待改进[29],且所需参数较多并难获取,导致6S模型精度达不到理想情况;ACOLITE算法一开始针对Landsat系列卫星影像,后期慢慢应用到Sentinel-2 MSI影像,在参数设置上的不同造成校正精度偏低[30],后期需进行进一步的改进以提高校正精度;QUAC算法虽运行速度较快,但默认选取水体像元作为暗像元进行校正,这就不可避免的导致校正精度较差[31]。
本文采用Sen2cor、6S、FLAASH、ACOLITE和QUAC等五种算法对鄱阳湖Sentinel-2B卫星MSI影像进行大气校正处理。并以实测水体遥感反射率光谱曲线为参考,对校正结果进行了验证。结果表明:
1)在Sen2cor、6S、FLAASH、ACOLITE和QUAC五种大气校正算法,Sen2cor整体校正效果最好,平均相对误差为29.37% ,然后是FLAASH算法,平均相对误差为30.11%,6S算法的平均相对误差为48.14%,ACOLITE算法的平均相对误差为52.91%,其中QUAC算法最差,平均相对误差为53.03%。在不考虑简单的情况下,五种算法都能在一定程度上消除大气的影响,其中 Sen2cor算法优于其他四种算法。
2)6S和ACOLITE算法在可见光波段表现相似,在665nm和705nm表现较好,但在其他波段表现较差,FLAASH算法在783nm和842nm的性能优于其他四种算法,Sen2cor算法在705nm和740nm波段表现优于其他四种算法,QUAC算法整体校正效果都较差。根据对比结果,建议根据需要的波段来选择合适的大气校正方法。
我国内陆湖泊较多,水环境也较为复杂,这五种方法虽在一定程度上校正了大气的影响,但由于鄱阳湖的高动态特征以及实测时间与卫星过境时间不完全一致导致校正过程还存在一定误差。为了进一步提高在内陆湖泊的适用性,有必要在校正方法上加强分析和改进,提高校正精度。在后续的研究中,需对鄱阳湖根据浑浊程度以及光学深浅水进行详细分析,并着重分析大气校正参数对其结果的影响,未来将会选择不同传感器卫星、不同季节的数据进行验证分析。