刘丹丹 黄印博 孙宇松2) 卢兴吉 曹振松†

1) (中国科学院安徽光学精密机械研究所，中国科学院大气光学重点实验室，合肥 230031)

2) (中国科学技术大学研究生院，科学岛分院，合肥 230026)

3) (皖西学院，电气与光电工程学院，六安 237012)

对流层顶作为对流层和平流层之间的过渡层,对大气痕量气体浓度反演有非常重要的影响.理论分析对流层顶高对大气分子含量垂直分布的影响,并结合2018年8月(6—16日)拉萨观测数据,定量分析对流层顶高的变化对气体柱-平均摩尔分数Xgas反演的影响.结果表明: 对流层顶高的变化改变气体分子的垂直分布,XCO2,XCH4与对流层顶高度变化呈正相关,相关系数分别为0.998和0.780;XCO与对流层顶高度变化呈负相关,相关系数为0.994;而XH2O与对流层顶高度变化的相关性非常小.对流层顶高变化3 km,XCO2,XCH4 及 XCO 误差范围在 8.640%,0.035% 及 0.049% 以内.观测期间,XH2O,XCO2,XCH4 及 XCO 日平均值分别在 3432—4287,406.1~408.2,1.673—1.720 及 0.082—0.095 ppmv 之间变化.观测数据显示,该地区CO2,CH4气体柱-平均摩尔分数日变化趋势相似,XCO2与XCH4的相关系数大部分高于0.5.观测结果可为我国研究高原温带地区温室气体浓度的时空分布及其变化规律提供参考和第一手的直接观测数据.

1 引 言

对流层顶是对流层与平流层之间的过渡层,热力学上定义为垂直温度廓线上温度最低点所对应的层高度[1].顶高的变化与众多因素有关,例如纬度、季节以及天气活动等,通常对流层顶高度随纬度增高而降低.一般来说,热带地区对流层顶高度约 15—20 km,中纬度地区 8—14 km;夏季对流层顶的高度高于冬季,例如青藏高原夏季对流层顶平均在 17 km[2−5].1949 年,Goody 等[6−9]在研究中发现,对流层顶的形成与 CO2,H2O,O3浓度有重要的关系,臭氧总量也与对流层顶高度有关.反之,对流层顶高度的变化也改变气体的垂直分布,从而影响其在大气中的总含量.

利用地基傅里叶变换光谱仪遥测温室气体含量时,需要利用反演算法获得温室气体浓度结果.国际上地基遥测的反演算法多数是基于最优估算法或最小二乘法,例如总碳柱观测网(total carbon column observing network,TCCON)利用标准的反演算法GFIT[10]对CO2柱浓度进行反演;大气成分观测网 (network for the detection of atmospheric composition change,NDACC)利用中红外波段的反演算法SFIT获得O3,HF,CO及CH4等气体浓度[11];德国卡尔斯鲁厄理工学院开发的PROFFIT反演算法[11]反演了柏林等地区CO2,CH4等气体的柱浓度,其反演精度与GFIT反演算法获得的结果相当.三种反演算法的辐射传输模型、后向校正存在一定的差异,例如SFIT反演算法采用FSCATM辐射传输模型计算模拟光谱,而PROFFIT采用KOPRA辐射传输模型.上述反演算法从流程上可以分为两部分: 一是基于大气前向模型计算大气透过率光谱,二是通过优化反演参数使得计算值与测量值达到最小或最优,最终反演获得待测大气分子的浓度.在反演过程中,仪器线型函数、先验廓线、地表压力、太阳天顶角等因素均影响气体浓度反演的精度[12−17].例如Wunch等[17]研究了先验廓线对CO2等气体浓度反演的影响,并利用机载廓线准确获得 CO2,CO,CH4等气体分子的校准因子;Frey等[18]和Hase 等[19]研究了仪器线型函数(instrument line shape,ILS)变化对 CO2浓度反演的影响,当 ILS振幅调制效率 (amplitude modulation efficiency,ME)改变 4%,XCO2变化 0.035%.先验廓线作为前向模型中的重要参数,决定着反演结果的精度,而对流层顶高度的变化直接改变分子浓度的垂直分布特征,并最终影响反演结果,但目前有关对流层顶高度变化对温室气体浓度反演的影响尚未见报道.

拉萨市位于青藏高原中部,平均海拔3650 m,是海拔最高的城市之一,属于高原温带半干旱季风气候,全年多晴朗天气,太阳辐射强,空气稀薄,气温偏低.本文利用傅里叶变换光谱仪EM27/SUN观测了拉萨地区温室气体的柱-平均摩尔分数,针对该地区对流层顶高度变化对温室气体浓度反演影响展开了理论及定量分析,并获得了拉萨地区H2O,CO2,CH4及 CO 四种温室气体的柱浓度信息,测量结果可为高原温带半干旱气候温室气体源与汇的研究提供数据支撑和重要参考.

2 基本原理

温室气体浓度反演包括垂直柱浓度(vertical column density,VCD)和柱-平均摩尔分数（Dryair Mole Fraction，DMF）,反演方法是基于德国卡尔斯鲁厄理工学院 (Karlsruhe Institute of Technology,KIT)气象和气候研究所开发的PROFFIT (PROFILE FIT) 反演算法,将最优估算法和非线性逐次迭代结合,利用Tikhonov-Phillips约束条件,在对数尺度上对温室气体柱浓度进行反演,反演算法主要包含前向模型和后向反演两个部分.前向模型可以表述为

其中x是n维参数,包括未知的独立参量,例如垂直柱浓度,u代表固定参数,例如压力、温度等;y是m维测量值,F是非线性模型,可线性化为下列形式:

其中K是m×n维的雅可比矩阵,计算公式K=DF/Dx,x0是线性化参考点,且y0=F(x0).利用最小二乘法,使得测量信号(ymeas)与模拟信号y的差值Dy最小:

S是反演状态参量的协方差矩阵,采用高斯-牛顿算法解决非线性并执行多次迭代,第i+1迭代:

在大气遥感时,由于未知量多于实测量,此时需要约束条件求解方程.在反演过程中,常采取确定的廓线或者对先验廓线进行缩放的方法.常采用Tikhonov-Phillips约束条件,规则化方程为

其中,xa是变化参量的先验值集合,B是规则化矩阵,g是规则化参数,在约束条件下的第i+1次迭代为

由于测量仪器是低分辨率光谱仪,因此,采用缩放先验廓线的方式进行反演,并采用氧气作为内部标准反演气体柱-平均摩尔分数DMFs以减小仪误差:

其中Columngas是气体的柱总量,ColumnO2是氧气的柱总量,压强和温度廓线均来自美国国家环境预报中心 (National Centers Environmental Prediction,NCEP)的分析数据.

3 观测与结果分析

EM27/SUN光谱仪包含干涉仪和太阳跟踪器(图1(a)),探测器光谱响应范围 0.83—2.5 µm,最大光程 1.8 cm,光谱分辨率 0.5 cm–1,可实现 CO2,CO,H2O,O2及 CH4气体的同时观测.EM27/SUN光谱仪观测点位于拉萨市气象局内(经度91.135°E,纬度 29.659°N),站点海拔高度 3.625 km(见图1(b)).光谱记录时间与UTC标准时间相差8 h.选择晴朗无云的天气进行观测采集,光谱分辨率设置为0.5 cm–1.

图1 (a) 观测站点位置 (拉萨市气象局);(b) 观测设备 (傅里叶变换光谱仪,EM27/SUN)Fig.1.(a) Observing site (Lhasa meteorological bureau);(b) FTIR spectrometer (EM27/SUN).

3.1 数据质量控制

在反演之前,为了保证数据质量,需要对上述特殊情形进行预处理,包括干涉图的直流分量校准、傅里叶变换、相位校准以及数据重新采样等.引入两个数据质量控制标准,第一个标准是最大信号强度的绝对值需大于5%,即采集光谱的最大信号强度的绝对值幅度低于5%,剔除该数据;第二个标准是干涉图的变化不能超过10%.与此同时,为了保证获得高精度的温室气体浓度,还需要记录观测站点的海拔高度、经纬度、温、湿、压等参数,以及观测站点的实时气象参数,例如: 地表温度、地表压强等.此外,观测设备的定标误差对反演结果也有影响,故在观测前后需要测量仪器的线型函数消除此影响,本文利用LINEFIT软件[19]计算H2O分子谱线获得线型函数的调制效率和相位偏差,从而降低仪器线型函数对观测数据的影响.

3.2 对流层顶高的变化对反演的影响

3.2.1 理论分析

对流层顶高的变化改变大气分子含量的垂直分布,即大气分子垂直廓线.为了分析对流层顶对痕量气体浓度反演的影响,将气体分子浓度南北梯度的差异及对流层顶高等参数引入到先验廓线中,首先引入下列三个参量:

其中下标ref表示参考点,obs表示观测点,grad_lat(jgas)是第j种气体分子的南北半球梯度浓度,lat_ref参考廓线的纬度,lat_obs观测点的纬度,z(i)表示大气分层的第i层高度,ztrop是对流层顶高.第j种气体(jgas)在第i层上(ilev)的浓度与参考廓线浓度之间的关系为

公式(10)显示,对流层顶高对气体廓线的调制程度取决于气体的南北半球梯度、层高及当前层上气体的参考廓线浓度.因此,对流层顶高对每一种气体廓线的调制程度不同.由于50%的水汽分子集中在大约1.5 km(气压约850 hPa)高度以下,且 90%以 上 的 水 汽 限 制 在 5.5 km(气 压 约500 hPa)高度以下的大气层中,因而对流层顶高的改变对H2O廓线几乎没有影响.而对集中分布在对流层附近的分子气体廓线调制较大,例如CO2、CH4分子等.

3.2.2 定量分析

将对流层顶高由 10 km 增加到 18 km,获得不同对流层顶高度的先验廓线,并作为后续反演的先验值.CO2,CO,CH4,O2以及 H2O 的浓度反演波段如表1所列.目标气体反演波段存在一种或多种干扰分子,采用缩放干扰分子先验廓线的方式扣除其影响.

由于目前缺少我国的大气模式,在辐射传输和遥感探测领域,常用6种标准大气模式(1976年美国标准大气,中纬度夏季、中纬度冬季、副极地夏季、副极地冬季、热带)进行相关计算与反演,本文主要基于TCCON网先验廓线,利用3.2.1节的廓线理论,获得反演过程中涉及的8种气体分子廓线,如图2所示.

表1 反演波段Table 1.Inversion band.

由于分子在大气层中的分布不同,对流层顶高的变化对不同分子的垂直分布影响不同.例如H2O,HDO分子 (图2(a),图2(b)是 HDO的先验廓线乘以丰度比后的结果)主要分布在10 km以下,对流层顶高的变化对其廓线改变非常小,但CO2分子廓线受对流层顶高变化较大.从图2可以看出,对流层顶从 10 km 增加到 18 km,CO2分子10 km 以 上 的 廓 线 上 移 (图2(c)),CH4(图2(d))及N2O(图2(f))分子含量分布(廓线)随对流层顶高变化趋势类似,均在10—30 km高度范围内有变化,而对流层顶高变化对CO廓线的影响主要集中在 7—20 km 高度上.通常情况下,O2分子在大气中均匀混合,其含量是一个定值(0.2095 ppmv).分析对流层顶高的变化对温室气体浓度的影响,选取2018年8月16日的观测数据进行分析,当天天气晴朗、少云,且采集数据时间长,如图3所示.

图2 不同对流层顶高度下的 H2O (a),HDO (b),CO2 (c),CH4 (d),CO (e),N2O (f),HF (g) 及 O2 (h) 廓线Fig.2.The profiles of H2O (a),HDO (b),CO2 (c),CH4 (d),CO (e),N2O (f),HF (g) and O2 (h) at different heights of top troposphere.

观测结果显示,上午 11:30到下午 17:30,H2O气体分子含量逐渐降低,日变化与对流层顶变化无关,如图3(a)和图3(e)所示.XCO2,XCH4在下午14:00达到峰值,其含量随着对流层顶的增高而增大,且两者的变化趋势相似(见图3(b)和图3(c)),绝对偏差在中午14:00最小(见图3(f)和图3(g));XCO随着时间逐渐减小,XCO与对流层顶变化的关系与CO2,CH4正好相反,其含量随着对流层顶的增高而减小,进一步分析对流层顶变化对XCO2,XCH4及XCO的影响,将平均摩尔分数与对流层顶高进行线性拟合,如图4所示.

XCO2和XCH4与对流层顶高变化呈正相关,XCO与对流层顶高变化呈负相关,CO2及CO日平均摩尔分数与对流层顶的相关性较高,相关系数分别为 0.998、0.994,CH4的相关性稍低,相关系数为0.78.依据线性拟合公式,对流层顶高变化1 km,XCO2变化 2.880%,XCH4变化 0.012%,XCO变化0.016%,由于一天的对流层顶高变化不超过3 km,由对流层顶高引起的XCO2,XCH4及XCO 误差 分别为 8.640%,0.035% 及 0.049%,对流层顶高变化对不同气体分子的影响程度不同,由于H2O主要分布在5 km以下,所以对流层顶高的变化对 H2O几乎没有影响,但对 CO2,CH4及CO影响较大,尤其是CO2与CO气体分子,其原因是对流层顶高的增加,致使CO2和CO的廓线上移大于CH4廓线(如图2).

3.3 浓度反演结果分析

综合以上分析,在地基反演温室气体浓度时,需要考虑对流层顶高度的变化,以下结果均考虑当天对流层顶高,以NCEP分析数据为基础,获得实验期间的温度廓线,通过温度廓线获得当日对流层顶高的信息.观测时间2018年8月,选取天气条件 较 好 的 观 测 日 (8 月 6,7,8,10,12,13,14,16日),10次扫描平均获得 1组光谱,切趾采用Norton-Beer函数,共 2300 组光谱数据.由于观测日的天气不稳定,EM27/SUN需要在晴天薄云条件下进行观测,每天的观测时间段也不同,因此对日观测结果取平均得到观测期间的日平均值时间序列,如图5所示.

在观测期间,XH2O日平均值范围在3432—4287 ppmv内变化,在 2018 年 8 月 10 日 (180810)变化较大,其标准偏差为162.8;XCO2日平均值为 406.1—408.2 ppmv,XCH4日平均值为 1.673—1.720 ppmv,但XCO2、,XCH4日变化幅度较大,尤其是7日、8日两天,XCO2标准偏差分别为1.62和 1.35;XCH4标准偏差分别为 0.0071和0.0084;XCO 日平均值为 0.082—0.095 ppmv,日变化较为平缓,标准偏差小于0.0027.在观测期间,XH2O,XCO2,XCH4及XCO 的 平 均 值 分 别 为3919.70,406.87,1.689 和 0.091 ppmv.与同台设备在敦煌地区观测的数据差异较大,尤其是CO2和CH4气体浓度,原因可能是拉萨地区高海拔,空气含量低,导致其 CO2,CH4浓度均偏低,再者,在敦煌地区将EM27/SUN观测的XH2O数据与探空气球探测的进行了对比分析,发现两者之间相差13%,可能是由测量时间差、反演误差、探空气球观测误差等原因引起.

图3 XH2O (a),(e),XCO2 (b),(f),XCH4 (c),(g)及 XCO (d),(h)日变化随对流层顶高度的变化Fig.3.The diurnal variation of XH2O (a),(e),XCO2 (b),(f),XCH4 (c),(g) and XCO (d),(h) with tropopause height.

图4 (a) CO2,(b) CH4 及 (c) CO 平均摩尔分数与对流层顶的线性拟合Fig.4.Linear fit of the average mole fraction of CO2 (a),CH4 (b) and CO (c) to the top of the troposphere.

图5 XH2O (a),XCO2 (b),XCH4 (c),XCO (d) 的日平均序列Fig.5.Time series of XH2O (a),XCO2 (b),XCH4 (c)and XCO (d).

3.4 XCO2与XCH4相关性分析

拉萨地区XCO2,XCH4两种气体分子的日平均变化趋势相似(图5),为了分析两者是否具有相同的源,比较了两天中XCO2与XCH4随时间的变化趋势,并计算了观测期间两者的相关性,如图6和图7所示.

2018年8 月7 日在11:00~14:00(8日13:00~15:00)时段内,天空中有云,导致信号强度幅值变化大于10%,在预处理过程中将其剔除,因此图6在该时间段出现了数据缺失.从图6中可以看出XCO2,XCH4变化趋势相似,随着时间先增大后减小.将观测期间XCO2与XCH4结果进行了相关性分析 (见图7),其结果发现 8月 6日、7日、8日、13日、16日相关系数高于0.5,尤其是7日、8日和13日相关系数达到0.86左右,高相关系数说明 CO2与 CH4分子具有相同的源.图7(c)、图7(d)、图7(f)和图7(h)图呈现斜“V”字型,可能在观测期间存在源的输送,例如中午时段的汽车尾气,周边地区外来源的输送或者大气的垂直输送等,且拉萨地区8月份常出现反气旋天气现象,可能对气体的输送也存在一定的影响,后续将进一步研究反气旋天气与气体的输送关系.本文初步结合观测期间拉萨地区空气流输送轨迹分析CO2与CH4的源,利用NOAA开发的轨迹模型HYSPLIT分析后向72 h气体分子输送轨迹,如图8所示.

从气体分子输送轨迹来看,3类轨迹的运输速度相当,8月6日大气气流的输送轨迹主要来自从东南方向,3类轨迹的起源均来自林芝地区,途经山南地区向西北方向输送到达拉萨;而山南市河流众多,水资源丰富,冰川蓄积量约 10亿立方米,山南市有山羊、绵羊等动物,且盛产青稞、小麦、玉米等作物,主要林木有杨、柳、落叶松等,这些均是CH4及CO2的源.8月7日和8月14日均存在垂直方向上的气流对流,但到达拉萨观测点已经沉降到地面,成为了局地源,8月7—16日气体的输送源主要来自当地的局地源(除了10—12日期间高层大气运输来自北边那曲地区).

图6 2018 年 8 月 7 日、8 日 XCO2 (a),(b),XCH4 (c),(d) 时间序列Fig.6.Time series of XCO2 (a),(b)and XCH4 (c),(d)at August 7,8,2018..

图7 XCO2 与 XCH4 相关性 (a) 2018-08-6;(b) 2018-08-7;(c) 2018-08-8;(d) 2018-08-10;(e) 2018-08-12;(f) 2018-08-13;(g) 2018-08-14;(h) 2018-08-16Fig.7.The correlation between XCO2 and XCH4: (a) August 6,2018;(b) August 7,2018;(c) August 8,2018;(d) August 10,2018.

3.5 地基观测结果与WACCM数据对比分析

美国大气研究中心发展的全球大气模型WACCM(whole atmosphere community climate model)可预测从地面到热层顶大气时空演化特性.在研究痕量气体浓度方面,通常用于分析其季节变化趋势[20],由于地基EM27/SUN的观测时间短,在观测期间没有直观的观测数据(探空气球、卫星数据)进行对比,只能将其与WACCM模型数据进行对比.

本文数据以拉萨为中心(经度91.135°E,纬度29.659°N),取经纬度 ± 0.5°内的数据作为拉萨地区样本点.图9是 WACCM模式数据与地基EM27观测结果的对比图,从图9可以看出,WACCM模式获得的拉萨地区XCO2,XCH4平均值分别为 407.76 ppmv 和 1.883 ppmv;相比较WACCM模式数据,拉萨地基观测的XCO2值平均低 1.68 ppmv,XCH4平均低 0.193 ppmv,两者之间的相对偏差为0.4%(XCO2)和10.24%(XCH4).推测两者差异较大的原因主要有两点: 一是地基观测是通过反演算法获得柱浓度,存在反演误差;二是WACCM模式基于GEOS-5模型获得的模拟数据,也会导致两者结果具有差异.

图8 2018 年 8 月 6 日—16 日 72 h 后向轨迹图 (a) 2018 年 8 月 4—6 日;(b) 2018 年 8 月 7—9 日;(c) 2018 年 8 月 10—12 日;(d) 2018年 8月 14—16日Fig.8.72-hour back trajectories of Lhasa during August 6—16,2018: (a) August 4-6,2018;(b) August 7-9,2018;(c) August 10-12,2018;(d) August 14-16,2018.

图9 地基观测 XCO2 (a),XCH4 (b) 日平均值与 WACCM 数据对比Fig.9.Comparison of XCO2 (a) and XCH4 (b) based on ground-based observations and WACCM data.

4 结 论

基于地基便携式傅里叶变换光谱仪EM27/SUN,获得了拉萨地区2018年8月份(观测时间15天)的太阳吸收光谱,反演得到 2018年8月6 日到 8 月 17 日大气中的 H2O,CO2,CH4及 CO浓度信息.理论分析了对流层顶高对气体分子垂直分布的影响,并定量分析了对流层顶高度变化对气体柱平均摩尔分数的影响,分析了拉萨地区XH2O,XCO2,XCH4及XCO 的变化趋势,并将其结果与WACCM模式数据进行了对比分析,以期为研究我国高原温带地区温室气体浓度的时空分布及其变化规律提供数据支持及理论依据,主要结果如下:

1)对流层顶高度的变化对不同气体浓度反演的影响程度不同.H2O浓度与对流层顶高度变化没有明显关系;对流层顶高度变化对CO2,CH4及CO 浓度影响较大,XCO2,XCH4与对流层顶高度变化呈正相关,XCO与对流层顶高度变化呈负相关,XCO2,XCH4及XCO 与对流层顶高度变化的相关性系数分别为0.998,0.78及0.994.当对流层顶高度变化 3 km，DXCO2、DXCH4及 DXCO 分别为8.640%、0.035%及0.049%。

2)在观测期间,拉萨地区XH2O,XCO2及XCH4日变化幅度较大,且XCO2,XCH4变化趋势相似,多数天相关系数高于 0.5,说明 CO2,CH4具有相似的源,观测期间空气流主要来自山南地区的输送;XH2O,XCO2,XCH4及XCO 日平均值分别在 3432—4287 ppmv,406.1—408.2 ppmv,1.673—1.720 ppmv以及 0.082—0.095 ppmv之间变化;XH2O,XCO2,XCH4及XCO 的 平 均 值 分 别 为3919.70,406.87,1.689 及 0.091 ppmv.

3)拉萨地区WACCM模式获得的XCO2,XCH4与地基观测值存在偏差,相对偏差分别为0.4%(XCO2) 和 10.24%(XCH4).相比较 WACCM 模拟值,地基观测的XCO2,XCH4结果均偏小,地基XCO2平均低 1.68 ppmv,XCH4平 均 低 0.193 ppmv.WACCM模式数据可能需要考虑海拔高度的影响.