大气气溶胶是悬浮在大气中的固体和液体粒子的总称，其主要来源于自然界产生的矿物质气溶胶和人为排放的细粒子气溶胶，大气气溶胶对地球辐射平衡、气候变化、大气环境、人类健康等有重要影响^[1]。气溶胶光学厚度（Aerosol Optical Depth, AOD）定义为气溶胶的消光系数在垂直方向上的积分，是气溶胶最基本的光学特性，也是评估气候变化中气溶胶辐射强迫的关键因子。地面站点太阳光度计可获取详细的气溶胶光学参数，但由于大气气溶胶的化学组成、尺度分布、生命周期等的时空多变性，地面站点不足以准确观测和评估大范围的大气气溶胶的分布特征。卫星遥感反演气溶胶产品的出现，为认识和研究气溶胶的气候效应和区域大气污染监测提供了强有力的数据支撑。

AOD反演的关键在于如何将大气信号从卫星接收到的表观信号中分离出来，Kaufman等^[2-3]研究发现在浓密植被等反射率较低的暗目标区，受大气影响较弱的近红外波段与红、蓝波段的反射率有很好的线性关系，据此关系可实现暗目标区域的地气解耦并能很好地反演得到较准确的气溶胶光学厚度，暗像元法被广泛应用于各类卫星的气溶胶光学厚度反演，其中MODIS利用2.1μm波段反射率推算可见光波段地表反射率从而实现AOD反演，且已经实现业务化运行并对外发布产品^[4]。王宏斌等^[5]对中国地区MODIS气溶胶产品进行精度验证，结果表明MODIS产品在我国东部地区的精度较高，同时还对误差来源进行了分析，得出MODIS误差来源主要是地表反射率关系、气溶胶模型的不适应。基于暗像元法的气溶胶反演算法近年来在不断改进^[6-7]，但其基本思想并没有改变，暗像元法依然存在两个主要问题：①暗像元法的波段比值关系在城市、裸地等高地表反射率的亮地表区域不适用;②由于不同时段、不同地物的光谱反射率存在差异，短波红外、红、蓝通道的地表反射率关系并不是固定的。为了弥补暗像元法在亮地表区域的AOD反演问题，Hsu等^[8-10]提出的深蓝算法可以应用于城市、沙漠和裸地等区域的AOD反演，但其反演精度低于暗像元法。Huang等^[6]结合最小反射率法和双星协同反演算法对MODIS进行AOD反演，在暗像元区域反演结果与MODIS AOD结果接近，但在亮地表区域的反演能力有很大提升。

极轨卫星的低轨道高度和全球观测的能力在全球尺度的气溶胶研究中有较大优势。但极轨卫星较长的重访周期对于持续监测区域性大气污染变化显得无能为力，且极轨卫星对不同地区反演AOD的时间不同，导致同一天不同地点的AOD缺乏可比性。大气气溶胶的分布可能与边界层日变化，天气系统影响等因素有关，大气气溶胶含量和分布在一天内变化较大，故需要高时空分辨率的观测才能较准确地获得其变化趋势及较细致的分布特征，对于同一区域的观测而言，静止轨道气象卫星的重访周期相对于极轨卫星来说有非常大的优势，但受限于静止轨道卫星光谱波段过少、光谱范围较宽、空间分辨率和定标精度较低等因素，以往的研究中静止轨道卫星一般用于海面或湖面等的气溶胶反演，在陆地上反演效果较差^[11-13]。

Himawari-8搭载的AHI传感器可见光波段达到0.5~1 km,中短波红外达到2 km。Yu等^[14]将AHI和极轨卫星NPP VIIRS数据进行交叉对比，发现两者在用于反演AOD的波段上吻合度较高，说明AHI传感器拥有较高的空间分辨率和定标精度。由此可见，Himawari-8卫星在区域大气污染监测领域有很大的应用潜力。Yoshida等^[15]基于最小表观反射率和暗像元法对Himawari-8进行了气溶胶光学厚度反演，同时利用韩国的2个AERONET站点（Baeksa、Taehwa）的观测资料对反演结果进行验证，反演结果（R=0.56，RMSE=0.174）表明该算法仍有待改进，但目前Himawari-8对外发布的AOD产品仍采用该算法。本文基于Himawari-8 AHI卫星数据和MOD09A1地表反射率产品，利用改进的红蓝波段比值反演方法，反演了中国区域2018年4月~2019年4月逐10分钟的AOD。利用中国地区西北、华北、华东、华南等6个AERONET地基观测资料对AOD反演结果进行精度验证，同时将反演结果与Himawari-8官方发布的AOD产品以及MODIS AOD产品进行时间序列、空间分布一致性对比，检验反演结果的可靠性。以期为了解中国地区大气气溶胶的时空变化规律和全天时的大气污染监测方法提供新的参考。

AERONET（Aerosol Robotic Network）是由NASA组织发起的全球气溶胶监测网络^[32]，AERONET网站（https://aeronet.gsfc.nasa.gov）提供了全球范围内多个站点的气溶胶光学厚度、水汽等数据。其发布的数据有3个质量等级：level 1.0为原始数据，level 1.5为云掩膜处理后的数据，level 2.0经过了云掩膜和质量控制，地基观测资料一般每小时有4~5次观测。由于中国地区提供level 2.0数据的站点较少，因此采用了经过云掩膜处理的level 1.5级数据。表1列出了用于AOD产品验证的站点的详细信息，站点主要分布在西北、华北、华东、华南等地。

AERONET AOD观测误差小于0.02^[16]，可以用于验证和评估卫星反演AOD的精度。但由于AERONET缺少550nm波段的观测，无法与550nm AOD反演结果进行直接对比，因此需要利用AERONET其他观测波段的AOD值推算550nm处的AOD值。有研究表明，气溶胶光学厚度的对数与波长的对数之间有非常好的二次多项式拟合关系，拟合误差小于0.02^[17]，本文利用AERONET 440、500、675和870nm 4个波段观测的AOD求得拟合方程进而推算出550nm处的AOD值。由于Himawari-8反演的AOD为500 nm波处的AOD，因此可直接与AERONET 500nm AOD进行对比。

Himawari-8号卫星发射于2014年10月7日，其搭载的AHI（Advanced Himawari Imager）传感器相比上一代静止气象卫星有了很大的改进^[18]，其时间分辨率达到了分钟级，全圆盘观测仅需10 min，区域加密观测时间分辨率达到2.5 min。空间分辨率方面，其中红波段为500 m×500 m，其他可见光波段为1 km×1 km，近红外和红外波段为2 km×2 km。Himawari-8观测数据和产品由日本宇宙航空研究开发机构（JAXA）发布（https://www.eorc.jaxa.jp/ptree）。

Himawari-8 L1圆盘数据时间分辨率为10 min，空间分辨率已重采样至5 km×5 km。L1数据存储格式为NetCDF格式，其中包含了可见光和近红外波段的反照率、中红外波段的亮度温度、卫星天顶角、卫星方位角、太阳天顶角、太阳方位角、观测时间等数据。

Himawari-8 CLP产品时间分辨率为10 min，空间分辨率为5 km×5 km。数据集包含有云光学厚度、云顶温度、云顶高度、云属性（ISCCP Definition）等产品。在反演AOD过程中，本文利用CLP产品中的云属性进行云掩膜。

Himawari-8 AOD产品分为L2级和L3级，L2级AOD产品实际分辨率为10 min，L3级AOD产品有1小时、1天和1个月3种时间分辨率，L2和L3级AOD产品空间分辨率均为5 km×5 km。数据集包含有Angstrom指数、500 nm处气溶胶光学厚度、QA质量文件等数据。本文采用L2级产品与反演结果进行对比。

地表反射率比值库基于MOD09A1地表反射率产品，MOD09A1产品提供经过大气校正后的MODIS Band1-Band7的地表反射率产品以及QA质量文件，产品空间分辨率500 米，时间分辨率为8天。MOD09A1 8天合成地表反射率产品考虑了高观测覆盖、低视角、无云及云的阴影以及气溶胶浓度等多种因素，选取8天内最合适的地表反射率。本研究基于Google Earth Engine平台对MOD09产品进行月最小值合成、重采样等处理，获得每个月的红蓝波段的最小地表反射率的比值，空间分辨率重采样至5 km×5 km。

MOD/MYD04_3K是NASA发布的Level 2级第6版气溶胶产品^[4]，其空间分辨率为3 km，其较第5版的10 km分辨率有较大提升，MOD04_3K原始数据储存格式为HDF-4，其中包含了气溶胶光学厚度、Angstrom指数、质量控制等51个波段信息。MODIS AOD产品是目前应用最为广泛的AOD数据之一，其AOD反演精度在浓密植被区较高。

卫星传感器观测的辐射值主要来自于地表贡献和大气贡献，在假设地表为朗伯体且大气均一的前提下，卫星观测到的表观反射率可表示为：

其中：μs与μv分别为太阳天顶角θs的余弦值、观测天顶角θv的余弦值，φ为相对方位角，T表示大气透过率，ρs为地表反射率，ρ0为路径程辐射项，S为大气下界的半球反照率。

AOD反演的关键在于如何精确计算出地表反射率，从而去除地表贡献获取大气贡献进而反演得到AOD。在反演过程中为加快反演速度，一般采用查找表方式进行反演。在构造查找表时，本文中气溶胶模式被假定为普适性较强的大陆型气溶胶模式，借助6S大气辐射传输模型构造不同卫星观测条件、不同波段、不同大气模式、不同AOD条件下的大气参数查找表。传统的基于地表反射率库反演AOD的方法^[19-20]，是假设短时间内地表反射率没有发生太大改变，将MOD09A1地表反射率代表真实地表反射率，利用生成的查找表和Himawari-8表观反射率，迭代计算不同AOD条件下蓝波段地表反射率ρ_blue，当ρ_blue与MODIS对应蓝波段地表反射率最为接近时，当前迭代的AOD值即为反演的AOD值(图1)。

反演过程中，先对Himawari-8 L1数据集进行云、冰、水、雪掩膜，读取筛选后的晴空像元的表观反射率和太阳天顶角、卫星天顶角等角度信息。将MOD09A1地表反射率产品重采样至和Himawari-8 L1数据一致的空间分辨率，读取MOD09A1红蓝地表反射率并合成红蓝比值。初始AOD条件设置为0.001，AOD增加步长为0.05，将Himawari-8表观反射率和角度信息以及AOD条件代入式（1）计算得到ρ₁和ρ₃，由于卫星接收到的红蓝表观反射率会随着AOD的增大而增大，波长越短的波段对气溶胶的增加越敏感，因此随着AOD的增加，ρ1/ρ3会逐渐减小，在迭代过程中当ρ₁/ρ₃-R<0时停止计算，此时迭代的AOD和计算得到的ρ₁和ρ₃的比值分别记为aod⁺和rb⁺，前一个迭代AOD值和ρ₁和ρ₃的比值分别记为aod^-和rb^-，根据权重公式（5）插值得到最终AOD结果。

基于红蓝地表反射率比值库的反演方法反演了中国地区2018年4月~2019年4月逐10 min的气溶胶光学厚度，共计19 840景数据。将验证站点所处位置对应的AERONET AOD（地基站点实测AOD）、RB AOD（红蓝比值法反演的AOD）、H8AOD（JAXA发布的Himawari-8 AOD产品）进行提取和匹配，在匹配过程中要求三者观测时间间隔小于5min且三者均有观测值或反演值。

图2为6个AERONET站点的地基观测AOD与RB AOD和H8 AOD的线性拟合散点图。由于卫星可见光波段受云的影响较大，因此北方站点较南方站点匹配的样本点多。

对比RB AOD与AERONER AOD结果（图2(a)、(c)、(e)、(g)、(i)、(k））表明：包头站（R=0.11、RMSE=0.57）的反演效果较差，其余5个站点反演效果较好。除香河站外，RB AOD在其余站点都有不同程度的低估，尤其是南方站点香港、台北站低估较为严重，分别有51%、31.7%的样本点被低估；北京、徐州、低估样本点数量占比分别为22%和22.1%。北京、徐州、香河站拟合直线斜率接近于1，而台北、香港站拟合直线斜率分别为0.643和0.463。另一方面，对比H8 AOD与AERONER AOD（图2(b）、（d）、（f）、（h）、（j）、（l））结果表明：包头站（R=0.07、RMSE=1.34）和台北站（R=0.22、RMSE=0.36）的反演效果非常差，其余4个站点反演结果相对较好。与RB AOD不同的是，北京站有较为严重的高估（42.5%），而南方站点徐州、台北、香港站低估较为严重，低估样本点占比分别为55.8%、44.2%、51%，香河站高估和低估的样本点占比接近。北京站和徐州站拟合直线斜率接近于1效果较好，而台北、香港站拟合直线斜率分别为0.232和0.386效果较差。

总体而言，RB AOD的的反演精度较H8 AOD产品精度高。虽然RB AOD在包头站整体精度很低，当AOD较大时，反演效果差，但在AOD较小时，RB AOD与地基站点AOD较为接近。

本研究选取北京站AERONET AOD验证RB AOD和H8 AOD数据在不同季节的时间变化一致性。图3分别给出了2018年4月（图3(a)，春季）、8月（图3(b)，夏季）、10月（图3(c)，秋季）、12月（图3(d)，冬季）地面观测AERONET AOD与RB AOD和H8 AOD时间序列变化曲线。春季的RB AOD与地基站点AERONET AOD的走势较为一致，在AOD>1时，小部分时段的RB AOD较AERONET AOD出现了明显高估；H8 AOD与AERONET AOD的走势也较为一致，但H8 AOD总体相对AERONET AOD偏高，且AOD值跳动较大。图3（b）为8月份的时间序列对比，与4月份相似的，当AOD>1时，RB AOD较AERONET AOD依然会出现明显高估，当AOD< 1时，RB AOD较AERONET AOD走势基本重合；8月份的H8 AOD总体相对AERONET AOD依然偏高，H8 AOD时间序列趋势较地基AERONET AOD较为一致，但走势不平稳，AOD值跳动较大。在秋冬季的10月（图3(c)）和12月（图3(d)），AERONET AOD显示北京站的AOD值总体为较低水平，大部分时段AOD 值在0.1水平以下，部分时间段AOD值有一个急剧升高的过程，可能对应有过程性雾霾天气。RB AOD与AERONET AOD的走势基本重合，但RB AOD在AOD > 0.7时会出现明显高估。而H8 AOD与AERONET AOD的走势偏离较大，当AERONET AOD值较低时，H8 AOD会出现很明显的高估，且有很明显的日变化特征，日际变化大，其总体AOD水平较AERONET AOD出现了明显的高估。

MODIS AOD产品因为其有较高的精度，常用于检验其他卫星的AOD产品。选取华北区域MODIS过境时刻的H8遥感影像，利用红蓝比值法反演AOD，将AOD反演结果与MOD04_3k AOD产品进行空间分布一致性检验。图4为2018年11月1日、2019年4月3日、1月2日MODIS过境时刻的Himawari-8真彩色合成、RB AOD、MOD04_3K AOD的空间分布对比图。

从图4(a)的H8 RGB真彩色合成图可见，2018年11月1日华北地区有区域性雾霾天气过程，特别是河北南部、河南北部一带灰霾较为严重，华北北部和华中地区能见度较高，图4(b)和图4(c)分别是对应时段的RB AOD反演结果和MODIS AOD产品的空间分布，从图中可见，两者AOD空间分布一致性很高，河北南部的AOD高值区（AOD >1.2）分布基本一致，华北北部和华中地区AOD值均较小，当AOD > 0.7时，RB AOD相对 MOD04_3K AOD值略偏小。MOD04_3K AOD在河北南部AOD高值区有部分缺值，可能是由于强雾霾天气下的大气散射作用导致2.1 μm处反射率较大，被误判为亮地表，因此无法使用暗像元法反演AOD。

从春季2019年4月3日反演结果来看，RB AOD与MOD04_3K AOD的空间、分布一致性也较好，AOD高值区主要分布在河南安徽等地，AOD低值区主要分布在山东半岛和华北北部。由于MODIS幅宽有限，因此在同一时刻下，处在另一轨道的河北西北部以及内蒙古一带没有AOD值，其他地域MOD04_3K AOD产品缺失值较多，AOD空间分布连续性较差。得益于Himawari-8卫星的大范围观测能力，RB AOD在这些地区反演能力较强，空间分布连续性较高，范围范围能覆盖更多区域。

从冬季2019年1月2日反演结果可以看出，由于冬季北方地区植被稀少、地表反射率高，暗像元法的MOD04_3K AOD产品在华北区域几乎没有反演值，而基于地表反射率比值库反演的RB AOD依然有能力反演得到大范围的AOD分布，对比图3(h)的真彩色影像，其反演得到的AOD空间分布能反映出华北地区的灰霾分布特征。

对4个站点按季节分别进行检验，季节划分标准为：春季（3~5月）、夏季（6~8月）、秋季（9~11月）、冬季（12~2月）。利用AERONET AOD按季节划分分别对RB AOD和H8 AOD进行对比检验，表2列出了各站点季节检验统计参数。

从相关系数R来看，除徐州站秋季的H8 AOD结果较好外，其余站点各季节RB AOD的相关系数R均较H8 AOD的相关系数高。包头站的H8 AOD在4个季节的反演结果都很差（-0.36<R<0.14），但RB AOD在包头站的夏季和秋季的相关系数（0.74<R<0.8）较H8 AOD的相关系数（-0.36<R<-0.06）有很大提升，春冬季RB AOD反演结果和H8 AOD一样差。位于华北的北京、香河站的RB AOD在秋冬季的反演结果（0.77<R<0.94）大大优于H8 AOD结果（0.32<R<0.82），北方地区冬季反演效果较差的原因可能是冬季地表反射率发生了较大改变，根据气溶胶反演公式（式1），卫星接收到的大气层顶表观反射率ρTOA与路径程辐射项ρ0和地表反射率ρs有关。对于短光波段，在ρs很小的情况下，路径辐射项对ρTOA的贡献比较大，当ρs很高时，地表反射率成为主要贡献，冬季北方地区植被较少，裸露地表的反射率较高，使用传统暗像元法反演的H8 AOD对亮地表区域的地表反射率估计能力不足，导致AOD反演误差较大。徐州站在秋季时H8 AOD反演结果较好，其他季节RB AOD反演结果较好。南方站点香港和台北站4个季节的RB AOD反演结果均优于H8 AOD，由于春夏季华南地区大气水汽含量较高，因此H8 AOD反演能力不足效果较差，秋冬季华南地区植被覆盖率依然较高且大气水汽含量较少，暗像元法反演的H8 AOD结果与RB AOD结果相差不大。

从RMSE来看，北京、香河站除夏季时的H8 AOD的RMSE较小外，其他季节均较RB AOD的RMSE高，尤其是春季和冬季大幅偏高。香港、台北站春季RB AOD和H8 AOD的RMSE均较高，表明南方地区站点春季可能容易受残云的影响，其反演结果与地基站点相比离散较大。徐州站在春夏季时 RB AOD与H8 AOD的RMSE相近且较高，在秋冬季时 RB AOD的RMSE（0.15~0.16）优于H8 AOD（0.23~0.28）。包头站的H8 AOD的RMSE全年都很大（0.87~1.49），RB AOD春季RMSE也很大（0.78），但在其他季节的RMSE较小（0.13~0.14），表明RB AOD在西北干旱区的反演结果较H8 AOD产品有很大的优势。

为减小单波段地表反射率的误差，在反演过程中，与直接采用蓝波段地表反射率与地表反射率库对比不同的是，本文采用MODIS红蓝波段地表反射率的比值R₁/R₃与Himawari-8红蓝波段地表反射率的比值ρ_red/ρ_blue做对比，以减小单波段反演可能造成的误差，该方法也被用于极轨卫星NPP的VIIRS传感器的气溶胶光学厚度反演^[21]，相比于目前NPP-VIIRS业务化运行的AOD产品结果有一定改进，反演结果尤其是在亮地表区域改善明显。图5为亮地表区域2018年2月~12月的逐8天的红波段（RED）地表反射率、蓝波段（BLUE）地表反射率、红/蓝（Red/Blue）地表反射率比值曲线，从图中可以看出，单波段的地表反射率年内月内变化较波段比值的变化大。研究表明，0.01的地表反射率误差会带来0.1的AOD反演误差^[3]，AOD反演精度高度依赖地表反射率的精确度。因此红/蓝地表反射率比值库可减小单波段地表反射率库的误差，有利于提高AOD的反演精度和稳定性。

图6为Himawari-8 和MODIS 红蓝波段光谱响应函数曲线对比，从图中可以看出MODIS红蓝波段的光谱响应函数曲线与Himawari-8 AHI 红蓝波段的光谱响应函数曲线较为接近，中心波长也较为接近，因此将MOD09A1地表反射率产品作为先验地表反射率库反演Himawari-8 AOD是可行的。

本文基于Himawari-8 AHI卫星数据和MOD09A1地表反射率产品，基于改进的地表反射率库的红蓝波段比值库反演方法，反演了中国区域2018年4月~2019年4月逐10 min的AOD。利用位于中国地区西北、华北、华东、华南等地的6个AERONET地基观测资料对AOD反演结果进行精度验证，同时将反演结果与Himawari-8官方发布的AOD产品以及MODIS AOD产品进行时间序列、空间分布一致性对比,得出以下结论:

（1）红蓝比值库反演方法反演的RB AOD在6个站点整体效果均优于H8 AOD产品。冬春季时RB AOD在北方反演效果较H8 AOD 有很大改进，弥补了暗像元法在高地表反射率区域反演能力不足的缺点，提升了冬春季AOD反演精度和反演区域范围。夏秋季节RB AOD和H8 AOD的反演效果都较好，RB AOD略优于H8 AOD结果。

（2）RB AOD与MOD04_3K AOD产品在空间分布上有很好的一致性，且RB AOD在空间上的连续性更好，弥补了MODIS 暗像元法在高反射率地表无法反演的缺点，RB AOD在冬季时的反演能力较MOD04_3K AOD产品有很大提升，扩大了反演范围，同时相较于MODIS 1 d的重访周期，RB AOD10分钟的高时间分辨率在监测区域性污染天气的变化有很大优势。

（3）在构建查找表时气溶胶类型仅选择了大陆型气溶胶类型，虽然大陆型气溶胶类型有较好的普适性，但不同地区不同季节的气溶胶类型是不同的，气溶胶类型的选择对AOD反演精度的影响较大，后期研究可考虑建立多种气溶胶类型查找表，在反演过程中根据多种气溶胶类型反演结果选择误差最小结果。