图1 中国陆地范围内346个地面气象观测站点的分布情况
Fig. 1 Distribution of 346 surface meteorological observation stations in China
摘要 为了获得大范围的能见度空间分布信息, 提出一种把 MODIS 卫星气溶胶光学厚度陆面产品转化为地面能见度分布的反演算法。该算法以 CALIPSO 卫星的气溶胶标高为更新观测值, 以 GEOS-Chem 模拟的气溶胶标高场为背景场, 通过最优插值的资料同化方法, 将两者整合成误差更小的气溶胶标高分析场, 然后通过气溶胶标高分析场, 把 MODIS 卫星气溶胶光学厚度转化为地面能见度。与中国地区地面气象多年观测资料的对比结果表明, 反演能见度与观测能见度之间, 点对点的月相关系数可达 0.5 以上, 两者多年的逐月变化趋势与地域分布形势较为一致。
关键词 地面能见度; 卫星反演; 气溶胶光学厚度; 气溶胶标高
能见度是具有正常视力的人在当时的天气条件下能够看清楚目标轮廓的最大距离[1], 是表征大气稳定度的天气指标, 也是保护交通运输安全的重要参考因素[2‒4]。目前, 在实际气象业务中, 地面能见度观测主要有人工观测和仪器观测两种方式[5]。这两种方式都有无法避免的局限性, 即无法得到能见度的空间分布信息。卫星遥感技术可以弥补这方面的不足。
卫星遥感反演地面能见度是基于气溶胶光学厚度(aerosol optical depth, AOD)、大气消光系数以及地面能见度三者的关系, 将卫星观测到的气溶胶浓度及垂直分布信息转化成地面能见度的过程。近年来, 随着人们对大气污染的关注, 有关 AOD 与能见度之间相关性的研究引起众多学者的兴趣。Retalis 等[6]对比塞浦路斯多个地面测站的 MODIS AOD、手持太阳光度计 AOD 以及地面能见度数据, 发现MODIS AOD 与由地面能见度估计的 AOD 有显著的相关性。陈静静等[7]定性地分析青岛 2006 年春季22 个晴天的 AOD 与实测水平能见度的关系, 并从浑浊度系数和消光系数两个角度计算水平能见度估算值。赵秀娟等[8]和李霞等[9]证明 AOD 与能见度之间近似地呈指数递减关系。张浩等[10]利用 MODIS气溶胶产品和气象观测资料, 分析能见度与 MODIS气溶胶产品和气象要素的相关性, 并运用多元回归方法, 建立安徽省不同季节能见度的反演模型。闫晓庆等[11]利用 2009 年 5—9 月的 MODIS 资料以及SACOL 站实测资料, 讨论兰州地区 AOD 与水平能见度的 4 种常见转换模型的适用性。Lee 等[12]提出一种基于 MODIS 卫星观测资料的分析反演方法, 用来估算气溶胶浓度、消光系数和地面能见度的空间分布, 并对韩国首尔地区 2007—2009 年 8 个集中观测时段内的上述 3 个物理量进行观测和统计, 发现相对于 Koschmieder 的关系式, 利用 MODIS 估算的能见度数值更接近观测值, 但这种反演方法的缺陷在于仅适用于无云的大气状况, 并且需要根据指数递减廓线预先计算气溶胶标高。Kessner 等[13]对美国东海岸 12 年的 ASOS 站点地面消光系数数据进行质量控制, 然后用 GOES-5 模式的气溶胶垂直廓线对标高进行订正, 建立两种由 AOD 获取能见度的模型, 并侧重讨论边界层高度内外的气溶胶分布情况对地面能见度反演的影响, 印证了气溶胶垂直分布对标高确定的重要性。
综上所述, 多数研究对气溶胶标高的考虑比较单一, 没有将实际观测与模拟结果结合起来提高其方法的准确性。同时, 这些研究往往只着眼于局地范围, 缺乏对大尺度区域(如中国陆地范围)能见度反演的研究。目前, 国内还没有发现通过卫星资料反演地面能见度的产品, 因此我们设计了本文的算法, 即通过资料同化手段, 将卫星观测与模式模拟得到的气溶胶标高进行整合, 形成误差更小的气溶胶标高区域分布, 然后用同化后的标高场对卫星观测的 AOD 进行转化, 反演得到中国陆地范围地面能见度的分布情况。
本研究使用的数据包括中国陆地范围地面气象观测站的水平能见度数据、MODIS 卫星 AOD 数据、CALIPSO 卫星气溶胶消光系数廓线数据以及微脉冲激光雷达的气溶胶消光系数廓线数据, 使用的模式为 GEOS-Chem 大气化学传输模式。其中, 气溶胶标高由 CALIPSO 观测廓线和 GEOS-Chem模式模拟获取, 二者同化后的分析场结合 MODIS AOD 进行能见度转化。用微脉冲激光雷达数据获得的气溶胶标高验证 CALIPSO 和 GEOS-Chem 的标高结果, 用地面气象观测站的水平能见度数据验证能见度反演结果。
中国陆地范围地面气象观测站的水平能见度数据来自美国国家海洋和大气管理局(National Oceanic and Atmospheric Administration, NOAA)资料库(http: //www.cdc.noaa.gov), 共计 346 个站点(图 1), 时间跨度为 2010 年 1 月至 2014 年 12 月。根据测站不同, 观测时次密度有逐小时观测结果和逐 3 小时观测结果。由于太阳同步轨道卫星过境时间相对固定, 可以与其匹配的观测时次集中在 UTC 3—6 时。
MODIS (MODerate-resolution Imaging Spectro-radiometer, 中分辨率成像光谱仪)光学厚度产品为Level2 MYD04 气溶胶产品中的 550nm 陆面气溶胶光学厚度数据, 空间分辨率为 10km×10km, 时间跨度为 2010 年 1 月至 2014 年 12 月。
CALIPSO (Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations, 云‒气溶胶激光雷达与红外探测者卫星)气溶胶廓线产品为 Level2 5kmAPro(CAL_LID_L2_05kmAPro-Standard-V4-10)中的 532nm 消光系数廓线数据, 星下点距离水平分辨率为 5km, 时间跨度为 2006 年 6 月 13 日至 2017年 12 月 31 日。对垂直方向上的消光系数进行积分, 可以得到光学厚度, 进而可以确定气溶胶标高。
MPL(Micropulse Lidar, 微脉冲激光雷达)是一种后向散射激光雷达系统, 理论上最低探测高度距地面 95 m, 垂直分辨率为 15 m, 可自动地进行全天连续的观测。对 MPL 接收的后向散射信号进行一系列的订正和处理, 可以得到气溶胶消光系数的垂直廓线。本研究所用的 MPL 安装在北京大学物理楼顶(39°59.5′N, 116°18.6′E), 数据为 30 分钟均值的廓线时间序列, 时间跨度为 2016 年 7 月至 2017 年12 月, 其标高的获取原理与 CALIPSO 相同。
图1 中国陆地范围内346个地面气象观测站点的分布情况
Fig. 1 Distribution of 346 surface meteorological observation stations in China
GEOS-Chem (http://acmg.seas.harvard.edu/geos/)是一个全球性的三维大气化学传输模式, 可以模拟大气中的各种物理、化学变化过程, 得到各种大气成分的浓度分布[14]。GEOS-Chem 由哈佛大学大气化学模式研究组开发并提供技术支持[15], 广泛地应用于大气化学和空气质量研究中。驱动模式的气象产品来自美国国家航空航天局(National Aeronautics and Space Administration, NASA)全球模式与资料同化办公室(Global Modeling and Assimilation Office, GMAO)提供的从 Goddard地球观测系统(Goddard Earth Observing System, GEOS)同化而得的再分析数据, 人为设定地面排放清单。该模式能够模拟大气中沙尘、黑炭、有机碳、粗海盐、细海盐和硫酸盐 6 种主要气溶胶光学厚度的垂直分布信息, 进而计算得到气溶胶标高。本研究中用到的模式结果均为水平分辨率 2.5°(经度)×2°(纬度), 垂直方向 47 层的全球模拟结果。提供背景协方差矩阵多年月平均时间序列的模式版本是 v10-01, 时间跨度为 1998 年 1 月至 2009 年 12 月, 气象数据产品为 MERRA。提供每日背景场资料的模式版本为 V11-02rc, 时间跨度为 2010 年 1 月至 2014 年 12 月, 根据模拟年份不同, 所用气象数据主要有 MERRA2 和 GEOS-FP两个产品。有关排放的设置均为模式默认值。
由于 CALIPSO 卫星只能探测其星下点的气溶胶垂直分布, 且脉冲信号进入地球大气层后的衰减很强, 因此最终得到的气溶胶廓线数量十分有限。鉴于其观测特点, 以月为时间单位, 对 2006—2017年的 CALIPSO 昼间观测数据进行统计, 对落在中国陆地范围内的经纬度分辨率为 1°×1°的网格内所有 532nm 消光系数廓线的标高逐月进行平均, 最终得到逐月的气溶胶标高分布情况(图 2)。
从图 2 可以看出, CALIPSO 的气溶胶标高分布呈现较明显的地域和季节差异。地域方面, 标高的变化与人口密度分布差异呈现相反的趋势, 青藏高原、新疆北部及东北北部等地区人口稀少, 空气比较洁净, 标高较高, 华北平原、中国东部和东南部沿海等人口和城市相对密集地区的地面气溶胶排放量大, 主要分布在边界层内, 因此标高相对较低。季节方面, 夏季标高明显高于冬季, 主要原因在于夏季地表加热明显, 边界层较高, 这也符合气溶胶垂直廓线的季节变化实际情况。
本文对 1998—2011 年 MERRA 气象数据产品驱动下的 GEOS-Chem 2.5°(经度)×2°(纬度)分辨率的全球模拟结果进行统计, 选取包含中国陆地区域的模式格点, 对上述 6 种气溶胶光学厚度的垂直分层模拟数据求和, 得到总光学厚度的垂直分布, 然后计算多年标高的月平均值, 得到逐月的气溶胶标高分布情况(图 3)。
从图 3 看出, 模式的结果也在一定程度上呈现地域和季节分布差异。地域方面, 总体上西高东低, 北高南低, 高值区主要分布在新疆、青海、西藏以及东北地区的西北部。季节方面, 从冬至夏, 高值区逐渐从新疆、内蒙古和东北地区的西北部转移至青藏高原的中东部, 然后从夏至冬再以相反的趋势转移回去。中东部以及东南部地区为标高的低值区, 季节差异比较小, 区域内的地域差异也不明显。与 CALIPSO 的结果相比, 模式标高的数值分布无论在时间还是空间上的变化都更平滑, 局地性差异也大大地小于前者。这种分布与选定的模式格点精度有关, 也与模式本身模拟的效果有关。由于本文选择的网格点较粗, 模式对区域的细节处理也比较粗, 但这样的精度可以满足一定范围的区域内能见度反演的需求。此外, GEOS-Chem 在时间和空间上的输出结果比较灵活, 能够设置输出与卫星过境时空匹配的结果, 可为后续研究提供方便。
使用 MPL 数据来验证 CALIPSO 气溶胶标高的准确性。MPL 和 CALIPSO 都是通过气溶胶消光系数廓线确定气溶胶标高。对 2016 年 7 月至 2017 年12 月的 MPL 与 CALIPSO 观测结果进行空间和时间上的匹配, 对与 MPL 所在地点直线距离不超过 50 km, 观测时间前后半小时内的 CALIPSO 廓线标高取平均值, 作为对应的卫星观测结果。因为 CALI-PSO 每日过境时间是固定的, 只有下午 13 时整和13 时 30 分的 MPL 观测时次才有机会与之匹配, 而两个时次比较接近, 因此取两个时次的平均值为MPL 观测结果, 最终得到的统计结果如图 4 所示。
图2 CALIPSO观测的中国陆地范围内气溶胶标高多年月平均分布情况
Fig. 2 Monthly distribution of aerosol scale height in China based on many years’ profile data observed by CALIPSO
图3 GEOS-Chem模拟的中国陆地范围内气溶胶标高的多年月平均分布情况
Fig. 3 Monthly distribution of aerosol scale height in China based on many years’ simulation of GEOS-Chem
2016 年 7 月至 2017 年 12 月, 与 CALIPSO 廓线相匹配的 MPL 廓线共有 15 条。从图 4 可以看出, 绝大多数情况下, 由 CALIPSO 观测得到的标高比GEOS-Chem 模式模拟的标高更接近 MPL 实地观测的标高。尤其是在春季和夏季, 用模式的格点标高来代表 MPL 所在位置的局地标高会产生很大的误差。这对能见度的反演非常不利, 因为能见度是一个局地性很强的指标, CALIPSO 的观测点虽然稀疏, 却可以保留气溶胶廓线的局地特征。因此, 模式模拟的标高对局地廓线模拟的准确性是有限的, 将CALIPSO 的观测结果同化进模式结果对提高反演能见度的准确性十分必要。
能见度反演算法中用到的资料同化方法是常用的最优插值法(optimal interpolation, OI)[16]。就本文所要解决的问题而言, 其基本思路是通过对 GEOS-Chem 模拟的背景值 xb 和 CALIPSO 观测值 xo 进行误差分析, 给出一个权重系数, 使得两种资料加权平均后的分析值 xa的误差尽可能地小。在实际问题中改写成矩阵形式, 通过式(1)~(3)计算得到分析场向量xa及其协方差矩阵 Pa。
(1)
(2)
, (3)
其中, xb为背景场向量, K 为权重系数矩阵, z 为观测值构成的向量, H为前向算子, H 为其换算矩阵, I为单位对角矩阵。背景协方差矩阵Pb由 GEOS-Chem模式每个格点 1998—2011 年背景值时间序列依次减去多年月平均值构建。观测协方差矩阵 R 中的代表性误差是通过统计 2006—2017 年落在每个模式格点范围内的所有 CALIPSO 廓线的标高数值, 然后取标准差得到。经统计, 将卫星测量误差假定为0.1 km。对每个 CALIPSO 观测标高, 减掉与其距离最近的 4 个模式格点值的面积加权平均值后, 得到观测值对背景值的更新增量。
此外, 考虑到 CALIPSO 廓线的局地性, 在同化观测结果的时候, 选定高斯函数
(4)为模式格点中心与观测点间水平距离
的相关系数函数[17], 以格点所在纬度的两倍纬向格点间距为高斯半径 L, 对背景协方差矩阵的元素乘以距离权重系数, 使卫星观测的结果对背景场的改善局地化。
将 MODIS 卫星观测到的 AOD(用
表示)视为气溶胶消光系数 kext 从地面到大气层顶的积分[18]:
。 (5)根据著名的Koschmieder公式:
, (6)
图4 MPL, CALIPSO和GEOS-Chem的气溶胶标高对比
Fig. 4 Comparison among MPL, CALIPSO and GEOS-Chem aerosol scale height
地面水平能见度
与地面消光系数 kext(0)成反比关系[19], 因此 AOD 与能见度通过消光系数建立联系。因为气溶胶主要由地面排放, 经湍流输送和扩散, 均匀地分布在边界层, 在边界层之上随高度迅速递减[20], 因此可以假设气溶胶的水平分布是均一的, 且光学性质不变, 大气密度随高度呈指数级递减, 气溶胶消光系数也随高度呈指数级递减:
。 (7)这里, H为大气标高。仿照其他大气物理量的标高定义, 我们定义从某一高度起, 向上至大气层顶的光学厚度为总光学厚度 1/e的高度为气溶胶标高。将 kext 带入 τ 的表达式中, 再借用 V 与 kext 的反比例关系, 可得
。 (8)
其中, K*为Koschmieder 常数, 取值通常为 3.912。本文只借用 Koschmieder 公式的反比例关系作为算法设计的理论支持, 并没有借用 3.912 这一数值, 这样做的原因会在后面进行误差分析时加以说明。
利用同化后的分析场对 MODIS AOD 进行转化, 根据式(8)得到反演的能见度数值。转化过程中还需考虑实际情况下整层大气的光学厚度不仅包括 MODIS AOD, 还包括瑞利散射光学厚度, 后者通过式(9)[21]计算得到:
(9)其中, τR(λ)代表对某一波长光的瑞利散射光学厚度, σ表示单个粒子散射截面, P表示大气压强, A是阿伏伽德罗常数, ma表示空气的平均相对分子质量, g表示重力加速度。
经过 2.1~2.5 节对反演过程中各个技术细节的设计, 本文以方便业务化应用为目的, 形成如图 5所示的反演算法流程。
1)卫星观测和模式数据准备。在进行算法转化前, 参照式(2), 先计算 OI 权重系数矩阵。反演某日地面能见度时, 应下载当日 MODIS AOD 数据和CALIPSO 观测廓线, 运行 GEOS-Chem 得到卫星过境时段的气溶胶空间分布。提取可作为更新增量的CALIPSO 卫星廓线并计算其标高, 同时对 GEOS-Chem 模拟结果进行计算, 得到标高的背景场。
图5 反演算法流程示意图
Fig. 5 Flow diagram of the proposed retrieval approach
2)将 CALIPSO 标高同化进 GEOS-Chem 模式结果。按照式(1), 将 CALIPSO 廓线标高同化进GEOS-Chem 模式背景场中。考虑到气溶胶分布的日变化和卫星过境的时段, 本文选择北京时间每天下午 13—16 点间的模式结果作为气溶胶标高的背景场。同化处理是将相对准确的 CALIPSO 气溶胶廓线信息, 通过背景协方差矩阵和距离权重系数扩展到整个背景场, 得到包含局地化信息的分析场。
3)MODIS AOD 转化为能见度。按照 2.5 节的转化方法, 计算 MODIS AOD 像素点处对应的瑞利散射光学厚度并与之相加, 再利用前两个步骤操作得到的标高分析场来转化整层大气光学厚度, 得到MODIS AOD 像素点对应的地面能见度分布。转化过程中, 设定 K*=1。最后, 对反演结果进行验证。
为了验证本文算法的效果, 对中国陆地范围内地面气象观测站的水平能见度数据与通过转化MODIS AOD 得到的能见度反演数据进行点对点的比较验证。为了排除水汽对能见度的影响, 设定以下 3 个条件对能见度资料进行筛选: 1)观测时未见雾、降水或降雪等天气现象; 2)露点温度差大于2°C; 3)云量记录不可为“无法辨认(obscure)”。
对符合筛选条件的能见度站点观测时次进行MODIS AOD 资料的时空匹配, 要求卫星观测时间在能见度观测时间的前后半小时内, 卫星观测的像素点与能见度观测站点的距离小于 20km。取所有符合匹配条件的 MODIS AOD 像素点的平均值作为该能见度站点观测时次的AOD转化基准。
对符合上述筛选条件的能见度站点观测值, 通过匹配的 AOD 数值, 计算同一时间、同一测站的反演能见度值, 结果如图 6 所示。图 6(a)中的反演能见度值由同化前的背景场转化而来, 图 6(b)中的反演能见度值由同化后的分析场转化而来, 图 6(a)和(b)中散点数同为 4422 个, 黑色斜线由所有散点的线性拟合得到。
从图 6 可以看出, 地面气象观测的水平能见度数据存在严重的离散化现象, 制约了点对点验证的相关性, 只能通过对大量数据的统计来消除这种离散化带来的误差。经过对大量数据的拟合, 本文提出的算法能够克服观测数据离散化的缺点, 使线性拟合的斜率接近 1, 说明通过 MODIS AOD 转化气溶胶标高所得的地面能见度与观测能见度在一定程度上具有可比性。分析场转化而得的能见度的拟合斜率系数比背景场更接近 1, 拟合度指数 R2 有所提高, 均方根误差RMSE有所降低, 说明分析场转化能见度提升了数据整体的相关性, 减小了偏差, 使反演值更接近观测值。由于观测数据本身的缺陷, 数据整体的相关性仍然较弱, 因此基于目前仅有的地面气象观测水平能见度的数据, 还无法得到一个准确而合理的能见度反演模型。
图6 对2010年数据反演‒观测能见度散点图的线性拟合效果
Fig. 6 Linear fitting between retrieval and observation visibility based on data of 2010
对 2010—2014 年的反演‒观测能见度的点对点月相关系数结果进行统计, 得到多年点对点月相关系数折线图(图 7)。从点对点的相关性来看, 观测散点个数相同时, 分析场转化所得能见度的相关性优于背景场, 但由于 CALIPSO 观测的稀疏性以及同化过程中考虑了距离权重函数和 GEOS-Chem 模式与 CALIPSO 廓线的误差比, 所以同化处理对背景场的改善效果有限, 相关系数的提高并不显著。此外, 算法对不同月份和季节的反演效果不同, 1月、6 月和 10 月份的反演效果较好, 相关系数可超过 0.5。从不同季节的平均相关性来看, 除春季相关性较差外, 冬季、夏季和秋季的相关系数依次为0.4743, 0.4691和0.4242, 均超过0.4。
图7 2010—2014年反演‒观测能见度的多年点对点月相关系数折线图
Fig. 7 Monthly point-to-point correlation coefficients between retrieval and observation visibility based on data from 2010 to 2014
以季节为时间单位, 将 2010 年同化标高场转化所得的反演能见度和观测能见度分别格点化成 1°× 1°的中国陆地区域分布图, 对同一季节落在每个格点内的所有能见度取平均值, 得到区域分布形势(图 8)。受测站分布的地理位置所限, 很多格点没有能见度观测值, 卫星数据反演可以在一定程度上弥补这种地域限制。
从图 8 可以看出, 就整个中国陆地范围来讲, 反演能见度与观测能见度具有相同的区域分布形势。格点对格点的相关系数在冬季、春季、夏季和秋季依次为 0.8468, 0.8217, 0.7396 和 0.7885, 相关性均较强。对比上下两行的颜色深浅程度可以看出, 反演算法对能见度的估计在观测能见度的高值区更高, 低值区更低。这种差异主要是因为观测能见度的记录习惯导致其取值范围有限, 而反演能见度是 MODIS AOD 和分析场标高两种数据的运算结果, 数值变化范围远远大于前者。这也是本文算法的一个优势, 即能够反演得到变化范围更大、取值更精确的能见度数值。此外, 反演结果明显地保持着 CALIPSO 和 GEOS-Chem 两种标高良好的地域差异性和季节差异性。相对于观测能见度不够显著的季节性差异, 算法反演能够更显著地突出不同季节以及不同区域的能见度变化差异。
图8 2010年中国陆地范围反演能见度与观测能见度区域分布形势
Fig. 8 Regional distribution of retrieval and observation visibility based on data of 2010
图9 对中国陆地范围内6个不同验证区域的划分
Fig. 9 Regional division of six different regions for validation
为了检验算法对不同区域的适用性, 在中国陆地范围内划分 6 个区域(图 9 中 A 区、B 区、C 区、D 区、E 区和 F 区)对算法进行验证。以月为时间单位, 分别对 6 个验证区域内 2010—2014 年的所有反演和观测能见度取平均值, 得到两种能见度的多年月平均变化折线, 并计算得到两种能见度逐月时间序列之间的相关系数(图 10)。多年资料的统计平均可以将很多噪声平滑掉, 最终得到比较清晰的逐月均值及变率。就逐月时间序列变化趋势而言, 由本文算法得到的反演能见度与观测能见度的相关性很强, 除 F 区域的夏季偏差较大以外, 其他区域反演能见度与观测能见度的的逐月时间序列之间相关系数均大于 0.7。
从数值上看, 本文算法在 A 区和 B 区容易对能见度有比较大的“高估”, 这种现象与地面观测中对水平能见度的记录习惯有关。在本文研究的多年观测数据中, 地面水平能见度的最大值都没有超过 15 km, 而在能见度良好的情况下, 这些区域的水平能见度实际上远远大于 15 km这个观测记录上限。因此, 对这种“高估”不能完全理解成算法的偏差, 还需要进一步的验证。对于常年气溶胶浓度较高的 C区, 与观测能见度相比, 反演能见度逐月变率更大, 季节性变化更显著。在 D 区, 反演能见度出现 6 个区域中少有的低估现象。结合两种气溶胶标高的多年逐月分布情况(图 2 和 3)来看, D 区的标高常年偏低, 反演能见度的多年月平均值也包含这样的信息。对比 F 区与其他 5 个区域反演能见度与观测能见度的相关性可以看出, 算法对沿海地区夏季的适用性较差。
表 1 给出对 2010—2014 年反演‒观测能见度进行点对点统计得到的各月份误差系数。可以看出, 本文算法反演能见度与观测能见度相关性最好的是1 月, 相关系数达到 0.5348, 最差的 3 月相关系数只有 0.1854; RMSE 变化呈现夏高冬低的趋势, 变化范围在 3.61~8.55 之间, 平均值为 6.3017。与 Lee 等[12]在首尔地区的定点观测结果(R=0.882, RMSE= 2.252)以及 Kessner 等[13]对美国东海岸 ASOS 站点的统计结果(R=0.70, RMSE=7.82)相比, 在相关性方面有一定的差距。主要原因在于, 本文算法验证中用到的地面气象观测资料大多是基于能见度目标物图的人工观测结果, 其数值本身的离散化就造成一定的误差, 再加上能见度观测记录上限的设定, 其不确定性就更不言而喻了。Lee 等[12]和 Kessner 等[13]研究的地面能见度或来自连续而准确的定点观测, 或进行了严格而细致的质量控制, 能够积累大量准确的实测数据, 最终可以作为统计模型的验证基准, 对统计模型的误差有显著的改善作用。本文算法从卫星观测角度出发, 着眼于整个中国陆地范围的区域能见度反演, 结果的精度无法做到与定点观测相提并论, 只能在逐月的时间尺度上得到比较显著的相关性及变率。在后续的工作中, 能够保证数据质量的地面能见度观测资料的积累是对算法进行验证和改进的重要前提。
图10 2010—2014年6个验证区域的反演‒观测能见度多年月平均变化折线及逐月时间序列相关系数
Fig. 10 Monthly mean plots and time series correlation coefficients of retrieval and observation visibility in the six given regions based on data from 2010 to 2014
表1 2010—2014年反演‒观测能见度多年点对点月误差指数
Table 1 Monthly point-to-point error indexes between retrieval and observation visibility based on data from 2010 to 2014
月份RRMSE月份RRMSE 10.53483.6146 70.40628.5481 20.41374.4051 80.49927.9988 30.18545.1755 90.36938.5568 40.20686.2488100.50995.6143 50.26807.7613110.39335.3470 60.50187.7851120.33534.5644
此外, 本文算法的误差还来自以下几个方面。
1)卫星观测数据的误差。除 MODIS 反演 AOD产生的误差外, CALIPSO 廓线每层的气溶胶消光系数也有一定的不确定性, 导致由其积分而得的气溶胶标高的观测误差难以量化, 在进行同化处理时只能人为假定。
2)GEOS-Chem 模式背景场的误差。本文选取的全球模拟格点分辨率对能见度这类局地物理量的刻画显然不够细致, 在后续的研究中可以通过区域嵌套模拟来提高局地能见度的反演精度。此外, 随着模式使用的排放清单和相关气象数据产品的日益完善, 其模拟效果也会越来越好, 这是本文算法在后续研究中可以继续依托 GEOS-Chem 模式建立标高背景场的重要保证。
3)对气溶胶垂直廓线的假设产生的误差。本文算法是基于气溶胶浓度指数递减垂直廓线的假定建立的, 而实际大气中气溶胶的垂直分布情况远比这种假定复杂, 这给标高的确定带来误差。此外, 本文算法未考虑边界层、高层输送和标高日变化等因素对气溶胶垂直分布的影响, 需要在后续工作中对这些特殊的气溶胶垂直分布做进一步的讨论。
4)Koschmieder 常数设定产生的误差。首先, Koschmieder 常数 3.912 的适用条件苛刻, 受众多因素制约, 无法直接作为算法的相关参数使用。再有, 本文设计的算法不完全等同于由地面消光系数换算的水平能见度, 是从卫星观测的角度, 通过气溶胶标高转化 AOD 得到地面能见度, 还需要经过大量的统计验证才能给出一个定量的反演系数。根据气溶胶的组成和种类、地表状况以及季节等因素的变化, 该反演系数可能有很大的变化和误差, 本文在对反演结果的讨论中暂未考虑这些差异。
本文对目前可以获取气溶胶垂直分布廓线的 3种手段进行资料统计, 将 CALIPSO 卫星观测到的532nm 气溶胶消光系数廓线和 GEOS-Chem 模式模拟的 6 种常见大气气溶胶光学厚度垂直分布换算成气溶胶标高, 通过多年资料的积累, 形成两者在中国陆地范围内的气溶胶标高逐月分布形势, 表明气溶胶标高可以在一定程度上反映气溶胶浓度的垂直分布。在利用 MPL 雷达观测数据对上述两种标高的验证中, 发现 CALIPSO 廓线比模式的模拟结果更接近局地的气溶胶垂直分布特征。
基于 Koschmieder 公式以及气溶胶光学厚度与消光系数之间的关系, 在假设大气气溶胶浓度在垂直方向呈指数递减分布并排除云和降水影响的前提下, 本文提出利用 MODIS AOD 反演地面能见度的算法。在将气溶胶标高转化成 AOD 的过程中, 加入最优插值的数据同化方法, 将 CALIPSO 卫星廓线资料与 GEOS-Chem 模式模拟结果整合成误差更小的气溶胶标高分析场, 再用分析场将 MODIS AOD 转化成地面能见度分布。与多年地面气象观测资料的对比验证结果表明, 本文提出的算法对地面能见度的反演效果较为理想, 月时间尺度上的点对点相关系数可在 0.5 以上, 逐月变化趋势以及区域分布形势也与地面观测十分接近, 有望在针对能见度的实际观测和预报业务中推广使用。
关于算法的准确性和适用性还有很大的讨论空间, 算法设计还存在很多值得探究的问题。云、降水等影响相对湿度的天气现象也是影响地面能见度的关键因素, 本文的算法未能考虑。由于缺乏准确而充足的实测能见度资料, 本文未能对算法的误差进行定量的分析。在对算法进行理论推导的过程中, 对特殊气溶胶垂直廓线的标高确定未做考虑, 对 Koschmieder 常数的适用性讨论不足。后续工作中, 在完备算法理论设计的基础上, 基于准确的地面能见度观测资料的算法验证非常重要。反演能见度模型的建立, 需要对不同区域进行更多的统计验证才能实现。
参考文献
[1] 盛裴轩, 毛节泰, 李建国, 等. 大气物理学. 北京: 北京大学出版社, 2003
[2] Ying Q, Mysliwiec M, Kleeman M J. Source ap-portionment of visibility impairment using a three-dimensional source-oriented air quality model. En-vironmental Science & Technology, 2004, 38(4): 1089‒1101
[3] Huang W, Tan J, Kan H, et al. Visibility, air quality and daily mortality in Shanghai, China. Science of the Total Environment, 2009, 407(10): 3295‒3300
[4] Thach T Q, Wong C M, Chan K, et al. Daily visibility and mortality: assessment of health benefits from im-proved visibility in Hong Kong. Environmental Re-search, 2010, 110(6): 617‒623
[5] 张文煜, 袁久毅. 大气探测原理与方法. 北京: 气象出版社, 2007: 24‒30
[6] Retalis A, Hadjimitsis D G, Michaelides S, et al. Comparison of aerosol optical thickness with in-situ visibility data over cyprus. Natural Hazards and Earth System Sciences, 2010, 10: 421–428
[7] 陈静静, 盛立芳, 郑元鑫. 2006 年春季青岛气溶胶光学厚度与气溶胶指数、能见度的关系分析. 中国海洋大学学报(自然科学版), 2010, 40(1): 10‒16
[8] 赵秀娟, 陈长和, 袁铁, 等. 兰州冬季大气气溶胶光学厚度及其与能见度的关系. 高原气象, 2005, 24(4): 617‒622
[9] 李霞, 刘艳, 王胜利. 新疆部分地区大气气溶胶光学厚度及其与能见度的关系//中国气象学会 2008年会大气环境监测、预报与污染物控制分会场论文集. 北京, 2008: 536‒542
[10] 张浩, 邓学良, 石春娥. MODIS 气溶胶产品在大气能见度监测中的应用. 环境科学与技术, 2015, 38 (8): 49‒55
[11] 闫晓庆, 张兴山. 气溶胶光学厚度与水平能见度转换模型分析[C/OL]//第 31 届中国气象学会年会 S6大气成分与天气、气候变化论文集. (2014)[2019‒ 04‒01]. http://cpfd.cnki.com.cn/Area/CPFDCONFArticle List-ZGQX201411006-6.htm
[12] Lee K H, Wong M S, Kim K, et al. Analytical ap-proach to estimating aerosol extinction and visibility from satellite observations. Atmospheric Environ-ment, 2014, 91: 127‒136
[13] Kessner A L, Wang J, Levy R C, et al. Remote sensing of surface visibility from space: a look at the United States East Coast. Atmospheric Environment, 2013, 81: 136‒147
[14] 燕莹莹, 林金泰, 旷烨, 等. 大气化学传输模型GEOS-Chem 全球‒区域双向耦合//创新驱动发展提高气象灾害防御能力——第 30 届中国气象学会年会论文集. 南京, 2013: 1‒28
[15] Bey I, Jacob D J, Yantosca R M, et al. Global modeling of tropospheric chemistry with assimilated meteorology:Model description and evaluation. J Geophys Res, 2001, 106: 23073‒23095
[16] Daley R. Atmospheric Data Analysis. Cambridge: Cambridge University Press, 1993: 763‒764
[17] Saroja P. Atmospheric data assimilation [EB/OL]. Toronto: University of Toronto. (2004‒09‒07) [2019‒04‒01]. http://www.atmosp.physics.utoronto.ca/ PHY2509/
[18] Coakley J A, Yang P. Atmospheric radiation: a primer with illustrative solutions. Berlin: Wiley-VCH, 2014: 81‒82
[19] 秦瑜, 赵春生. 大气化学基础. 北京: 气象出版社, 2003
[20] Koschmieder H. Theorie der horizontalen sichtweite II: kontrast und sicht-weite beitrage zur physik der freien. Atmosphere, 1925, 12: 171‒181
[21] Bodhaine B A, Wood N B, Dutton E G, et al. On rayleigh optical depth calculations. Journal of Atmospheric & Oceanic Technology, 1999, 16(11): 1854‒1861
Retrieval of Surface Visibility Using Satellite-Based Aerosol Measurements
Abstract In order to obtain the large scale spatial distribution of surface visibility, a novel technique is developed to convert satellite-retrieved aerosol optical depth (AOD) to surface visibility through aerosol scale height. Specifically, the aerosol scale height is calculated using the advanced optimal interpolation (OI) method by assimilating CALIPSO (Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations) measurements into GEOS-Chem stimulation results. The assimilated analysis field is then combined with MODIS AOD to calculate surface visibility. Validation against ground observation in China shows that best monthly correlation between collocated satellite-retrieved visibility and surface visibility data is above 0.5. Reasonable agreement is also found in seasonal and spatial variability.
Key words surface visibility; satellite retrieval; aerosol optical depth; aerosol scale height
doi: 10.13209/j.0479-8023.2019.124
国家自然科学基金(41575018)资助
收稿日期: 2019‒04‒02;
修回日期: 2019‒05‒15