2015, Vol. 44
[PDF全文]
2. 解放军94259部队 山东 蓬莱 265660;
3. 空军工程大学信息与导航学院 西安 710077
2. Unit 94259 Penglai Shandong 265660;
3. Information & Navigation College, Air Force Engineering University Xi'an 710077
对流层位于大气底层，一般将距地面10～15 km的大气层称为对流层。对流层靠近地面，直接决定了地面天气气象情况，同时也对穿过该层的电磁波产生折射、反射和延迟等效应，影响卫星通信和散射通信等通信效果。针对对流层延迟效应，国内外学者做了大量研究，核心思想是将对流层天顶延迟(zenith troposphere delays，ZTD)通过映射函数的映射得到斜延迟(SPD)。近年来，随着参数气象测试记录和相关探空技术的发展，映射函数考虑的因素越来越充分，精度越来越高，如文献[1, 2]根据26个全球探空气球站数据和200 hPa气压层高度的相关数据提出了Niell映射函数(Niell mapping function，NMF)和等压映射函数(isobaric mapping functions，IMF)[2]。文献[3]根据欧洲中期天气预报中心(European centre for medium-range weather forecasts, ECMWF)40年气象数据资料于2006年提出了维也纳映射函数1(Vienna mapping functions 1，VMF1)和全球映射函(global mapping function，GMF)[4]等。
本文将设计基于对流层散射信道双向时间比对系统和对流层链路延迟进行建模分析，并结合映射函数和东北亚中纬度地区部分城市气象探空数据对对流层双向时间比对中的斜路径延迟进行研究。
1 对流层双向时间比对斜延迟分析 1.1 对流层双向时间比对原理借鉴卫星双向时间比对技术[5, 6]和流星余迹双向比对技术[7]理论，结合对流层散射通信原理，对流层双向时间比对(two way troposphere time transfer，TWT3)系统如图 1所示。
 |
| 图1 TWT3组成框图 |
A、B两站配置相同，由本地原子钟源(一般采用铷钟或者铯钟)、时间间隔计数器(time interval counter，TIC)、相关数据处理单元等部分组成。比对时刻A、B两站同时向对方站发送秒脉冲(one pulse per second，1pps)信号，该信号在发送前经过调制至中频后混频到散射频段(4.4～5.0 GHz)上，经发射机发射后经过对流层散射信道的传输被对方站接收，接收端经过与发射端相逆的信号处理过程对接收到的信号进行处理，该过程时间关系如图 2所示。
 |
| 图2 TWT3原理几何示意图 |
根据图 2的时间几何关系，A、B两站间的原子钟差△T可表示为：
| ${\rm{\Delta }}T = \frac{1}{2}\left\{ \begin{array}{l} ({\rm{TI}}{{\rm{C}}_{{\rm{AB}}}} - {\rm{TI}}{{\rm{C}}_{{\rm{BA}}}}) + ({t_{{\rm{TA}}}} - {t_{{\rm{RA}}}}) - ({t_{{\rm{TB}}}} - {t_{{\rm{RB}}}})\\ + ({t_{{\rm{ATropo}}}} - {t_{{\rm{TropoA}}}}) - ({t_{{\rm{BTropo}}}} - {t_{{\rm{TropoB}}}}) \end{array} \right\}$ | (1) |
式中，TICAB 和TICBA 分别表示A、B两站同时向对方站发送1pps信号时间间隔计数器开始计数至收到对方站1pps信号后停止计数时的计数间隔；tTA和 tRA、tTB和tRB 分别表示A、B两站的发射时延和接收时延；tTAropo、tBTropo和tTropoA、tTropoB分别表示A、B两站发射信号上行到对流层散射体间的时延和信号经对流层散射体后下行到A、B两站间的时延[8]。
1.2 对流层延迟映射函数分析关于对流层路径延迟误差的分析，文献[9]提出了一种连分式映射函数方法，映射函数为：
| ${\rm{M}}{{\rm{F}}_i} = \frac{{1 + \frac{{{a_i}}}{{1 + \frac{{{b_i}}}{{1 + {c_i}}}}}}}{{\sin \theta + \frac{{{a_i}}}{{\sin \theta + \frac{{{b_i}}}{{\sin \theta + {c_i}}}}}}}{\rm{ }}\begin{array}{*{20}{c}} \begin{array}{l} i = h{\rm{ 静力学映射}}\\ i = w\;湿映射 \end{array} \end{array}$ | (2) |
式中，MFi表示静力学映射函数或湿映射函数；θ为高度角；参数ai，bi，ci表示i条件下的参数值，其值因映射函数不同而不同。
1) NMF模型
Niell根据26个全球探空气球站数据，建立了全球大气延迟映射函数，并根据探空数据总结得出了部分纬度下参数a，b，c的经验值[1]，其他纬度条件下的参数值根据线性内插获得，参数a如式(3)和式(4)所示，求解参数b，c方法与求解参数a的方法类似[10- 11]：
| $a(\varphi ,{\rm{doy}}) = \left\{ \begin{array}{l} {a_{{\rm{avg}}}}(\varphi ) - {a_{{\rm{amp}}}}(\varphi )\cos \left( {2{\rm{\pi }}\frac{{{\rm{doy}} - {D_0}}}{{365.25}}} \right)\begin{array}{*{20}{c}} {}&{\left| \varphi \right|} \end{array} \le 15^\circ \vee \left| \varphi \right| \ge 75^\circ \\ {a_{{\rm{avg}}}}(\varphi ) + [{a_{{\rm{avg}}}}({\varphi _{i + 1}}) - {a_{{\rm{avg}}}}({\varphi _i})]\frac{{\varphi - {\varphi _i}}}{{{\varphi _{i + 1}} - {\varphi _i}}} - \\ \left( {{a_{{\rm{amp}}}}(\varphi ) + [{a_{{\rm{amp}}}}({\varphi _{i + 1}}) - {a_{{\rm{amp}}}}({\varphi _i})]\frac{{\varphi - {\varphi _i}}}{{{\varphi _{i + 1}} - {\varphi _i}}}} \right)\cos \left( {2{\rm{\pi }}\frac{{{\rm{doy}} - {D_0}}}{{365.25}}} \right)\begin{array}{*{20}{c}} {}&{15^\circ < \left| \varphi \right|} \end{array} < 75^\circ \end{array} \right.$ | (3) |
| $\begin{array}{l} {c_{\rm{h}}} = {c_0} + \left[ {\left( {\cos (\frac{{{\rm{doy}} - 28}}{{365}}2{\rm{\pi }} + \psi ) + 1} \right)\frac{{{c_{11}}}}{2} + {c_{10}}} \right] \times (1 - \cos \varphi )\\ \end{array}$ | (4) |
2) VMF1模型
文献[3, 12, 13]在Niell的工作基础上，采用欧洲中期天气预报中心40年的观测数据，修改完善了连分式参数，形成了VMF1函数。其参数ah，aw可从http://ggosatm.hg.tuwien.ac.at网站下载，并利用全球栅格网点实测数据内插得到测站对应的ah，aw值，bh=0.002 905，ch值如式(4)所得，对应参数如表 1所示。
| 表1 VMF1静力学系数相关参数表 |
 |
湿映射函数继承了NMF函数的参数，直接取45°下的NMF湿映射函数参数值：bw=0.001 46，cw=0.043 91。
3) GMF模型
VMF1精度和可靠性较高，但参数a的获得有34 h时延，为解决此问题，文献[4, 14]在ECMWF数值天气模型的基础上提出了GMF模型。
GMF模型在VMF1的基础上发展而来，其映射函数模型采用VMF1模型的系数b，c系数a则修改为：
| $a = {a_0} + A\cos \left( {\frac{{{\rm{doy}} - 28}}{{365}}2{\rm{\pi }}} \right)$ | (5) |
| ${a_0} = \sum\limits_{n = 0}^9 {\sum\limits_{m = 0}^n {{P_{nm}}(\sin \varphi )} } [{A_{nm}}\cos (m\lambda ) + {B_{nm}}\sin (m\lambda )]$ | (6) |
式中，a0为各站记录年数据的均值；年平均幅值A根据相位差估计得到；Pnm为n阶m次Legendre多项式；φ、λ分别为测站纬度和经度；Anm、Bnm为最小二乘计算得到的球谐系数。
1.3 对流层双向时间比对斜延迟分析斜延迟与对流层天顶延迟间的几何对应关系如图 3所示，为了便于分析和研究，将对流层天顶延迟分割为静力学天顶延迟(zenith hydrostatic delay，ZHD)和湿天顶延迟(zenith wet delay，ZWD)两部分，将斜延迟也对应分为静力学斜延迟(slant hydrostatic delay，SHD)和湿斜延迟(slant wet delay，SWD)两部分。
 |
| 图3 斜延迟与天顶延迟几何关系示意图 |
| $\left\{ \begin{array}{l} {\rm{SPD = SHD + SWD}}\\ {\rm{SHD = ZHD}} \times {\rm{M}}{{\rm{F}}_{\rm{h}}}\\ {\rm{SWD = ZWD}} \times {\rm{M}}{{\rm{F}}_{\rm{w}}} \end{array} \right.$ | (7) |
目前对于对流层散射延迟机理学术界没有统一结论，针对其中对流层斜延迟，本文以对流层散射模型为机理，结合映射函数分层理论，对对流层双向时间比对中对流层斜延迟进行建模分析，做如下假设：
1) 根据映射函数的分层推理理论思想，将对流层大气简化成分层大气，各层内折射率和气象条件相同；
2) 对流层高度与具体气象气候条件密切相关，本文假设为15 km；
3) 在对流层散射通信中，尽管电磁信号在对流层大气中传输路径各异，也不会全经过整个对流层，为了能够得到延迟最大值，假设上下行信号路径均包括整个对流层大气。
基于上述假设，根据对流层双向时间比对原理，结合图 3所示的斜延迟与天顶延迟几何关系，简化建立对流层延迟反射机理图如图 4所示(根据映射函数的简化原理，电波实际路径用直线简化标示)，设A、B两站的高度角为θ。
 |
| 图4 对流层延迟反射机理示意图 |
根据对流层双向时间比对原理，A站高精度1pps信号通过相关信号处理方法从发射机发射后，经过对流层延迟到达对流层层顶O，经过大量发射和折射后到达B站，根据对流层散射机理，只有与B站接收天线方向一致的反射信号才能到达B站，可用如下理论解释：散射体内的大量有极分子相当于一个个小偶极子天线，它们受到发射电波的激励后在偶极子天线方向上产生二次辐射，只有少量与发射波的波前相切的偶极子天线才能受到最大的激励，它们的二次辐射波强度最大，形成前向散射波。另外也只有这少量偶极子天线的轴向与接收天线的线极化方向一致，才能被接收天线有效接收。同理，B站信号到A站的过程与此类似[15]。
根据Snell定理和电波传输可逆性原理，在简化分层化的对流层大气中，图 4中对流层A→B传输过程可以分解为A→O→B过程，结合图 3所示的对流层斜延迟和天顶延迟的几何关系，将A→O和O→B过程分别利用A、B两站的气象地理条件进行延迟计算，得到：
| $\left\{ \begin{array}{l} {d_{{\rm{AB}}}} = {d_{{\rm{AO}}}} + {d_{{\rm{OB}}}}\\ {d_{{\rm{AO}}}} = {\rm{SH}}{{\rm{D}}_A} + {\rm{SW}}{{\rm{D}}_A}\\ {d_{{\rm{OB}}}} = {\rm{SH}}{{\rm{D}}_B} + {\rm{SW}}{{\rm{D}}_B} \end{array} \right.$ | (8) |
本文以卫星信号经过对流层所产生的延迟误差相关成熟理论为基础，将对流层大气做简化分层后，建立了对流层双向时间比对斜延迟计算模型，主要考虑的是两种情况所经过的介质均为对流层大气，延迟误差均为对流层散射体对信号的折射、散射等效应产生的，在对流层相关假设的前提下，对流层双向时间比对斜延迟模型中dOB可近似认为是信号经过对流层最大高度的对流层大气后产生的对流层传播延迟。
下文将结合东北亚中纬度地区部分测站气象和探空数据对对流层斜延迟进行分析研究。
2 算例运算分析考虑到测站情况和对流层散射通信的实际距离，并考虑到我国国内GPS测站分布较为分散的情况，本文选取与我国气象特征类似的东北亚中纬度日本地区4个GPS测站进行运算分析，所选测站基本地理情况如表 2所示，根据测站地理位置信息和1.3节中对流层大气相关假设，计算比对站之间的距离L和高度角θ如表 3所示。
| 表2 测站情况表 |
|
| 表3 所选站间的参数表 |
|
本文采用GMF函数分析表 3中前3组测站间的总斜延迟(即式(8)中的dAB)分布曲线，并根据分布情况对其进行拟合，如图 5a～图 5c所示；选择最后一组测站分析不同的映射函数下的总斜延迟分布情况，如图 5d～图 5f所示，散点图为实际计算得到的延迟分布情况，曲线为散点图的二阶傅里叶拟合曲线图。
 |
| 图5 不同映射函数下各站间SPD分布曲线图 |
从图 5a～图 5c可知，在GMF和上述条件设置的前提下，随着通信传输距离的增加，高度角相应减小，使得所选测站对的SPD增大；SPD呈年周期余弦函数分布，且大致呈夏季最大冬季最小的分布特征，最大最小值差值约2～3 m。从图 5d～图 5f可知，所选的KSMV和KGNI站之间在上文所述的参数设置条件下，不同映射函数下的SPD大致相同，而GMF根据数值天气模型计算得到的SPD解决了VMF1在参数获得方面34 h时延的缺陷，且其精度与VMF1相当，这也是GMF的优势所在。
为了计算对流层双向时间比对的斜延迟，需要分析链路延迟误差，根据对流层双向时间比对原理，在比对时刻近似认为对流层散射体稳定，则两站间信号传输路径上的误差可能来自于二者比对时刻高度角的微小变化，当△θ时，GMF下，选择表 3中前3组测站间双向时间比对SPD误差分布情况如图 6所示，并选择A&C测站对的误差进行月统计和年统计，其结果分别如表 4和表 5所示。
 |
| 图6 GMF下SPD误差分布曲线图 |
| 表4 A&C测站间GMF下SPD月误差统计表(2012年) |
|
| 表5 A&C测站间GMF下SPD年误差统计表 |
|
从图 6可知，由于高度角的微小变化导致SPD变化分布近似呈余弦函数分布，呈现夏季最小冬季最大的分布特征，且波峰持续时间较长，波谷持续时间较短。从图 6和表 4、表 5分析可知，单向误差最大值3 cm左右，带来的双向时延误差最大约为0.2 ns。
由于此分析在假设比对时刻对流层散射体近似稳定的前提下进行的，而实际对流层散射体是随机变化的，故实际测量值估计将大于0.2 ns。
3 结 束 语本文设计了对流层双向时间比对方法，对比分析了NMF、VMF1和GMF共3种映射函数特征和参数选择方法，对对流层大气进行了简化假设并建模分析，提出了一种基于映射函数的对流层双向时间比对斜延迟分析方法，通过选择与我国气候条件类似的东北亚中纬度日本地区4个GPS观测站，根据实测探空数据进行了对比分析，所得结论对研究我国东部沿海地区对流层双向时间比对延迟分析具有借鉴意义。
研究表明，所选测站在本文所示的条件下，SPD呈年周期余弦函数分布，且夏季最小冬季最大；随着测站间通信距离的增大和高度角的减小，SPD相应增大；所选测站下3种映射函数SPD分布大致相同，GMF较VMF1有实时优势；当△θ=-0.01，对流层双向时间比对链路时延误差理论估计约为0.2 ns。
| [1] | NIELL A E. Global mapping functions for the atmospheric delay at radio wavelengths[J]. J Geophys Res, 1996, 101(B2): 3227-3246. |
| [2] | NIELL A E. Improved atmospheric mapping functions for VLBI and GPS[J]. Earth Planets Space, 2000(52): 699-702. |
| [3] | BOEHM J, WERL B, SCHUH H. Troposphere mapping functions for GPS and very long baseline interferometry from European centre for medium-range weather forecasts operational analysis data[J]. J Geophys Res, 2006, 111(B2): B02406. |
| [4] | 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): L07304. |
| [5] | SONG H, IMAE M. Improvement of precision on two-way satellite time and frequency transfer[J]. Electronics Letters, 2009, 45(18): 945-946. |
| [6] | LIN H T, LIAO C S, CHU F D, et al. Full utilization of TWSTT network data for the short-term stability and uncertainty improvement[J]. IEEE Transactions on Instrumentation and Measurement, 2011, 60(7): 2564-2569. |
| [7] | KORNEYEV V A, EPICTETOV L A, SIDOROV V V. Time & frequency coordination using unsteady, variable-precision measurements on meteor burst synchronization and communication equipment[C]// Proceedings of the 2003 IEEE International Frequency Control Symposium and PDA Exhibition Jointly with the 17th European Frequency and Time Forum. [S.l]: IEEE, 2003. |
| [8] | 刘强, 孙际哲, 陈西宏, 等. 对流层双向时间比对及其时延误差分析[J]. 测绘学报, 2014, 43(4): 341-347.LIU Qiang, SUN Ji-zhe, CHEN Xi-hong, et al. Analysis of two way troposphere time transfer and its delay errors[J]. Acta Geodaetica Cartographica Sinica, 2014, 43(4): 341-347. |
| [9] | MARINI J W. Correction of satellite tracking data for an arbitrary tropospheric profile[J]. Rad Sci, 1972, 7(2): 223-231. |
| [10] | MENDES V B, LANGLEY R B. An analysis of high-accuracy tropospheric delay mapping functions[J]. Phys Chem Earth (A), 2000, 25(12): 809-812. |
| [11] | SAHA K, RAJU C S, PARAMESWARAN K. A new hydrostatic mapping function for tropospheric delay estimation[J]. Journal of Atmospheric and Solar-Terrestrial Physics, 2010(72): 125-134. |
| [12] | KOUBA J. Implementation and testing of the gridded Vienna mapping function 1 (VMF1)[J]. J Geod, 2008(82): 193-205. |
| [13] | BOEHM J, KOUBA J, SCHUH H. Forecast Vienna mapping functions 1 for real-time analysis of space geodetic observations[J]. J Geod, 2009(83): 397-401. |
| [14] | KOUBA J. Testing of global pressure temperature (GPT) model and global mapping function (GMF) in GPS analyses [J]. J Geod, 2009(83): 199-208. |
| [15] | 张明高. 对流层散射传播[M]. 北京: 电子工业出版社, 2004: 8-19.ZHANG Min-gao. Troposcatter propogation[M]. Beijing: Electronic Industry Press, 2004: 8-19. |