2019, Vol. 44文章信息
- 陈俊平, 王君刚, 王解先, 谭伟杰
- CHEN Junping, WANG Jungang, WANG Jiexian, TAN Weijie
- SHAtrop：基于陆态网GNSS数据的中国大陆区域ZTD模型
- SHAtrop: Empirical ZTD Model Based on CMONOC GNSS Network
- 武汉大学学报·信息科学版, 2019, 44(11): 1588-1595
- Geomatics and Information Science of Wuhan University, 2019, 44(11): 1588-1595
- http://dx.doi.org/10.13203/j.whugis20170384
-
文章历史
收稿日期: 2018-06-07
2. 中国科学院大学天文与空间科学学院, 北京, 100049;
3. 同济大学测绘与地理信息学院, 上海, 200092
2. School of Astronomy and Space Science, Chinese Academy of Sciences, Beijing 100049, China;
3. College of Surveying and Geo-Informatics, Tongji University, Shanghai 200092, China
大气延迟是全球导航卫星系统(Global Navigation Satellite System, GNSS)数据处理的主要误差源之一，主要包括对流层延迟与电离层延迟，对流层延迟通常作为待估参数解算或通过模型改正。对流层天顶延迟(zenith total delay, ZTD)约为2.3 m，在5°高度角时由投影函数映射到斜路径上的延迟量可达20 m[1]，GNSS常用的映射函数有NMF[2]、GMF[3]以及VMF1[4]等。
采用实测气象数值模型积分计算ZTD精度较高，但是计算复杂；Saastamoinen、Hopfield等模型在采用实测气象数据时也可以取得较高精度，但采用标准大气时其精度有限。EGNOS模型[5]是欧盟星基广域增强系统所采用的对流层延迟改正模型。UNB3是针对美国的广域增强导航系统(wide area augmentation system, WAAS)建立的对流层天顶延迟模型，在北美地区精度较好[6]；UNB3m改进了UNB3模型，其参数中采用相对湿度代替UNB3模型的水汽压, 改善了湿延迟的精度, 其偏差（Bias）为－0.5 cm[7]。GPT、GPT2[8-9]模型是基于欧洲中期天气预报中心(European Centre for Medium-Range Weather Forecasts, ECMWF)建立的简化气象模型，利用其气象参数与Saastamoinen模型计算的对流层天顶延迟精度优于EGNOS、UNB3模型。文献[10]根据美国国家环境预报中心（National Centers for Environmental Prediction，NCEP）再分析资料，采用三维参数表构建了IGGtrop全球对流层天顶延迟改正模型，全球平均偏差与均方根误差(root mean square, RMS)分别为－0.8 cm、4.0 cm，优于EGNOS、UNB3、UNB3m模型结果。
近年来，部分学者利用GNSS实测对流层延迟数据，建立了与气象参数无关的对流层天顶延迟经验改正模型，取得了良好的效果。文献[11]采用中国区域的ECMWF分层资料建立了区域对流层格网模型SHAO-C，其精度(RMS)达到4.5 cm，分辨率高于EGNOS。文献[12]利用NCEP全球再分析资料，把SHAO模型扩展到全球范围，初步建立了全球对流层延迟改正模型SHAO-G，其平均Bias为-1.2 cm，RMS为4.1 cm。文献[13]分析了GGOS Atmosphere提供的全球天顶对流层延迟格网数据，采用球谐函数构建了全球非气象参数对流层天顶延迟模型GZTD，其精度(Bias: -0.02 cm, RMS: 4.24 cm)优于EGNOS和UNB3系列模型。文献[14]采用国际GPS服务（International GPS Service，IGS)提供的对流层天顶延迟产品，建立了一种简单的全球对流层天顶延迟模型，该模型在全球范围内的RMS为4.9 cm。
本文利用中国大陆构造环境监测网络(Crustal Movement Observation Network of China, CMONOC)的GNSS测站的对流层天顶延迟，分析了中国区域223个GNSS测站的对流层天顶延迟时空特性，根据其周期变化特征的参数，建立了基于区域格网的适用于中国区域的天顶对流层延迟模型SHAtrop。SHAtrop充分考虑了不同经纬度、不同高程地区的对流层特征，建模方法简单，使用方便，且精度良好。新模型精度(RMS)为3.4 cm，能满足实时导航与定位的精度要求，可以作为中国区域GNSS用户的对流层参考模型。
1 建模数据上海天文台GNSS分析中心(GNSS Analysis Center at Shanghai Astronomical Observatory, SHA)对IGS测站以及陆态网约300个GNSS测站的观测值进行处理，提供精密轨道和钟差等产品[15-16]，其精密GNSS轨道和钟差与IGS一致。SHA在轨道钟差解算时也估计对流层天顶延迟[17]，其以GPT模型计算的天顶延迟为初值，湿延迟作为分段常数每小时估计一次。为评估SHA ZTD的精度，将SHA ZTD与IGS提供的ZTD进行对比。统计了2012年1月至2017年10月期间中国区域的7个IGS测站的SHA ZTD ZiIGS 与IGS ZTD ZiSHA 的平均偏差B和RMS（R），统计方法见式(1)。由于IGS ZTD采样间隔为5 min，统计时通过取均值的方法将其重采样至1 h。中国区域的7个IGS测站的平均偏差与RMS结果见表 1。
| 精度指标 | BJFS | SHAO | CHAN | LHAZ | URUM | KUNM | TWTF | 平均 |
| 平均偏差 | 1.2 | -0.8 | 1.6 | 1.5 | 1 | 1 | -0.7 | 0.7 |
| RMS | 2.2 | 2 | 2.5 | 2.7 | 1.6 | 3.8 | 2.4 | 2.5 |
| $ \left\{ {\begin{array}{*{20}{l}} {B = \frac{1}{N}\sum\limits_{i = 1}^N {\left( {Z_i^{{\rm{IGS}}} - Z_i^{{\rm{SHA}}}} \right)} }\\ {R = \sqrt {\frac{1}{N}\sum\limits_{i = 1}^N {{{\left( {Z_i^{{\rm{IGS}}} - Z_i^{{\rm{SHA}}}} \right)}^2}} } } \end{array}} \right. $ | (1) |
从表 1中可以看出，中国区域测站的SHA ZTD与IGS ZTD一致性很好，各测站平均偏差均值为0.7 mm，RMS均值为2.5 mm。以上分析表明，SHA解算的ZTD具有与IGS提供的ZTD相同的精度。
SHA提供的陆态网ZTD数据包含了全部260个GNSS测站。选用其中观测较长的223个测站2012年1月至2017年10月约6 a的结果，分析了其时空分布特性，并建立区域高精度经验模型。其中，152个测站的数据用于建模，71个测站的数据用于模型的精度检验，测站分布与数据可用性见图 1。
|
| 图 1 GNSS测站分布和每个GNSS测站ZTD数据的可用天数 Fig. 1 GNSS Site Distribution and Days When ZTD Is Available for Each Site |
图 1（a）中红色测站用来建模，蓝色测站用来检验精度；图 1（b）为各测站的数据时间可用性。可以看出，绝大多数测站的数据可用天数为1 500 d，最少可用天数也达2 a，保证了建模与精度检验的可靠性。
2 中国区域ZTD时空特性分析 2.1 ZTD随测站高程的变化首先计算223个测站2012—2017年6 a的ZTD均值，测站高程与ZTD均值的关系见图 2。本文中的测站高程都是基于椭球面的大地高。由图 2可见，ZTD与测站高程近似呈反比关系，随着测站高程的增加，ZTD不断减小。
|
| 图 2 GNSS测站的高程和每个测站的ZTD年均值 Fig. 2 GNSS Site Ellipsoid Height and Mean ZTD of Each Site |
ZTD Z(h)与测站高程的关系常用指数模型或线性模型表示，指数函数为：
| $ Z\left( h \right) = {Z_0} \times {{\rm{e}}^{\beta h}} $ | (2) |
式中，Z0为测站在椭球面的总延迟量；e为自然常数；h为测站高程（m）；β为ZTD随高程衰减系数。
考虑到ZTD随高程的衰减系数与测站纬度相关，本文在纬度方向以5°为间隔，将大陆区域分成5个纬度区间，不同纬度区间分别采用不同系数进行拟合，拟合结果见图 3，拟合的参数见表 2。由图 3可以看出，指数函数能够合理地反映ZTD与测站高程的变化，而在不同的纬度区间拟合的ZTD高程衰减系数也呈现一定的规律性：随着纬度增高，参数β呈递增趋势（见表 2）。
|
| 图 3 GNSS测站ZTD与测站椭球面高程的关系 Fig. 3 ZTD Fitting Result with Respect to Site Ellipsoid Height |
| 纬度 | β/10－4 |
| < 25°N | －1.411 |
| 25°N~30°N | －1.357 |
| 30°N~35°N | －1.329 |
| 35°N~40°N | －1.303 |
| > 40°N | －1.300 |
由式（2）可知，ZTD的大小与测站高程相关，但是其时间变化规律（年周期等）与测站高程无关。图 4为青海门源（QHME）、山东泰安（TAIN）两个测站的时间序列以及傅里叶频谱分析结果。从图 4中可以看出，两个测站的ZTD都存在显著的周期信号。频谱结果显示周年信号显著，振幅都超过了5 cm；半周年信号也较为明显，其振幅也有厘米量级。因此，采用周期函数对年周期和半年周期信号Z0(d)进行建模：
|
| 图 4 天顶延迟时间序列与傅里叶变化结果 Fig. 4 ZTD Time Series and Fourier Transformation |
| $ {Z_0}\left( d \right) = {A_0} + {A_1}{\rm{cos}}\left( {\frac{{2{\rm{\pi }}}}{{365.25}}\left( {d - {d_1}} \right)} \right)$ | (3) |
| $ \begin{array}{*{20}{c}} {{Z_0}\left( d \right) = {A_0} + {A_1}{\rm{cos}}\left( {\frac{{2{\rm{\pi }}}}{{365.25}}\left( {d - {d_1}} \right)} \right) + }\\ {{A_2}{\rm{cos}}\left( {\frac{{4{\rm{\pi }}}}{{365.25}}\left( {d - {d_2}} \right)} \right)} \end{array}$ | (4) |
式中，A0是ZTD年均值；(A1，d1)是年周期项振幅与相位；(A2，d2)是半年周期项的振幅与相位；d为年积日。
对QHME、TAIN的ZTD按照年周期以及年周期+半年周期进行拟合，结果也绘制在图 4中。从图 4中可看到，年周期+半年周期的效果优于只采用年周期拟合的结果。
进一步统计各个测站两种拟合方式的RMS差，结果见图 5。
|
| 图 5 年周期+半年周期拟合模型相对于仅考虑周年项模型RMS的改进量 Fig. 5 Fitting RMS Differences Between Annual+Semi‑annual Model and Annual Only Model |
图 5中RMS差值序列全为正, 说明所有测站附加半年周期的拟合模型的精度高于仅考虑周年项的模型；且在东北、华北区域改善尤为明显，可达5~8 mm；附加半年周期拟合的改善效果与测站经纬度、高程没有明显的关系，主要表现为地域特性。
3 区域对流层模型SHAtrop的建立利用式(4)对陆态网223个GNSS测站6 a的ZTD时间序列进行拟合。各测站的拟合中误差与拟合参数见图 6。
|
| 图 6 ZTD周期模型的拟合参数 Fig. 6 ZTD Fitting Parameters Using Annual+Semi-annual Model |
图 6中年周期与半年周期相位的单位是d。由图 6可以看出，在水汽较为活跃的地区（华东华南地区），拟合中误差远大于其他地区，拟合中误差最大可达7 cm，最小为1 cm；对于拟合得到的参数，年周期项的振幅表现为西北地区较小，东南地区较大；半年周期的振幅表现为东北地区大于其他区域；并且拟合函数的常数项、周期项的振幅与初相位分布呈地域特性。
利用以上拟合结果，建立适用于中国大陆区域的ZTD经验改正模型SHAtrop。建模步骤为：（1）采用指数改正模型以及表 2中的ZTD高度衰减参数，对中国区域152个测站（图 1）6 a的ZTD时间序列进行高程改正，并将其归化到椭球面；（2）对每个测站归算后的ZTD时间序列进行年周期+半年周期的三角函数拟合，得到每个测站ZTD的常数项、年周期项和半年周期项（A0，A1，d1，A2，d2）；（3）考虑到各拟合参数呈现地域分布的特性，将陆态网络GNSS测站覆盖范围（70°E~135°E, 18°N~54°N）划分为2.5°×2.0°的经纬格网，根据每个测站拟合的常数项与周期项参数，采用双线性内插方法计算每个格网点的常数项与周期项参数。
SHAtrop的格网参数提供给用户使用时，用户站计算ZTD的过程为：
1）由测站经纬度确定所在格网，根据4个格网点的参数值双线性内插出用户站参数（A0, A1, d1, A2, d2）；
2）采用式(4)，计算当前年积日用户站在椭球面上的ZTD参数；
3）根据所在纬度区间，查询表 2中的高程改正系数，采用式(2)将ZTD由椭球面上改正到测站实际高程。
4 模型精度评估为分析SHAtrop的适用性，采用参与建模的152个测站（图 1中红色测站）与未参与建模的71个测站（图 1中蓝色测站），利用SHAtrop模型计算了所有测站6 a内的天顶对流层延迟。精度评定结果统计时，本文同时计算了EGNOS、UNB3m、GPT2模型的对流层天顶延迟，并以SHA解算的ZTD为参考值，计算各个模型的Bias与RMS，计算公式同式(1)。各模型的RMS的分布见图 7，偏差分布见图 8。
|
| 图 7 每个测站不同ZTD模型的RMS Fig. 7 RMS of the Different Models for Each Site |
|
| 图 8 每个测站不同ZTD模型的系统偏差 Fig. 8 Bias of the Different Models for Each Site |
由图 7可知，EGNOS与UNB3m模型的精度相当，且地域分布也类似，在东南部精度低于西北地区；GPT2模型中误差整体小于EGNOS与UNB3m，但是在四川区域内RMS较大，最大可达9 cm；SHAtrop模型精度普遍好于前3种模型，并且SHAtrop模型的Bias分布较均匀，不存在偏差较大的区域。
不同ZTD模型RMS与Bias的统计见表 3。
| 测站 | 精度指标 | EGNOS | UNB3m | GPT2 | SHAtrop |
| 建模站 | RMS | 5.6(1.8, 8.6) | 5.2(1.7, 8.2) | 4.7(1.6, 9.8) | 3.4(1.6, 7.0) |
| Bias | －0.6(－5.9, 3.8) | －0.3(－5.1, 3.2) | 1.7(－3.9, 7.5) | 0.1(－3.2, 2.9) | |
| 验证站 | RMS | 5.5(1.8, 8.3) | 5.0(1.8, 8.3) | 4.5(1.4, 10.5) | 3.4(1.1, 6.6) |
| Bias | －0.3(－4.2, 4.8) | －0.1(－3.9, 3.6) | 1.5(－3.7, 8.3) | 0.0(－0.9, 1.0) |
表 3中，建模站的结果是对参与SHAtrop建模的152个测站进行统计，验证站的结果是对未参与建模的71个测站进行统计。可以看出，SHAtrop对于建模站和验证站的精度相当，RMS均为3.4 cm；建模站的ZTD精度比EGNOS、UNB3m、GPT2模型分别改善了39%、35%、28%；验证站的ZTD精度分别改善了38%、32%、24%。SHAtrop模型的Bias均值小于1 mm，优于EGNOS、UNB3m以及GPT2模型。在水汽变化活跃的长江三角洲区域内，由于利用了实际计算的ZTD进行建模，SHAtrop模型的拟合精度明显优于其他模型。同时在该区域模型的拟合精度低于其他区域，这是由于在水汽变化活跃的区域，ZTD存在模型没有考虑的高频信号。
5 结语本文利用上海天文台GNSS分析中心解算的陆态网223个测站2012年1月至2017年10月的对流层天顶延迟数据，分析了其时空特性，结果表明：中国区域对流层天顶延迟量随测站高程增加近似呈指数递减，随纬度增加近似呈线性递减；天顶延迟随时间的变化主要变现为年周期和半年周期，周期项的振幅与地域分布有关，年周期项与半年周期项的初相位均为10~35(年积日)，但在不同的地域分布不同。在以上分析的基础上，本文建立了适合中国大陆区域的对流层天顶延迟格网模型SHAtrop。SHAtrop按照陆态网络GNSS测站覆盖的范围（70°E~135°E, 18°N~54°N），将其划分为2.5°×2.0°的经纬格网，基于陆态网每个GNSS测站拟合得到的ZTD常数项、年周期项和半年周期项（A0，A1，d1，A2，d2）计算得到每个格网点的对应参数；将以上格网模型提供给用户使用，使用方便，与传统模型相比精度取得了良好的改善。相对于其他常用建模方法，本文采用了以GNSS实测数据作为数据源的建模方法，最大程度上保证了模型在GNSS应用中的自洽。此外，本文在高程改正方面也更为精细，因此能够保证模型的精度。
SHAtrop模型在中国区域内拟合的RMS为3.4 cm，相对于EGNOS、UNB3m、GPT2平均改善了39%、33%、26%；此外，SHAtrop模型平均系统偏差小于1 mm，优于EGNOS、UNB3m、GPT2等模型。该经验模型使用简单，改正效果良好，相对于传统模型有良好的精度改善，能够满足中国大陆区域GNSS导航定位用户的对流层延迟改正需求，可以作为中国大陆区域对流层延迟的参考模型。SHAtrop模型的格网数据文件以及使用程序在SHA分析中心官网(http://www.shao.ac.cn/shao_gnss_ac)可以获取。
| [1] |
Yin Haitao, Huang Dingfa, Xiong Yongliang, et al. New Model for Tropospheric Delay Estimation of GPS Signal[J]. Geomatics and Information Science of Wuhan University, 2007, 32(5): 454-457. (殷海涛, 黄丁发, 熊永良, 等. GPS信号对流层延迟改正新模型研究[J]. 武汉大学学报∙信息科学版, 2007, 32(5): 454-457. ) |
| [2] |
Niell A E. Global Mapping Functions for the Atmosphere Delay at Radio Wavelengths[J]. Journal of Geophysical Research Solid Earth, 1996, 101(B2): 3227-3246. DOI:10.1029/95JB03048 |
| [3] |
Boehm J, Niell A, Tregoning P, et al. Global Mapping Function (GMF): A New Empirical Mapping Function Based on Numerical Weather Model Data[J]. Geophysical Research Letters, 2006, 33(7). DOI:10.1029/2005gl025546 |
| [4] |
Kouba J. Implementation and Testing of the Gridded Vienna Mapping Function 1 (VMF1)[J]. Journal of Geodesy, 2008, 82: 193-205. DOI:10.1007/s00190-007-0170-0 |
| [5] |
Penna N, Dodson A, Chen W. Assessment of EGNOS Tropospheric Correction Model[J]. Journal of Navigation, 2000, 54(1): 37-55. |
| [6] |
Collins J P, Langley R B. The Residual Tropospheric Propagation Delay: How Bad Can It Get[C]. ION GPS, Nashville, Tennessee, USA, 1998
|
| [7] |
Leandro R F, Santos M C, Langley R B. UNB Neutral Atmosphere Models: Development and Performance[C]. ION NTM, Monterey, California, USA, 2006
|
| [8] |
Lagler K, Schindelegger M, BöHm J, et al. GPT2: Empirical Slant Delay Model for Radio Space Geodetic Techniques[J]. Geophysical Research Letters, 2013, 40(6): 1069-1073. DOI:10.1002/grl.50288 |
| [9] |
Wang Jungang, Chen Junping, Wang Jiexian. Analysis of Tropospheric Propagation Delay Mapping Function Models in GNSS[J]. Progress in Astronomy, 2014, 32(3): 383-394. (王君刚, 陈俊平, 王解先. GNSS对流层延迟映射模型分析[J]. 天文学进展, 2014, 32(3): 383-394. DOI:10.3969/j.issn.1000-8349.2014.03.07 ) |
| [10] |
Li Wei, Yuan Yunbin, Ou Jikun, et al. A New Global Zenith Tropospheric Delay Model IGGtrop for GNSS Applications[J]. Chin Sci Bull, 2012, 57(15): 1317-1325. (李薇, 袁运斌, 欧吉坤, 等. 全球对流层天顶延迟模型IGGtrop的建立于分析[J]. 科学通报, 2012, 57(15): 1317-1325. ) |
| [11] |
Song S L, Zhu W Y, Chen Q M, et al. Establishment of a New Tropospheric Delay Correction Model over China Area[J]. Science China(Physics, Mechanics and Astronomy), 2011, 54(12): 2271-2283. DOI:10.1007/s11433-011-4530-7 |
| [12] |
Chen Qinming, Song Shuli, Zhu Wenyao. The Establishment of the Global Model(SHAO-G) for the Tropospheric Delay[C]. China Satellite Navigation Foncerence, Shanghai, China, 2011 (陈钦明, 宋淑丽, 朱文耀.全球对流层延迟改正模型(SHAO-G)的初步建立[C].中国卫星导航年会, 上海, 2011)
|
| [13] |
Yao Yibin, He Changyong, Zhang Bao, et al. A New Global Zenith Tropospheric Delay Model GZTD[J]. Chinese Journal of Geophysics, 2013, 56(7): 2218-2227. (姚宜斌, 何畅勇, 张豹, 等. 一种新的全球对流层天顶延迟模型GZTD[J]. 地球物理学报, 2013, 56(7): 2218-2227. ) |
| [14] |
Mao Jian, Zhu Changqing, Guo Jifa. A New Global Zenith Tropospheric Delay Model[J]. Geomatics Information Science of Wuhan University, 2013, 38(6): 684-688. (毛健, 朱长青, 郭继发. 一种新的全球对流层天顶延迟模型[J]. 武汉大学学报∙信息科学 版, 2013, 38(6): 684-688. ) |
| [15] |
Chen Junping, Wu Bin, Hu Xiaogong, et al. SHA: The GNSS Analysis Center at SHAO[J]. Lecture Notes in Electrical Engineering, 2012, 160: 213. |
| [16] |
Chen Junping, Zhang Yize, Zhou Xuhua, et al. GNSS Clock Corrections Densification at SHAO: From 5 Minutes to 30 Seconds[J]. Science China(Physics, Mechanics and Astronomy), 2014, 57(1): 166-175. DOI:10.1007/s11433-013-5181-7 |
| [17] |
Wang Jungang, Chen Junping, Wang Jiexian, et al. Assement of Tropospheric Delay Correction Model over China[J]. Geomatics and Information Science of Wuhan University, 2016, 41(12): 1656-1663. (王君刚, 陈俊平, 王解先, 等. 对流层经验改正模型在中国区域的精度评估[J]. 武汉大学学报∙信息科学版, 2016, 41(12): 1656-1663. ) |