EN中文

地球信息科学学报 >

2020 , Vol. 22 >Issue 10: 2038 - 2050

DOI: https://doi.org/10.12082/dqxxkx.2020.190394

专栏：城乡生态环境综合监测

利用Sentinel-3A OLCI可见光通道反演台湾岛气溶胶光学厚度

  • 王峰 , 1, 2, 3 ,
  • 汪小钦 , 1, 2, 3, * ,
  • 丁宇 1, 2, 3
展开
  • 1.数字中国研究院（福建）,福州 350003
  • 2.福州大学卫星空间信息技术国家地方联合工程研究中心,福州 350118
  • 3.空间数据挖掘与信息共享教育部实验室,福州 350108
*汪小钦(1972— ),女,福建古田人,研究员,主要从事资源遥感环境应用。E-mail:

王峰(1994— ),男,河南开封人,硕士生,主要从事大气遥感研究。E-mail:

收稿日期: 2019-07-24

  要求修回日期: 2020-09-16

  网络出版日期: 2020-12-25

基金资助

国家重点研发计划项目(2017YFB0504203)

中央引导地方发展专项(2017L3012)

版权

版权所有，未经授权，不得转载、摘编本刊文章，不得使用本刊的版式设计。
收起

Retrieval of Aerosol Optical Depth over Taiwan Island Using Visible Channels of Sentinel-3A OLCI

  • WANG Feng , 1, 2, 3 ,
  • WANG Xiaoqin , 1, 2, 3, * ,
  • DING Yu 1, 2, 3
Expand
  • 1. The Academy of Digital China (Fujian), Fuzhou 350003, China
  • 2. National & Local Joint Engineering Research of satellite-spatial Information Technology, Fuzhou University, Fuzhou 350108, China
  • 3. Key Laboratory of Spatial Data Mining & Information Sharing of MOE, Fuzhou 350108, China
*WANG Xiaoqin, E-mail:

Received date: 2019-07-24

  Request revised date: 2020-09-16

  Online published: 2020-12-25

Supported by

National Key Research and Development Program of China(2017YFB0504203)

Central Guide Local Science and Technology Development Projects(2017L3012)

Copyright

Copyright reserved © 2020
Fold

摘要

OLCI（Ocean Land Colour Instrument）作为MERIS（Medium Resolution Imaging Spectrometer）的后继升级版传感器,在气溶胶反演中存在潜在优势,但是目前利用OLCI数据进行气溶胶监测的研究较少。因此,本文针对OLCI多通道反射特征开发了OLCI云检测算法,并对传统查找表构建方法进行改进,根据观测几何特征提出动态查找表法,并通过光谱卷积方式等效转换MODIS和OLCI红蓝通道地表反射率并获取OLCI红蓝通道地表反射率固定关系,进而实现台湾岛550 nm处的AOD反演。与550 nm处AERONET level 2.0 AOD验证结果首先表明不同季节、不同站点的精度表现存在一定差异,其次相对于同期MOD04_3K AOD产品,本文反演结果与全球气溶胶自动观测网（AERONET）站点实测值之间表现出更显著 的相关性（R2=0.8199）,均方根误差（RMSE）从0.175下降到0.113,相对平均误差（RME）从33.6%下降到26.7%,且67.5%的OLCI AOD落在预测误差（EE）区间内,明显大于MOD04_3K AOD落在预测误差区间的百分比（55.7%）。此外,误差分析表明,当实际AOD值较低时,红蓝通道地表反射率之间关系的误判会导致较为明显的AOD反演相对误差。

关键词： OLCI; 气溶胶光学厚度（AOD）; 云识别; 动态查找表; 台湾岛

本文引用格式

王峰 , 汪小钦 , 丁宇 . 利用Sentinel-3A OLCI可见光通道反演台湾岛气溶胶光学厚度[J]. 地球信息科学学报, 2020 , 22(10) : 2038 -2050 . DOI: 10.12082/dqxxkx.2020.190394

Abstract

OLCI (Ocean Land Colour Instrument) sensor, as a successor of MERIS (Medium Resolution Imaging Spectrometer) sensor, has better temporal resolution, spacial resolution and image width, so it has potential advantages in aerosol retrieval, but there are few studies on aerosol monitoring using OLCI data. In order to expand the application field and scope of Sentinel-3AOLCI image data, a high-precision OLCI cloud detection algorithm based on the multi-channel advantage and observation geometry characteristics of OLCI was first proposed, which can effectively and accurately detect thin and thick clouds. Secondly, dynamic lookup tables based on observed geometric features was constructed. Compared with the traditional lookup table, the efficiency of AOD retrieval using dynamic lookup table is improved significantly. The surface reflectance of MODIS and OLCI red and blue channels were also transformed in an equivalent manner through spectral convolution and obtain the fixed relation of surface reflectance of OLCI red and blue channels. Finally, OLCI Aerosol Optical Depth (AOD) retrieval at 550 nm over The Taiwan Island was carried out, and the spatial distribution of AOD we estimated is highly similar to MODIS (Moderate-resolution Imaging Spectroradiometer) AOD. The validation results of OLCI AOD and AERONET (AErosol Robotic Network) level 2.0 AOD at 550 nm show that the accuracy performance of different seasons and stations is different, and then compared with the same period MOD04_3K AOD at 550 nm products, the correlation between the retrieved OLCI AOD and the measured AOD at AERONET stations is more significant (R2=0.8199), the Root Mean Squared Error (RMSE) decreases from 0.175 to 0.113, the Relative Mean Error (RME) decreases from 33.6% to 26.7%, and 67.5% of OLCI AODS fall within Expected Error (EE), which is significantly greater than that of MOD04_3K AOD falling within the prediction error range (55.7%). Error analysis showed that when the actual AOD value is low, the misjudgment of the relationship between surface reflectance of red and blue channels will lead to relatively obvious relative error of AOD retrieval, which means that there are some potential difficulties in retrieving high-precision AOD in areas with good air quality, but for areas with severe pollution, OLCI AOD has significant precision advantage.

Key words： OLCI; Aerosol Optical Depth (AOD); cloud detection; dynamic look up table; Taiwan Island

1 引言

星载气溶胶遥感研究已经进行了近40年[1]。气溶胶是辐射平衡、降水形成、能见度降低以及酸沉降等大气现象形成的关键组分之一,同时能通过吸收散射太阳辐射和地面辐射直接或者间接影响太阳能在地球表面的分布以及地球本身能量的分布,并且改变云的性质和存在时间,进而影响天气、空气质量、大气化学和生物地球化学循环[2,3]。此外,暴露在高浓度的气溶胶环境中—尤其是高浓度细颗粒物（PM2.5）—会对人体健康造成严重危害[4]。因此,准确、快速、大范围地监测气溶胶微物理特性以及时空分布具有重要意义。但是气溶胶的化学成分、来源较为复杂多样,而且在大气中生命周期较短,很难对其直接测量。气溶胶光学厚度（Aerosol Optical Depth, AOD）可以被定义为气溶胶粒子对太阳辐射的大气消光（吸收或散射）在垂直方向上的积分,可用来表示大气浑浊度,是描述气溶胶特性的重要参数之一[5],相较于笼统的气溶胶,AOD更容易获取。获取方式包括地面监测和卫星遥感监测,传统的地面监测方式虽然具有高精度、高频率的优点,但是受限于地形分布和仪器成本等问题,难以做到面状监测,而卫星遥感监测技术则弥补了这些缺点。
基于卫星标量信息监测气溶胶的遥感理论始于20世纪70年代,Griggs等[6]利用搭载在ERTS-1上的MSS（Multispectral Scanner）传感器对海洋上空的AOD进行了反演。随着美国宇航局NASA研制的Terra（1999年）和Aqua（2002年）2颗卫星携带的MODIS（Moderate-resolution Imaging Spectroradiometer）和MISR（Multi-angle Imaging Spectroradiometer）传感器的在轨运行,更是使通过卫星光学信息监测气溶胶技术逐渐迈向业务化阶段。由欧空局（ESA）发射的ENVISAT（Environmental Satellite）卫星携带的MERIS（Medium Resolution Imaging Spectrometer）由于提供较多的通道,可以提供不同应用的信息[7],ESA和部分学者针对MERIS的气溶胶反演开展了很多研究并开发了对应的算法[8]。同样由ESA研制并发射升空的哨兵-3卫星携带了OLCI（Ocean Land Colour Instrument）传感器,该传感器包含了MERIS的所有通道,还增加了新通道,其较高的信噪比、辐射分辨率和双星（A和B）组合小于1 d的时间分辨率特性使其在气溶胶监测方面具有巨大的潜力。Mei等[9]使用XBAER算法从OLCI传感器反演了2016年10月全球尺度的AOD,展示了OLCI在气溶胶反演方面良好的应用潜力,但是构建地表反射率数据库的前提在一定程度上限制了XBAER算法的发展。除此之外,目前利用OLCI进行气溶胶反演的研究几乎没有,极大限制了此传感器在气溶监测方面的应用。
红蓝通道固定地表反射率比值法进行气溶胶反演具有简单、稳定的特点,而且能适用于大多卫星传感器：王中挺等[10]针对GF-1 WFV传感器数据,利用地面观测的植被光谱数据模拟了其红蓝通道的地表反射率比值,对天津地区和北京地区进行了AOD反演,较好地观测了气溶胶分布;郭强等[11]同样利用该方法基于MODIS数据进行了气溶胶反演,并取得了良好地效果。考虑到OLCI的通道范围涵盖红蓝波长而且通道的信噪比较高,因此本文开展利用红蓝通道地表反射率固定关系法进行OLCI气溶胶光学厚度反演实验。
目前云像元的存在将会直接影响地表参数及大气参数定量反演的精度,因此准确且快速的云检测方法将极大的促进卫星资料的实际应用。迄今为止,已经提出了许多云检测算法,主要可分为统计法和阈值法。阈值法利用云与地面物体的光谱辐射特性差异和纹理特征差异设置合适的阈值对二者进行区分,在实际遥感定量产品生产时,阈值法由于计算成本较小、速度较快受到广泛欢迎。OLCI从可见光到近红外多达21个通道的光谱信息更是使利用阈值法进行云检测具有精度保证。
气溶胶反演需要考虑多种大气气溶胶参数,逐像元利用辐射传输模型计算量大,不具有实际意义[13]。因此利用查找表法存储不同模式下大气计算参数,计算气溶胶光学厚度是一套极为可行的方案。目前气溶胶光学厚度反演研究绝大部分是利用查找表法进行查算：陈健等[14]基于查找表法利用GOCI数据对长江三角洲地区进行了AOD反演; 赵志强等[15]同样利用辐射传输模型计算查找表反演HJ CCD影像的气溶胶光学厚度。但是传统的查找表构建方法存在计算成本与精细参数设置无法兼顾的缺点,美国宇航局利用MODIS数据反演中国地区1 km分辨率的大气污染参数需要耗费数个小时甚至更长的时间,因此如何快速准确的构建查找表,是气溶胶光学厚度反演的一个重要研究方向,对于实现气溶胶产品的定量生产具有重要意义。
为拓展OLCI数据在气溶胶光学厚度监测方面的研究,本文以台湾岛为研究区,在针对OLCI特性研发云识别算法、动态查找表构建技术的基础上,利用OLCI红蓝通道地表反射率固定关系式进行了气溶胶光学厚度反演实验。并采用AERONET实测AOD值作为精度验证数据、MOD04_3K AOD值作为交叉验证数据,对比分析本文AOD反演结果。

2 研究区概况、数据来源与研究方法

2.1 研究区概况

台湾岛（120°08'E—122°01'E,21°53'N—25°18'N）,位于中国东海南部,总面积约36 193 km2,70%为山区地形,平原主要分布在西部沿海（图1）。北回归线横穿台湾岛中部地区,致使北部地区表现出明显的亚热带季风气候,而南部地区则主要受热带季风气候支配。
图1 台湾岛海拔及AERONET站点分布

Fig. 1 AERONET Sites distribution and altitude of Taiwan island

在我们的研究时间范围内,台湾岛内可用的地面AERONET监测站点有EPA-NCU、斗六、嘉义、高雄等站（图1）。

2.2 数据来源

（1）Sentinel-3A/OLCI数据
Sentinel-3A卫星于2016年2月16日成功发射,其搭载的OLCI（Ocean and Land Colour Instrument）传感器具有从可见光到近红外的21个通道,空间分辨率为300 m,波段范围从400~1020 nm。OLCI level-1B产品包含的辅助数据集有：经纬度文件、观测几何文件（太阳天顶角、太阳方位角、卫星天顶角、卫星方位角）、大气参量信息文件（水平风速、海平面压力、臭氧总量、相对湿度、参考压力水平、大气温度廓线、水汽总量）,在本文中,臭氧总量和水汽总量信息被用来在辐射传输模型确定大气模式。OLCI数据由欧空局（European Space Agency）免费向公众发布,用户可通过以下网址下载：https://scihub.copernicus.eu/s3/。本文利用台湾岛天气状况较为稳定且云量相对较少的2017年11月到2018年10月的影像数据进行AOD反演实验。由于云污染和地基监测数据缺失等问题存在,首先针对选取时间范围内的影像进行匹配筛选。
（2）MODIS气溶胶产品数据
MODIS的气溶胶产品具有可靠的验证精度和稳定的监测周期（1~2 d）,因此被大量的用于气溶胶特性研究。在本文中,选用与Sentinel-3A过境时间较为接近的Terra 星MOD04_3K AOD（550 nm）产品作为OLCI AOD的交叉对比分析数据,选取影像数据时间范围为2017年11月到2018年10月。可通过美国航天航空局数据下载网站进行下载,数据下载网址如下：https://ladsweb.modaps.eosdis.nasa.gov/search/。
（3）AERONET数据
全球气溶胶自动观测网AERONET（AErosol Robotic Network）是由NASA、LOA-PHOTONS、研究机构及个人联合在全球建立,至今是世界上站点数量最多的气溶胶地基观测网,在全球各地分布超过七百个站点（V3算法AOD数据下载网址https://aeronet.gsfc.nasa.gov/new_web/aerosols.html）。其AOD观测精度为±0.01[16]。AERONET AOD数据可分为3个等级：level 1.0（未筛选）、level 1.5（云筛选）、level 2.0（云筛选与质量保证）。本文选取的AERONET数据时间范围为2017年11月到2018年10月,此外,为了保证验证结果的精度,选用level 2.0级AOD通过二次多项式法插值出550 nm的AOD作为OLCI AOD精度验证数据,相对于Angstrom波长指数法,二次多项式插值法插值得出的550 nm通道的AOD数据精度更高,这主要是因为Angstrom波长指数依赖于荣格尺度分布,而一般情况下气溶胶粒径分布不遵循这一特征[17]。该方法利用3种已知通道（本文利用通道为440、500、675 nm）的AOD值,通过拟合得出公式参数,进而获得550 nm通道的值。
ln τ λ = a 0 + a 1 ln λ + a 2 ln λ 2
式中： τ λ 表示 λ nm处的AOD值; a 0 、 a 1 、 a 2 为相应项未知系数。

2.3 研究方法

2.3.1 AOD反演原理及流程
假设陆地为均匀的朗伯面,不考虑大气的吸收,则卫星传感器观测到的表观反射率（ ρ TOA ）可被表示为：
ρ TOA μ s , μ υ , ϕ = T g ρ o μ s , μ υ , ϕ + T μ s × T μ υ × ρ s μ s , μ υ , ϕ 1 - ρ s μ s , μ υ , ϕ × S
式中： ρ o 为大气的路径辐射项等效反射率; ρ s 表示地表二向反射率; ϕ 为相对方位角; μ s = cos θ s , μ υ = cos θ υ , θ s 与 θ υ 为太阳天顶角和卫星天顶角; S 为大气下界的半球反射率; T μ s 、 T μ υ 分别为上下行散射透过率（大气分子和气溶胶）; T g 是由水汽、臭氧等其他气体引起的总的气体分子透过率（包括上行和下行）[18]。
卫星传感器观测的表观反射率既包含了地表信息,又包含了大气信息。AOD反演的关键就是从卫星信号中将地表反射率信息去除进而通过先验的大气模式及气溶胶类型知识获取气溶胶光学厚度信息。而且,由于较小的地表反射率使得 ρ TOA 对于气溶胶光学特性的变化更加敏感,因此在陆地地表上空气溶胶遥感研究中,选择浓密植被红蓝通道地表反射率进行反演受到广泛应用[19]。此外,红蓝通道地表反射率之间的相关性及稳定性还要显著优于2.1 μm [11],且不受散射角等的影响。因此假设浓密植被区红蓝通道地表反射率关系如下：
ρ s r = k × ρ s b + c
将红蓝通道表观反射率代入式（2）,结合可得以下方程组：
ρ s r μ s , μ υ , ϕ = ρ TOA r μ s , μ υ , ϕ T g - ρ o r μ s , μ υ , ϕ T r μ s × T r μ υ + ρ TOA r μ s , μ υ , ϕ T g - ρ o r μ s , μ υ , ϕ × S ρ s b μ s , μ υ , ϕ = ρ TOA b μ s , μ υ , ϕ T g - ρ o b μ s , μ υ , ϕ T b μ s × T b μ υ + ρ TOA b μ s , μ υ , ϕ T g - ρ o b μ s , μ υ , ϕ × S ρ s r μ s , μ υ , ϕ = k × ρ s b μ s , μ υ , ϕ + c
式中：上标r、b分别代表红、蓝通道;k为浓密植被区红蓝通道地表反射率之间一元线性拟合等式的斜率;c为截距。
将相关先验知识输入6SV2.1辐射传输模型,并通过迭代的方式循环输入AOD值获取相应大气光学参数（ ρ o 、 S 、 T ）,通过上式计算得到相应的红蓝通道地表反射率,然后根据二者先验一元线性拟合方程得到合适的AOD值。
利用OLCI数据进行AOD反演流程示意图如图2所示,主要过程分为OLCI数据输入、数据预处理和AOD反演。
图2 AOD反演流程

Fig. 2 The flow chart of AOD retrieval

（1）OLCI数据输入
本文所使用数据为OLCI level-1.0B级产品的Oa2、Oa3、Oa4、Oa8、Oa17、Oa19辐射亮度数据;太阳天顶角（SZA）、太阳方位角（SAA）、卫星天顶角（OZA）和卫星方位角（OAA）等几何数据文件;以及臭氧总量（TO）、水汽总量（TCWV）等大气参量文件。
（2）数据预处理
利用ESA提供的哨兵应用平台（Sentinel Application Platform, SNAP）对输入的数据进行部分预处理。SNAP包含灵活强大的Sentinel-3处理工具,同时支持大批量数据处理。基于SNAP,对辐射亮度数据进行辐射定标及投影转换,对几何及大气数据进行投影转换。
利用编程语言对辐亮度数据进行数据裁剪、云水像元检测及暗像元挑选;对于几何数据,为了更加灵活高效的构建查找表,需要进行角度匹配。
（3）AOD反演
经预处理之后,分别得到已进行云水掩膜和暗像元挑选的OLCI红蓝通道表观反射率数据（本文中分别指Oa8、Oa4通道表观反射率）、已进行角度匹配的几何数据和已进行单位转换的大气参量数据。将几何数据和大气参量数据输入6SV 2.1辐射传输模型进行计算,得到不同预设AOD值下的大气参数。6SV（Second Simulation of the Satellite Signal）是目前应用较为广泛的辐射传输模型,与它的前一版本6S相比,增加了极化效应,提高了瑞利散射和气溶胶散射的计算精度[20],6SV 2.1是目前公开发行的最新版,该模型可从以下网站下载：http://6s.ltdri.org/pages/downloads.html。
利用计算得到的大气参数基于式（4）前两式计算得到红蓝通道地表反射率,并计算二者地表反射率与式（4）最后一式的距离,获取AOD值,计算公式如下：
D = k × ρ s b μ s , μ υ , ϕ - ρ s r μ s , μ υ , ϕ + c k 2 + - 1 2
在本文中,涉及关键内容有云检测、红蓝通道地表反射率比值确定、动态查找表构建、气溶胶类型确定、结果验证指标及匹配原则。
2.3.2 云、水检测及浓密植被像元挑选
（1）云检测算法
在被动遥感测量中,从卫星影像的太阳光谱区域反演地表和大气的地球物理参数时,充分且准确的云掩膜是必要的预处理步骤,云检测不充分会使AOD反演出现较大的误差[21]。因此本文基于OLCI多通道表观反射率信息和空间变化信息,提出了对应的快速云检测算法。具体流程如图3所示。
图3 云检测流程

Fig. 3 The flow chart of cloud detection

MODIS官方气溶胶反演算法是利用了通道的表观反射率和表观反射率的空间变化信息设置阈值进行云识别,具有相当高的可行性和准确度[22]。本研究在利用这2个特征的基础上,引入归一化植被指数（Normalized Difference Vegetation Index, NDVI）定义动态阈值进行薄云辅助识别。厚云在可见光通道具有较高的反射率,根据厚云在Oa3通道的高反射特性,根据样区统计,设置合理的阈值（0.35）进行初步云检测。Oa19通道处于水汽的吸收通道,且在该通道浓密植被区域的表观反射率高于城区的表观反射率,却低于薄云的表观反射率,根据其反射特性差异,设置阈值进行薄云检测,阈值为对应归一化植被指数（Normalized Difference Vegetation Index, NDVI）大于0.6的所有Oa19像元表观反射率的平均值。
为了进一步增加薄云检测的精度,利用薄云在Oa2通道的非均一性,计算3像元×3像元的标准差作为薄云检测的另一个阈值。此外,为消除云边界的“模糊效应”影响[23],对于被判别为云的像元,其周围5像元×5像元自动扩展为云像元。图4为部分日期影像的云识别结果,可以看出,云识别算法虽然在部分高亮城市地区存在云像元虚检问题,但是依然可以有效准确地检测出薄云与厚云,这为后续AOD反演奠定了坚实的基础。
图4 台湾岛OLCI云识别结果和真彩色影响对比

Fig. 4 Comparison of cloud recognition results and true color effects OLCI Taiwan island

（2）水像元检测及浓密植被像元挑选
用归一化水体指数（Normalized Difference Vegetation Index, NDWI）检测水像元[24],所用通道为Oa6、Oa17。当NDWI大于0时,待判断像元被标记为水像元。
归一化植被指数与植被浓密度具有较好的相关性,常用作浓密植被筛选指标[25]。为保证足够多的像元参与AOD反演,设置当像元NDVI大于0.45时为浓密植被像元。
2.3.3 浓密植被区红蓝通道地表反射率比值确定
在NASA的V 5.2气溶胶算法中,对于MODIS数据,浓密植被区红蓝通道地表反射率之间的关系可表示为：
ρ S M 3 = 0.49 × ρ S M 1 + 0.005
M3代表MODIS数据第3通道（蓝光通道）,M1代表MODIS数据第1通道（红光通道）。
OLCI数据与MODIS数据在成像方式、光谱响应函数、通道范围方面均存在较大差异,因此二者对相同目标地物观测得到的反射特征也不尽相同。因此将OLCI和MODIS相应红蓝通道的光谱响应函数与1 nm重采样后的USGS中实测植被（500多种）光谱进行卷积,从而得到OLCI和MODIS相应红蓝通道实测目标植被的地表反射率。通过一元回归拟合建立二者对应通道之间的关系（图5）,二者对应红蓝通道之间具有较好的相关性,相关系数均高于0.99。
图5 MODIS与OLCI红蓝通道地表反射率散点图

Fig. 5 Scatter plot of surface reflectance of red and blue channels for MODIS and OLCI

根据拟合得到的关系式,结合式（6）,从而获得OLCI红蓝通道地表反射率关系模型：
ρ S Oa 4 = 0.497 × ρ S Oa 8 + 0.008
2.3.4 动态查找表构建
传统查找表构建方式受到参数步长、传感器 观测几何条件等因素的影响,不具有很好的移植性,往往需要耗费大量的时间计算参数组合。而且由于参数步长设置较大,很容易导致结果AOD影像存在分块现象。针对上述存在的问题,利用动态构建查找表法逐影像计算获取大气参数数据并将其应用于该影像的AOD获取,该算法可成倍提高AOD查算效率。动态查找表构建算法示意图如图6所示。
图6 动态查找表构建示意

Fig. 6 Construction chart of dynamic look-up table

动态查找表构建算法核心步骤描述如下：
（1）几何角度数据文件预处理。首先将待计算OLCI影像对应的几何角度数据文件的像元值取整,然后统计记录取整之后每个几何角度数据文件中不重复的角度值数量（角度整值数量）;
（2）影像计算区域划分。首先按照角度整值数量对几何角度数据文件进行从大到小排序,角度整值数量最多的几何角度文件作为第一基准匹配文件,根据该文件的角度整值信息逐角度划分影像区域（A1）;数量排名第二的几何角度文件作为第二基准匹配文件,根据该文件的角度整值信息逐角度划分影像区域（A2）;然后对A1和A2区域求并集获得最终的影像划分区域A3;
（3）逐区域计算大气参数。根据步骤（2）得到的影像计算划分区域A3,逐区域对每种几何角度数据文件求平均作为该区域的几何角度数据输入,并几何预设AOD值（550 nm处AOD步长设置范围为0~0.9之间,步长为0.05）以及气溶胶模型（气溶胶模型的选择对于AOD反演精度有重要影响,王毅等[26]根据有效的太阳辐射计数据进行分析,发现中国东部沿海地区冬季较为符合6SV2.1模型内置的大陆型气溶胶模型;贾亮亮等[27]利用GF-1 卫星数据,将不同模型反演出的AOD与和地基监测AOD进行比较,同样发现将大陆型气溶胶模型作为输入精度较高。因此,本文选取大陆型气溶胶模型作为输入条件。）输入6SV2.1辐射传输模型得到大气参数数据：大气下界半球反射率（ S ）、路径辐射项等效反射率（ ρ o ）、总散射透过率（ T μ s × T μ υ ）、总气体分子透过率（ T g ）。
由构建过程可以看出,针对不同的影像（几何角度信息不同,甚至光谱响应函数都会发生变化）,动态查找表算法可以自动划分不同的计算区域并得到相应区域的大气参数,相对于传统的查找表法,适用性更强,而且不需要提前计算海量大步长参数组合对应的大气参数。
2.3.5 评价指标
本文采用影像AOD预测值与地面观测值之间的决定系数（R2）、均方根误差（Root Mean Squared Error, RMSE）、预测误差区间（Expected Error, EE）[28]、相对平均误差（Relative Mean Error, RME）,以及一元线性回归方程多种统计指标来对反演结果进行评价。其中EE、RMSE和RME的计算公式如下：
RMSE = ∑ i = 1 n τ ˆ - τ 2 n
EE = ± 0.05 ± 0.2 × τ
RME = ∑ i = 1 n τ ˆ - τ n × τ ̅ × 100 %
式中： τ ˆ 为AOD预测值; τ 为地基AOD值; τ ̅ 为地基AOD平均值;n为验证的点对数量。

3 结果及分析

为了全面地验证本文反演的OLCI AOD的精度及空间分布合理性,首先利用地基监测AOD值进行逐站点分析,然后与MOD04_3K AOD产品进行影像对比和AOD数值对比分析。

3.1 利用AERONET实测数据进行验证

为保证结果的合理性,AERONET level 2.0 AOD值取相应日期影像拍摄时刻前后半小时 （±0.5 h）数据的平均值。反演结果取地面观测站点对应像元AOD值,如该像元无值,则取其周围3 km范围内均值。为避免云像元的“模糊效应”对反演结果造成影响,如果地基站点周围3 km范围内存在云像元,则该匹配数据对不参与精度验证,可参与匹配的影像共计61景,有效匹配数据对83对。
图7为OLCIAOD与AERONET AOD的逐站 点对比结果,可以看出,在4个站点中,二者具有 相当一致的变化趋势。但是在不同的站点反演的结果和实际观测值分布趋势也存在微弱的差别, 其中EPA-NCU站点（图7（a））和Kaohsiung站点 （图7（d））OLCIAOD表现出了微弱的高估趋势,且当AOD值较低时这种趋势更为明显。对于Chiayi站点（图7（b））,除2018年3月31日、4月1日和4月5日部分日期的反演结果质量相对较差外,其余反演结果和地基AOD结果一致性较好。Douliu站由于在2018年1月15日之后便不再提供AERONET AOD数据,因此可参与匹配的站点数据极少,但是在有限的数据对比中可以看出（图7（c））,除2018年2月13日反演结果出现较大偏差外,其余日期反演的AOD值均与站点值较为接近。
图7 卫星OLCIAOD与地基监测AOD逐站点对比

Fig. 7 Comparison between retrieved AOD and in-situ AOD at different sites

表1为OLCI AOD在不同季节、不同站点的精度表现（横线代表无统计数据）。可以看出,对于全年数据而言,Chiayi站点和Kaohsiung站点的表现最好,前者的RMSE和RME分别低至0.111和19.0%,后者的RMSE和RME分别低至0.101和32.9%;EPA-NCU站点的精度表现稍次之,RMSE和RME分别为0.120和33.2%;Douliu站点则出现了较大的误差,RMSE和RME分别高达0.133和58.1%,很可能是由于匹配数据对较少（4对）,而且该站点匹配数据对集中于冬季,在一定程度上缺乏精度表现的评价意义。对于季节匹配数据,各站点表现出了明显不同的季节分异,其中EPA-NCU站点精度随着季节的更替发生了明显的变化,主要表现为在春夏秋冬季节误差逐渐上升的趋势,尤以冬季的表现最差,RMSE和RME分别高达0.153和59.5%;Chiayi站点则在冬季表现最优,RMSE和RME分别低至0.074和14.2%,春秋季表现稍差一些;Kaohsiung站点在春夏秋3个季节的精度表现极为接近,RMSE和RME分别稳定于0.1和33%左右,冬季表现稍差一些,RMSE和RME分别为0.137和29.8%。
表1 不同季节各站点OLCI AOD和地基监测AOD精度比较

Tab. 1 Accuracy comparison between retrieved AOD and in-situ AOD in all seasons

站点 指标 春 夏 秋 冬 全年
EPA-NCU RMSE
RME/%		 0.095
27.4		 0.120
28.9		 0.128
33.8		 0.153
59.5		 0.120
33.2
Chiayi RMSE
RME/%		 0.135
19.8		 -
-		 0.120
31.4		 0.074
14.2		 0.111
19.0
Douliu RMSE
RME/%		 -
-		 -
-		 -
-		 0.133
58.1		 0.133
58.1
Kaohsiung RMSE
RME/%		 0.101
33.9		 0.091
33.4		 0.108
33.6		 0.137
29.8		 0.101
32.9

3.2 与MODIS AOD产品对比

将本文OLCI AOD与MOD04_3K AOD 产品进行对比交叉验证,二者部分日期（2018年1月24日和2018年3月1日）空间分布对比结果如图7所示,可以看出,本文反演结果与MOD04_3K AOD在空间分布上高度一致,不仅高低值相互对应,且相同位置的AOD值也较为相近,不同的是,本文的反演结果更为精细,且空间覆盖有效范围更大。由图8可看出,台湾岛西部沿海地区以及南部高雄等地AOD值明显高于中部五大山脉区,可知人类活动对于大气质量变化具有显著影响。
图8 OLCI AOD与MOD04_3K AOD空间分布对比

Fig. 8 Spatial distribution Comparison between OLCI AOD and MOD04_3K AOD

进一步将OLCI AOD和MOD04_3K AOD分别与经时空匹配的地基监测AOD进行散点回归分析（图9）。图中红色虚线为预测误差区间（EE）。从图9（a）可看出,OLCI AOD与实测AOD值较为接近,决定系数R2达到了0.8199,均方根误差RMSE和相对平均误差RME分别低至0.113和26.7%;OLCI AOD与实测AOD一元线性拟合方程的斜率为0.8538,截距为0.1235,表明OLCI AOD出现微弱的低值高估和高值低估情况;此外,67.5%的OLCI AOD落在预测误差区间内（等于EE）,33.5%的OLCI AOD落在预测误差区间上方（大于EE）。图9（b）为MOD04_3K AOD产品与地基实测AOD值的同期精度图,从结果可以看出,二者决定系数R2为0.7049,RMSE和RME分别为0.175和33.6%;根据其一元线性拟合方程和结果图观察,其低值高估和高值低估趋势得到进一步放大;此外,55.7%的MOD04_3K AOD落在预测误差区间内（等于EE）,37.2%落在预测误差 区间上方（大于EE）,7.1%落在预测误差区间下方（小于EE）。通过对比分析可以发现,本文OLCI AOD在精度评价指标上均优于同期MOD043K_AOD产品,具有较高的精度和稳定性。
图9 OLCI AOD与MOD04_3K AOD同期精度验证对比

Fig. 9 Accuracy comparison between retrieved AOD and MOD04_3K AOD

3.3 误差分析及讨论

本文反演AOD的前提是：在浓密植被区域,红蓝通道地表反射率表现出固定的关系（式（3））,准确的斜率k和截距c是实现高质量AOD反演的决定性因素。因此,有必要对斜率或者截距改变带来的误差做进一步分析讨论。
根据第2.2.2节的计算结果,本文将OLCI红蓝通道（分别为Oa4和Oa8波段）的地表反射率之间的斜率设为0.497,截距设为0.008。将这2个参数值作为理论真实值,假设红通道的地表反射率为0.1,则蓝通道相应地表反射率为0.0577。在模拟大气状况下（AOD:0.1~0.9,步长为0.1;大气模式：中纬度夏季;气溶胶模型：大陆型）,改变红蓝通道地表反射率的斜率或者截距（斜率：0.470~0.520,步长为0.001;截距：0.006~0.010,步长0.0001）,然后利用6SV2.1模型进行辐射传输运算,获取斜率或者截距变化之后的AOD计算模拟值,并计算与模拟真实值的相对差异（RD,%）和绝对差异（AD）,二者的计算公式如下：
RD = β ̅ - β β × 100 %
AD = β ̅ - β
β 为AOD模拟真实值, β ̅ 为AOD计算模拟值,计算时几何角度条件设置：太阳天顶角为41°、太阳方位角为128°,卫星天顶角为19°,卫星方位角为102°。
由图10可以看出,斜率改变0.01或者截距改变0.001会带来最少0.01的AOD反演误差,而且斜率或者截距正向改变和负向改变带来的AOD绝对误差和相对误差基本成线性对称分布,随着AOD的增加,绝对误差的绝对值和相对误差的绝对值均有不同程度的增大,前者增加的幅度远小于后者,表明在AOD低值区间斜率和截距误差将会带来较小的绝对误差绝对值和较大的AOD反演相对误差绝对值,如图10所示,当AOD低于0.2时,斜率或者截距的误差会造成超过13%的相对误差绝对值。
图10 红蓝通道地表反射率公式改变对AOD反演结果的影响

Fig. 10 Effect of surface reflectance equation of red and blue channels on AOD retrieval results

在相同的大气状况和几何角度条件下,同样模拟了采用NASA V5.2气溶胶算法中MODIS浓密植被区红蓝通道地表反射率公式（式（6））和本文采用的OLCI浓密植被区红蓝通道地表反射率公式 （式（7））进行AOD反演带来的差异,如图11所示。可以看出,由于斜率和截距的同时改变,带来的绝对差异和相对差异要远大于斜率或者截距单独改变带来的绝对差异和相对差异。当AOD比较小时,差异非常大,随着AOD的增加,二者差异快速变小。当AOD小于0.4时,相对差异远大于10%,且变化较为剧烈;当AOD处于0.4~0.9时,相对差异变化较为平缓,逐渐稳定在4%左右。
图11 利用MODIS和OLCI红蓝通道地表反射率公式进行AOD反演带来的差异

Fig. 11 Difference of AOD results between using MODIS and OLCI surface reflectance equation of red and blue channels

4 结论

基于浓密植被区红蓝通道具有稳定的地表反射率关系,利用Sentinel-3A星的OLCI传感器对台湾岛陆地气溶胶进行了气溶胶反演。获得结论如下：
（1）构建了针对OLCI传感器数据的快速云识别算法,能有效优化薄云识别区域,一定程度上避免了云像元对气溶胶反演结果的影响。
（2）利用植被地表在OLCI红蓝波段的地表反射率之间的固定关系能有效较好地反演陆地气溶胶光学厚度,且误差分析表明,不准确的斜率和截距会造成AOD反演误差,尤其是当二者同时出现误差时,AOD反演精度会进一步大幅下降。
（3）传统的查找表构建方法耗时较长且步长较大,根据每幅影像几何角度特点划分区域构建的动态查找表,能有效提高AOD反演效率。
（4）精度验证结果表明不同站点的全年精度表现以及季节表现存在一定分异,此外,将反演结果与MOD04_3K气溶胶产品进行了比较,发现不仅二者空间分布较为一致,而且本文反演结果与地基实测值的精度验证指标均优于同期MOD04_3K气溶胶产品,结果表明,本文反演结果和地面实测值的决定系数R2达到了0.8199,RMSE和RME分别低至0.113和26.7%,具有较高质量。
本研究也存在一些不足之处。气溶胶模型只考虑大陆型,忽略了气溶胶组分变化的多样性与复杂性。浓密植被区地表反射率关系设为固定值,还需根据更多观测数据分季节分区域进行进一步研究,而且没有考虑BRDF的影响。

参考文献

原文顺序 | 文献年度倒序 | 文中引用次数倒序
[1]
蒋哲, 陈良富, 王中挺. 细粒子气溶胶光学厚度和谱分布偏振的反演[J]. 地球信息科学学报, 2012,14(4):460-464.

[ Jiang Z, Chen L F, Wang Z T. Study on the retrieval of fine mode aerosol optical depth and size distribution using polarized signal[J]. Journal of Geo-information science, 2012,14(4):460-464. ]

[2]
Chin M, Diehl T, Tan Q, et al. Multi-decadal aerosol variations from 1980 to 2009: A perspective from observations and a global model[J]. Atmospheric Chemistry and Physics, 2014,14(7):3657-3690.

[3]
Safai P D, Raju M P, Budhavant K B, et al. Long term studies on characteristics of black carbon aerosols over a tropical urban station Pune, India[J]. Atmospheric Research, 2013,132-133:173-184.

[4]
Jung C R, Hwang B F, Chen W T. Incorporating long-term satellite-based aerosol optical depth, localized land use data, and meteorological variables to estimate ground-level PM2.5 concentrations in Taiwan from 2005 to 2015[J]. Environ Pollut., 2018,237:1000-1010.

[5]
Diner D J, Beckert J C, Reilly T H, et al. Multi-angle Imaging SpectroRadiometer (MISR) instrument description and experiment overview[J]. IEEE Transactions on Geoscience and Remote Sensing, 1998,36(4):1072-1087.

[6]
Griggs M. Measurements of atmospheric aerosol optical thickness over water using ERTS-1 data[J]. Air Repair. 1975,25(6):622-626.

[7]
Verstraete M, Bernardpinty, Curran P. MERIS potential for land applications[J]. International Journal of Remote Sensing, 1999,20(9):1747-1756.

[8]
Mei L, Rozanov V, Vountas M, et al. Retrieval of aerosol optical properties using MERIS observations: Algorithm and some first results[J]. Remote Sensing Environment, 2017,197:125-140.

[9]
Mei L, Rozanov V, Vountas M, et al. XBAER-derived aerosol optical thickness from OLCI/Sentinel-3 observation[J]. Atmospheric Chemistry and Physics, 2018,18(4):2511-2523.

[10]
王中挺, 辛金元, 贾松林, 等. 利用暗目标法从高分一号卫星16m相机数据反演气溶胶光学厚度[J]. 遥感学报. 2015,19(3):530-538.

[ Wang Z T, Xin J Y, Jia S L, et al. Retrieval of AOD from GF-1 16m camera via DDV algorihm[J]. Journal of Remote Sensing, 2015,19(3):530-538. ]

[11]
郭强, 唐家奎, 何文通, 等. 利用MODIS可见光波段反演陆地气溶胶光学厚度[J]. 地理与地理信息科学, 2015,31(2):38-43.

[ Guo Q, Tang J K, He W T, et al. Aerosol optical depth retrieval over land using MODIS visible bands imagery[J]. Geography and Geo-information Science, 2015,31(2):38-43. ]

[12]
刘向培, 王毅, 石汉青, 等. 基于统计特征的中国东南沿海区域云检测[J]. 中国图象图形学报, 2010,15(12):1783-1789.

[ Liu X P, Wang Y, Shi H Q, et al. Cloud etection over the southeast China basing on statistical analysis[J]. Journal of Image and Graphics, 2010,15(12):1783-1789. ]

[13]
郑文武, 曾永年. 利用MODIS数据和查找表的城区TM影像大气校正[J]. 测绘科学, 2012,37(3):27-29.

[ Zheng W W, Zeng Y N. An atmospheric correction algorithm for urban area TM data based on MODIS data and look-up table[J]. Science of Surveying And Mapping, 2012,37(3):27-29. ]

[14]
陈健, 周杰, 李雅雯. 基于静止卫星GOCI数据的陆地上空气溶胶光学厚度遥感反演[J]. 遥感技术与应用. 2017,32(6):1040-1047.

[ Chen J, Zhou J, Li Y W. Retrieving aerosol optical depth over land based on GOCI data onboard geostationary satellite[J]. Remote Sensing Technology and Application, 2017,32(6):1040-1047. ]

[15]
赵志强, 李爱农, 边金虎, 等. 基于改进暗目标法山区HJ CCD影像气溶胶光学厚度反演[J]. 光谱学与光谱分析, 2015,35(6):1479-1487.

[ Zhao Z Q, Li A N, Bian J H, et al. An improved DDV method to retrieve AOT for HJ CCD image in typical mountains areas[J]. Spectroscopy and Spectral Analysis, 2015,35(6):1479-1487. ]

[16]
Holben B N, Tanre D, Smirnov A, et al. An emerging ground-based aerosol climatology: Aerosol optical depth from AERONET[J]. Journal of Geophysical Research Atmospheres, 2001,106(D11):12067-12097.

[17]
李龙, 施润和, 张璐, 等. 华东地区MODIS与OMI气溶胶光学厚度数据融合[J]. 地球信息科学学报, 2015,17(10):1224-1233.

[ Li L, Shi R H, Zhang L, et al. Data fusion of MODIS AOD and OMI AOD over east China using universal Kriging[J]. Journal of Geo-information Science, 2015,17(10):1224-1233. ]

[18]
Vermote E F, Saleous M N. Operational atmospheric correction of MODIS visible to middle infrared land surface data in the case of an infinite Lambertian target[J]. Earth Science Satellite Remote Sensing, 2006: 123-153.

[19]
田信鹏, 孙林, 刘强, 等. 北京地区Landsat 8 OLI高空间分辨率气溶胶光学厚度反演[J]. 遥感学报, 2018,22(1):51-63.

[ Tian X P, Sun L, Liu Q, et al. Retrieval of high-resolution aerosol optical depth using Landsat 8 OLI data over Beijing[J]. Journal of Remote Sensing, 2018,22(1):51-63. ]

[20]
Tanre D, Deuze J L, Herman M, et al. Second simulation of the satellite signal in the solar spectrum-6s code[J]. IEEE Transactions on Geoscience & Remote Sensing. 2002,35(3):675-686.

[21]
Yang W, Chen L, Shang H. A new cloud mask algorithm used in aerosol retrieval over land for Suo-NPP VIIRS[C]//. IGARSS 2016-2016 IEEE International Geosicence And Remote Sensing Symposium, 2016.

[22]
Remen A L, Chambless D L, Steketee G, et al. Second-generation operational algorithm: Retrieval of aerosol properties over land from inversion of Moderate Resolution Imaging Spectroradiometer spectral reflectance[J]. Journal of Geophysical Research Atmospheres, 2007,112(D13).

[23]
Mei L, Vountas M, Gómez-Chova L, et al. A Cloud masking algorithm for the XBAER aerosol retrieval using MERIS data[J]. Remote Sensing of Environment, 2016,197:141-160.

[24]
Mcfeeters S K. The use of the Normalized Difference Water Index (NDWI) in the delineation of open water features[J]. International Journal of Remote Sensing, 1996,17(7):1425-1432.

[25]
Liang S L, Fallahadl H, Kalluri S, et al. An operational atmospheric correction algorithm for Landsat Thematic Mapper imagery over the land[J]. Journal of Geophysical Research Atmospheres, 1997,102(D14):17173-17186.

[26]
王毅, 石汉青, 黄思训. 中国东南地区及近海海域气溶胶反演遥感研究[J]. 遥感技术与应用, 2009,24(1):13-21.

[ Wang Y, Shi H Q, Huang S X. The retrieval of atmospheric aerosol optical depth in southeast China[J]. Remote Sensing Technology and Application, 2009,24(1):13-21. ]

[27]
贾亮亮, 汪小钦, 苏华, 等. 台湾岛高分一号卫星WFV数据气溶胶反演与验证[J]. 环境科学学报, 2018,38(3):1117-1127.

[ Jia L L, Wang X Q, Su H, et al. Validation of retrieving aerosol over Taiwan Island using GF-1 satellite WFV data[J]. Acta Scientiae Circumstantiae, 2018,38(3):1117-1127. ]

[28]
Levy R C, Mattoo S, Munchak L A, et al. The Collection 6 MODIS aerosol products over land and ocean[J]. Atmospheric Measurement Techniques, 2013,6(11):2989-3034.

Options
文章导航

/

〈 〉