谢劭峰，王义杰，黄良珂，彭 华，黎峻宇，刘立龙
1.桂林理工大学测绘地理信息学院，广西 桂林 541004；
2.广西空间信息与测绘重点实验室，广西 桂林 541004
水汽是地球大气的重要组成成分，虽然水汽在大气中占比很小，但却是大气中最活跃的部分，它在大气辐射、水汽循环和气候变化等方面都起着至关重要的作用[1-4]。大气加权平均温度(Tm)是全球导航卫星系统(GNSS)对流层水汽反演过程中计算水汽转换系数的关键参数[5]。近年来，利用大气再分析资料(如ERA-Interim、ERA5、NCEP和MERRA-2等)积分计算的Tm信息通过空间插值应用于GNSS大气水汽反演获得了广泛关注[6-7]。因为水汽在不同大气高度范围变化趋势不一致，所以获取不同高度处的高精度Tm信息十分重要。由于Tm在高程上的变化远大于水平方向上的变化[8-9]，而大气再分析资料的格网点高程通常与GNSS站、探空站等用户所处高程不一致(排除不同数据源之间的高程基准差异)，尤其在地形起伏较大的中国区域，这种高程差异更为显著。因此，研究中国区域高精度的Tm垂直递减率模型能够将再分析资料格网点的Tm数据精确插值到用户位置，同时，它也是构建中国区域高精度Tm模型的重要基础。
目前，常用的Tm模型主要分为两大类：需要实测气象参数的模型和非气象参数模型。对于前者，最具代表性的是Bevis模型[5]，后续诸多学者对Bevis模型在区域或全球范围内进行了精化[10-13]，但这些Tm模型依赖实测气象参数，且未考虑Bevis模型在不同高度处的适用性，从而限制了它在实时GNSS水汽监测中的应用。针对上述模型的不足，诸多学者建立了许多区域或全球范围的非气象参数Tm模型[14-16]，极大地促进了GNSS气象学的发展。文献[17]利用全球135个探空站5 a的数据，建立了一个全球加权平均温度模型(GWMT模型)，可实现全球实时的GNSS-PWV反演。文献[18]构建了GPT2w模型，该模型可在全球范围内提供高精度的Tm信息，但是它未顾及Tm在高程上的改正，在地形起伏较大的区域使用时，存在显著的系统偏差。针对GPT2w模型的不足，文献[19]建立了中国区域顾及Tm高程变化的改进模型(IGPT2w)，IGPT2w模型在中国区域表现出明显的精度提升。文献[20]建立了中国区域统一的Tm垂直递减率模型，在此基础上提出了顾及垂直递减率函数的中国区域大气加权平均温度模型(CTm模型)，它在地形起伏较大的中国西部地区具有显著的精度改善。文献[21]进一步研究表明Tm垂直递减率是Tm高程改正的重要参数，并利用Tm垂直递减率将再分析资料格网点高度处的Tm值改正到目标高度处，最终构建了中国地区增强Tm模型(GH-Tm模型)。文献[22]在全球范围内利用同一个Tm垂直递减率进行高程改正，建立了全球Tm模型(GTm-Ⅲ模型)。文献[23]使用欧洲中尺度天气预报中心(European Centre for Medium-Range Weather Forecasts,ECMWF)的产品对Tm在垂直方向上的变化特性进行了分析,并给出每个格网点上非线性的Tm高程改正函数，建立了一种全球Tm模型(GTm-H模型)，在全球范围内都具有较高的精度。
尽管现有的Tm模型顾及了Tm在高程上的改正且表现出了各自的优越性，但仍存在构建Tm垂直递减率模型时仅使用单一格网点数据、建模数据使用月均剖面资料等不足。同时，已有的经验Tm模型提供的Tm信息难以精确捕获Tm的日周期变化，而中国区域具有地域辽阔、地形起伏较大和气候多样化等特点，尤其在西部地区地形地貌和气候条件更为复杂。因此，亟须建立一个高精度的中国区域Tm垂直递减率模型，以提高大气再分析资料的高时间分辨率Tm格网产品的插值精度，并提供实时、高精度的Tm高程改正值。为此，本文拟引入滑动窗口算法，并利用多年的MERRA-2(The second modern-era retrospective analysis for research and applications)大气再分析资料，在分析Tm垂直递减率精细季节变化的基础上，建立顾及时变垂直递减率的中国区域实时高精度Tm垂直递减率格网模型，并利用MERRA-2、GGOS大气格网数据和探空站资料检验模型的精度。
MERRA-2大气再分析资料(https:∥goldsmr4.gesdisc.eosdis.nasa.gov/data/MERRA2)是美国国家航空航天局(National Aeronautics and Space Administration,NASA)推出的最新一代高时空分辨率的大气再分析资料，可提供分层气象参数和地表气象参数，其水平分辨率为0.5°×0.625°(纬度×经度)，地表气象参数时间分辨率为1 h，分层气象参数时间分辨率为6 h，垂直分辨率为42层(顶层高度约为50 km)，根据气象参数进行积分计算可得到相应的分层和地表Tm格网数据。Global Geodetic Observing System(GGOS)Atmo-sphere(https:∥vmf.geo.tuwien.ac.at/)可提供时间分辨率为6 h、水平分辨率为2°×2.5°(纬度×经度)的Tm大气格网产品及相应的地表高程，其Tm大气格网产品是基于ERA-Interim大气再分析资料积分计算得到。无线电探空资料(http:∥weather.uwyo.edu/upperair/sounding.html)可提供从地面到约30km高度范围内实测的分层气象参数及相应地表气象参数，时间分辨率为12 h，因此，无线电探空资料常作为独立观测值用来评价其他再分析资料产品和模型的精度。其中，中国区域89个探空站点位分布如图1所示。
图1 中国区域89个探空站和4个代表性格网点分布
Tm可利用探空站的分层资料或者大气再分析分层资料积分计算获取，其积分计算表达式为
(1)
式中，e为水汽压；
T为温度；
ei为第i层大气的平均水汽压；
Ti为第i层大气的平均温度；
ΔHi为第i层大气的厚度；
ht、hb分别为积分的最顶层高度和最底层高度。
由于大气再分析资料通常在使用前需要对其进行精度评定，为此本文利用2017年中国区域89个探空站的数据对MERRA-2大气再分析资料计算Tm的精度进行评估，并统计MERRA-2资料积分计算的Tm值的偏差(bias)和均方根误差(root mean square,RMS)，统计结果见表1。
表1 2017年中国区域探空站数据检验MERRA-2再分析资料计算Tm的精度
由表1可知，MERRA-2再分析资料在中国区域计算Tm的平均bias为0.10 K，平均RMS为1.34 K。因此，MERRA-2再分析资料在中国区域计算Tm具有较高的精度，可作为构建中国区域Tm垂直递减率模型的数据源。
2.1 Tm垂直递减率时空特性分析
中国区域地势西高东低，总体呈三级阶梯状分布，地形起伏较大，导致再分析资料格网点高程与用户高程存在较大差异，因此，探究Tm垂直递减率在中国区域近地空间范围内(指0～10 km的高程范围)的时空变化特性具有重要意义[24]。文献[21—22]研究表明Tm与高度和纬度具有较大的相关性，为了更好地分析Tm随高程的变化关系，本文选取了中国区域4个具有代表性的MERRA-2再分析资料格网点(45°N,120°E；
42°N,105°E；
42°N,90°E；
30°N,105°E)2016年1月1日的分层Tm数据及相应的位势高参数进行分析，其中分层Tm数据是利用式(1)将气压从1000 hPa到250 hPa的气压层(高程范围约为0～10 km)作为底层往顶层积分得到的对应高度面的Tm，结果如图2所示。
由图2可知，中国区域4个格网点的Tm与高程均存在近似的线性变化关系，其变化关系可表示为
图2 中国区域2016年1月1日4个格网点Tm随高程的变化
Tm=β×φh+l
(2)
式中，β表示Tm垂直递减率，单位为K/km；
φh表示椭球高，单位为km。由此可知，β是Tm高程改正的重要参数。文献[20]研究表明β具有明显的季节性变化，为了进一步分析β的精细季节变化，根据式(2)计算出每个格网点2012—2016年的日均β时间序列，并选取上述4个格网点采用年周期和半年周期的余弦函数对β进行拟合，结果如图3所示。
图3 中国区域4个格网点Tm垂直递减率时序变化
由图3可知，β在中国区域表现出明显的年周期和半年周期变化，其周期性变化可能与大气垂直结构的季节性变化有关[25-26]，此外，还发现低纬度格网点β的季节性变化没有其他地区明显，可能是受亚热带季风气候和热带海洋气候等复杂气候的综合影响。由此可知，β在不同地区的变化规律不一致。为了进一步分析β的年均值、年周期和半年周期振幅在中国区域的分布特性，计算了2012—2016年中国区域MERRA×2再分析资料每个格网点β的年均值、年周期和半年周期振幅，结果如图4所示。
图4 中国区域β的年均值、年周期振幅和半年周期振幅分布
由图4可知，β在西部地区具有绝对值相对较大的年均值，主要原因可能是该地区地形起伏较大，Tm随着高度变化剧烈。而β在高纬度地区具有较大的年周期振幅，在西部高海拔地区具有相对较大的半年周期振幅，主要原因可能是这些地区Tm季节性变化更强。因此，为了保证Tm垂直递减率模型的精度，在模型构建时需同时顾及β的年周期和半年周期变化。
2.2 CTm-H模型建立
针对已有Tm垂直递减率模型建模时仅采用单一格网点数据以及使用月均剖面数据等不足，本文引入滑动窗口算法将中国区域剖分为大小一致的规则窗口。首先对中国区域进行规则格网剖分，得到与MERRA-2再分析资料大小相同的水平分辨率格网，即0.5°×0.625°(纬度×经度)。滑动窗口算法的关键需要确定其窗口的大小及步长，同时窗口大小的确定需顾及中国区域窗口剖分个数的整数性、窗口的连续性及窗口内模型参数的可求解性等原则。基于此，为了探究窗口大小对模型精度的影响，本文选取了几种不同的窗口分辨率进行建模，分别为1°×1.25°、2°×2.5°和4°×5°(纬度×经度)。滑动窗口算法的具体流程详见文献[27]。
本文利用中国区域各窗口内所有格网点数据,在分析Tm垂直递减率精细季节变化的基础上，建立顾及时变垂直递减率的实时高精度Tm垂直递减率格网模型，其表达式如下
(3)
(4)
针对中国区域每个窗口，利用窗口内9、25、81个MERRA-2格网点2012—2016年6 h分辨率的Tm分层数据及相应位势高参数，通过最小二乘法，估计出中国区域每个窗口的Tm垂直递减率模型的5个系数，并存储于分辨率为1°×1.25°、2°×2.5°和4°×5°窗口的几何中心，最终构建了中国区域水平分辨率为1°×1.25°、2°×2.5°和4°×5°的Tm垂直递减率格网模型(简称“CTm-H1、CTm-H2和CTm-H3模型”)。模型的使用非常便捷：①用户仅需提供年积日和目标点位置信息，根据目标点位置信息查找其所处窗口的模型参数；
②根据查询获得的模型参数，利用式(3)和式(4)即可将目标点在参考高程处的Tm值改正到目标高程处。
本文选取未参与建模的2017年的中国区域MERRA-2、GGOS大气格网数据和89个探空站资料为参考值，并以bias和RMS作为精度指标[19,27-28]，评价CTm-H模型在中国区域的精度。
3.1 CTm-H模型在垂直插值中的精度分析
为验证CTm-H模型在中国区域的垂直插值精度和适用性，同时探究不同格网分辨率数据对模型在垂直插值中精度的影响，为此，对MERRA-2格网数据进行提取，得到分辨率分别为1°×1.25°、2°×2.5°和4°×5°(纬度×经度)的MERRA-2格网数据。在本次模型垂直插值检验中，以2017年3种分辨率的MERRA-2格网数据为参考值，通过CTm-H模型，将MERRA-2地表格网数据改正到分层格网数据在近地空间范围内各层高度处进行检验，最后对每个格网点Tm的垂直插值的bias和RMS进行统计，并与中国区域统一的Tm垂直递减率模型[20](简称“统一模型”)进行比较，结果见表2。
表2 利用2017年MERRA-2格网数据检验CTm-H和统一模型的精度统计
由表2可知，从格网分辨率由低到高来看，各模型的bias和RMS相差不大，主要原因是CTm-H和统一模型是对格网点Tm进行垂直插值，没有进行水平方向的插值。因此，以格网分辨率为1°×1.25°为例，对各模型进行分析。其中，统一模型的bias和RMS的变化范围最大，其平均bias和RMS分别为-3.32 K和4.78 K；
而CTm-H 3个模型的bias和RMS的变化范围较小，表现相似，其平均bias和RMS分别在-0.40 K和3.32 K左右。相比于统一模型，CTm-H模型的精度(RMS)提高了30%，主要原因是CTm-H模型顾及了Tm垂直递减率在不同地区的变化规律，在使用时能够较好地改正Tm在垂直方向上的变化，且当两种Tm数据高程相差较大(高差最大可达到10 km左右)时，CTm-H模型改正效果更加明显，表现出显著的优势。
以格网分辨率为1°×1.25°为例，各模型的bias和RMS在中国区域的分布如图5所示。由图5可知，统一模型在西部地区表现较好，具有绝对值相对较小的负bias和RMS，而在东南及东北地区表现出绝对值较大的负bias和RMS，主要原因是在中国区域使用统一的Tm垂直递减率引入了较大的系统误差。相比之下，CTm-H1、CTm-H2和CTm-H3模型在中国区域均表现出较小的bias和RMS，CTm-H 3个模型的性能相当。CTm-H模型在东南地区表现出绝对值较小的负bias，其他地区表现出较小的正bias，说明在东南地区CTm-H模型提供的Tm改正值较参考值略小，其他地区提供的Tm改正值较参考值略大；
在RMS方面，CTm-H模型在中国区域均表现出较小的RMS，尤其在西部地区，表现出更为显著的优势。总体而言，引入滑动窗口算法进行建模可以有效克服复杂地形因素带来的系统性误差，能够显著提升模型的垂直插值精度及稳定性。
图5 利用2017年MERRA-2格网数据检验CTm-H和统一模型的bias和RMS分布(上为bias，下为RMS)
3.2 CTm-H模型在空间插值中的精度分析
3.2.1 模型在MERRA-2地表格网数据空间插值中的精度检验
探空站提供的Tm信息均为实测值，具有较高的精度和可靠性。为验证CTm-H模型在中国区域的空间插值精度和适用性，同时探究不同格网分辨率数据对模型在空间插值中精度的影响，为此，对MERRA-2地表格网数据进行提取，得到分辨率分别为1°×1.25°、2°×2.5°和4°×5°的MERRA-2地表格网数据。在本次模型空间插值检验中，以2017年中国区域89个探空站数据为参考值，通过CTm-H模型将3种分辨率MERRA-2地表格网数据中探空站周围的4个MERRA-2格网点Tm改正到探空站高度处，再将高程改正后的4个格网点Tm利用反距离加权法插值到探空站处进行检验，最后对每个探空站Tm的空间插值的bias和RMS进行统计，并与统一模型和未顾及垂直改正的情况进行比较，结果如表3、图6和图7所示。
图6 利用2017年探空站和MERRA-2地表格网数据检验不同垂直改正模型的bias分布
图7 利用2017年探空站和MERRA-2地表格网数据检验不同垂直改正模型的RMS分布
表3 利用2017年探空站和MERRA-2地表格网数据检验不同垂直改正模型的精度统计
由表3可知，从格网分辨率由低到高来看，各模型的bias和RMS的变化范围和平均绝对值逐渐减小，主要原因是在空间插值中，格网分辨率越高插值精度越高。因此，以格网分辨率为1°×1.25°为例进行分析。在未顾及垂直改正时，bias的变化范围最大，其大小变化范围为-7.56～1.12 K，表现出明显且绝对值较大的负bias，平均bias为-1.84 K，而经过CTm-H和统一模型改正后bias的变化范围较小，基本保持在±2 K以内，其平均bias均为较小的正bias；
在RMS方面，在未顾及垂直改正时RMS的变化范围和平均值最大，其平均RMS为2.69 K，而经过CTm-H和统一模型改正后RMS的变化范围和平均值较小，其中CTm-H模型的平均RMS在1.83 K左右，统一模型的平均RMS为1.88 K。与统一模型相比，CTm-H模型的性能提升较小，精度只提高了3%，主要原因是以探空站数据为参考值时，在多数地区MERRA-2地表数据格网点的高程与探空站高程相差不大，Tm高程改正的空间较小，只需要利用简单地Tm垂直递减率模型进行高程改正即可达到较高的精度。但CTm-H和统一模型相比于未顾及垂直改正性能提升较大，精度提高了30%左右，说明Tm垂直递减率模型能够显著提高Tm在空间插值中的精度。
结合图6、图7和表3可知，随着格网分辨率的提高，各模型在中国西部地区出现显著异常值的情况也随之减少，当格网分辨率为1°×1.25°时，各模型表现最佳。未顾及垂直改正时，在西部地区表现出绝对值较大的负bias和RMS，在其他地区表现出绝对值相对较小的bias和RMS，主要原因是西部地区MERRA-2地表格网数据的格网点高程与探空站高程相差较大，在空间插值时没有对格网点Tm进行高程改正引入了较大的误差；
而经过CTm-H和统一模型改正后，模型的bias和RMS在中国区域均表现出较小的绝对值，CTm-H和统一模型整体表现相似，稳定性较高。相较于未顾及垂直改正，CTm-H和统一模型在中国东部地区改善不够显著，但在西部地区表现出显著的优势。总体而言，CTm-H和统一模型能够显著提高Tm在空间插值中的精度和稳定性。
为了分析CTm-H模型与统一模型和未顾及垂直改正的bias和RMS随高程的变化关系，以格网分辨率为1°×1.25°时各模型的bias和RMS为例对其进行统计，结果如图8所示。由图8可知，在未顾及垂直改正时，随着探空站所处高程逐渐增大，表现出较大的负bias和RMS，与前述在西部地区表现出较大的负bias和RMS的结果相同(海拔较高的探空站基本分布在西部地区)。而经过CTm-H和统一模型改正后随着探空站所处高程逐渐增大，Tm高程改正的精度逐渐提高，其bias均在0 K左右，RMS均在2 K左右，且在不同探空站高度处，CTm-H和统一模型均表现出较高的精度和稳定性。
图8 不同垂直改正模型的bias和RMS随探空站高程变化的分布
3.2.2 模型在GGOS大气格网产品空间插值中的精度检验
为进一步验证CTm-H模型在中国区域的空间插值精度和适用性，利用探空站数据检验模型在GGOS大气格网产品空间插值中的精度。GGOS大气格网产品的最高分辨率为2°×2.5°，由于GGOS大气格网产品的高程是大地高，而探空站的高程是海拔高，两者之间存在高程基准差异，因此，本文采用EGM2008模型来实现两者高程基准统一[28]。为了探究不同分辨率的GGOS大气格网产品对模型在空间插值中精度的影响，对GGOS大气格网产品进行提取，得到分辨率分别为2°×2.5°和4°×5°的GGOS大气格网产品。以2017年中国区域89个探空站数据为参考值，利用与模型在MERRA-2地表格网数据空间插值中相同的方法进行检验，然后对每个探空站Tm的空间插值的bias和RMS进行统计，并与统一模型和未顾及垂直改正的情况进行比较，结果如表4、图9和图10所示。
表4 利用2017年探空站和GGOS大气格网产品检验不同垂直改正模型的精度统计
图9 利用2017年探空站和GGOS大气格网产品检验不同垂直改正模型的bias分布
图10 利用2017年探空站和GGOS大气格网产品检验不同垂直改正模型的RMS分布
由表4可知，从格网分辨率由低到高来看，各模型在空间插值中的精度逐渐提高。因此，以格网分辨率为2°×2.5°为例进行分析。在未顾及垂直改正时，bias和RMS的变化范围最大，其平均bias和RMS分别为-2.19 K和3.16 K；
而经过CTm-H和统一模型改正后，bias和RMS的变化范围较小，表现出较小的平均bias和RMS，其中统一模型的平均bias和RMS分别为0.85 K和2.16 K，CTm-H模型的平均bias和RMS在0.83 K和2.05 K左右。与统一模型相比，CTm-H模型精度提高了5%，性能提升较小；
与未顾及垂直改正相比，CTm-H和统一模型精度提高了约46%，性能提升较大。总体而言，与各模型在MERRA-2地表格网数据空间插值中检验的结果相似。与表3中相应格网分辨率的MERRA-2地表格网数据检验结果相比，CTm-H和统一模型在GGOS大气格网产品空间插值中的平均bias相对较大，而平均RMS相差不大。
结合图9、图10和表4可知，随着格网分辨率的提高，各模型精度逐渐提高，且出现显著异常值的情况也随之减少，当格网分辨率为2°×2.5°时，各模型表现最佳。未顾及垂直改正时，在西部地区表现出绝对值较大的负bias和RMS；
经过CTm-H和统一模型改正后，其bias和RMS在中国区域均表现出较小的绝对值，整体表现相似，与模型在MERRA-2地表格网数据空间插值中检验的结果相似。由此进一步说明，CTm-H和统一模型能够显著提高Tm在空间插值中的精度和稳定性。
利用大气再分析资料获取的Tm信息可应用于高时间分辨率GNSS-PWV监测，但需使用高精度的Tm垂直递减率模型对其进行高程改正。针对已有Tm垂直递减率模型建模仅使用单一格网点数据等不足，本文引入滑动窗口算法建立了顾及时变垂直递减率的中国区域不同水平分辨率Tm垂直递减率格网模型(CTm-H1、CTm-H2和CTm-H3模型)。联合2017年的MERRA-2、GGOS大气格网数据和探空站资料，对CTm-H模型进行精度检验，并与统一模型进行精度对比。结果表明：以MERRA-2格网数据为参考值，CTm-H 3个模型的性能相当，相比于统一模型，CTm-H模型性能提升明显且在中国区域内均表现出较高的精度和稳定性；
以探空站数据为参考值，CTm-H 3个模型与统一模型的性能相似，但与未顾及垂直改正相比，CTm-H和统一模型性能提升明显，尤其西部地区表现出显著的优势。CTm-H 3个模型在中国区域均具有较高的精度，不需要实测气象参数即可提供中国区域近地空间范围内任意位置实时高精度的Tm高程改正，因此，在中国区域的实时高精度GNSS水汽探测中具有重要的应用。
由于CTm-H 3个模型精度相当，其中CTm-H3模型参数最少，使用更加便捷，建议使用CTm-H3模型进行中国区域的Tm高程改正。本文在建模时只利用了MERRA-2再分析资料构建了中国区域近地空间范围内的Tm垂直递减率格网模型，后续研究将采用更高分辨率的ERA5再分析资料构建全球范围的Tm垂直递减率格网模型以及其在GNSS-PWV中的应用。
致谢：感谢NASA提供的MERRA-2再分析资料、GGOS Atmosphere提供的Tm格网数据和美国怀俄明大学提供的无线电探空资料。
猜你喜欢递减率探空格网用L波段探空测风雷达评估风廓线雷达测风准确性成都信息工程大学学报(2021年3期)2021-11-22实时电离层格网数据精度评估空间科学学报(2020年6期)2020-07-21水驱油藏产油量递减率影响因素的分解方法——以锦州9–3油田和渤海油田的应用为例石油地质与工程(2019年2期)2019-06-19一种水驱油藏产油量递减率评价新方法天然气与石油(2019年1期)2019-03-20上市储量评估中指数递减法合理递减率的选取特种油气藏(2018年6期)2019-01-11TK-2GPS人影火箭探空数据与L波段探空数据对比分析高原山地气象研究(2016年1期)2016-11-10延长油田递减率计算方法研究延安大学学报（自然科学版）(2015年3期)2015-06-07浅谈净举力对探空气球升速及施放高度的影响河南科技(2014年14期)2014-02-27郑州探空数据库的建设简介河南科技(2014年4期)2014-02-27平均Helmert空间重力异常格网构制方法测绘科学与工程(2014年4期)2014-02-27