通常不低于10 km的基线叫中长基线，伴随着基线长度的增加，与距离有关的误差如电离层、对流层误差等的相关性会越来越弱^[1]。对于中长基线解算，电离层延迟和对流层延迟不会像短基线那样通过观测间双差处理就得以消除或减弱，各因素的影响更为复杂。例如，太阳活动会影响电离层延迟，对流层延迟误差与观测频率无关，因此, 难以建立正确的模型或组合观测值对其进行改正^[2-4]。为了得到比较理想的结果，除了尽可能选择更准确的初始坐标之外，选择合适的精密星历和钟差、电离层延迟、对流层延迟改正模型来进行解算是很有必要的。

本文通过对常用的GPS原始观测数据组合方式的分析，提出伪距宽巷组合的数据处理方法，有效消除了电离层、对流层传播误差，形成抗差性强、大气误差自由的GPS组合观测值，并通过对固定模糊度方法实现精确定位等问题进行了综述。

1 中长基线处理中误差的改正 1.1 对流层天顶延迟改正

鉴于对流层折射对GPS信号传播的影响情况很复杂，一般采用改正模型进行减弱，常用的有霍普菲尔德(Hopfield)公式、萨斯塔莫宁(Saastamoinen)公式和UNB3天顶延迟模型^[5]3种模型。UNB3天顶延迟模型的修正精度在x、y平面上与萨斯塔莫宁模型相当，但在高程方向上优于萨斯塔莫宁，所以这里只对UNB (University of New Brunswick)系列中的UNB3天顶延迟模型作详细说明。

UNB3模型跟GPS对流层延迟和其他改正模型一致，分为干延迟改正和湿延迟改正两部分。对于任意方向上的GPS对流层延迟改正可以表示为:

$ {d_{{\rm{trop}}}} = d_{{\rm{hyd}}}^z{m_{{\rm{hdy}}}} + d_{{\rm{wet}}}^z{m_{{\rm{wet}}}} $

(1)

式中，d_trop为总的对流层延迟量; d_hyd^z为天顶方向的对流层干延迟; d_wet^z为天顶方向的对流层湿延迟; m_dry为对流层干延迟的映射函数; m_wet为对流层湿延迟的映射函数。

需要注意的是，UNB3天顶延迟模型中的气象参数值是从海平面处开始起算的，包括了大气压P₀(mbar)、温度T₀(K)、水气压e₀¹(mbar)、温度变化率β(K/ m)和水气压变化率λ(mbar/m)，这些气象参数值依据的是测站处的大地纬度和测量时间，按照UNB3对流层天顶方向干、湿延迟的气象参数格网值进行内插获得。

UNB3对流层天顶方向的干、湿延迟模型为:

$ \left\{ \begin{array}{l} d_{{\rm{hyd}}}^z = ({10^{ - 6}}{K_1}{R_d}/{g_m}){(1 - \beta H/{T_0})^{g/({R_d}\beta )}}{P_0}\\ d_{{\rm{wet}}}^z = [{10^{ - 6}}{{K'}_3}{R_d}/({g_m}\lambda ' - \beta {R_d})] \times \\ \;\;\;\;\;\;\;\;\;{(1 - \beta H/{T_0})^{\lambda 'g/({R_d}\beta )}}^{ - 1}({e_0}/{T_0}) \end{array} \right. $

(2)

式中, K₁=77.604 K/mbar；K′₃=382 000 K²/mbar；g=9.806 65 m/s²；R_d=287.054 J/(10^-3kg)；H为用户的高程(单位为m)；g_m=9.784×[1-2.66×10^-3cos(2Φ)-2.8×10^-7H](m/s²)；λ′=λ+1。

1.2 映射函数模型

对于任意方向上的对流层延迟改正需要建立映射函数。常用的映射函数有NMF(Neill mapping function)函数、VMF1(Vienna mapping function)函数、GMF(global mapping function)函数。因为VMF1函数是目前精度最好、可靠性最强的模型，且该模型求得的基线向量比NMF模型具有更好的重复精度^[6-8]，所以本实验主要使用VMF1函数。下面详细介绍VMF1映射函数。

VMF1映射函数是由奥地利维也纳理工大学建立的模型, 分为干分量投影函数和湿分量投影函数，该大学的大地测量研究所依据实测气象资料生成的经差为2.5°、纬差为2°、时间间隔为6 h的格网提供系数a_d和a_w，欧洲中尺度天气预报中心提供系数b_d和c_d^[8]。其中b_d为常数，取0.002 9，c_d由下式计算：

$ \begin{array}{*{20}{l}} {{c_d} = {c_0} + \left[ {\left( {\cos \left( {\frac{{d - 28}}{{365}} \times 2\pi + \psi } \right) + 1} \right) \times \frac{{{c_{11}}}}{2} + {c_{10}}} \right] \times }\\ {\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left( {1 - \cos \varphi } \right)} \end{array} $

(3)

式中，d表示测量时的年积日；φ为测站纬度，其他各项系数如表 1所示。

1.3 电离层延迟改正(利用双频改正)

消除电离层观测量的观测方程为：

$ {\mathit{\Phi }_{{\rm{IF}}}} = \frac{\rho }{{{\lambda _1}}} - (\frac{{f_1^2}}{{f_1^2 - f_2^2}}{N_1} - \frac{{{f_1} - {f_2}}}{{f_1^2 - f_2^2}}{N_2}) + {\varepsilon _{{\rm{IF}}}} $

(4)

式中，ρ为伪距观测量；λ为载波波长；f为载波频率；N为模糊度；ε_IF为误差项。

电离层残差观测量的观测方程为：

$ {\mathit{\Phi }_1} = \frac{{\lambda _2^2-\lambda _1^2}}{{{\lambda _1}}} \times \frac{{40.3{\rm{TEC}}}}{{{c^2}}}-({N_1}-\frac{{{f_1}}}{{{f_2}}}{N_2}) $

(5)

式中，c为光速；TEC(total electron content)为信号传播路径上的电子总数。

利用双频线性组合而形成的消电离层观测量可将电离层的影响予以消除。消电离层线性组合观测量可以变形为：

$ {\mathit{\Phi }_c} = (\frac{1}{{1-{K^2}}})({\mathit{\Phi }_1}-K{\mathit{\Phi }_2}) $

(6)

式中，$K = \frac{{{f_2}}}{{{f_1}}}$。

1.4 卫星星历误差改正

从GPS卫星导航电文中获取的广播星历只适用于短基线。对于中长基线，其准确性明显不能满足要求^[9]。因此, 中长基线解算应采用精密星历。GPS卫星精密星历采用的是全球跟踪站的观测资料而计算出来的，目前由IGS精密星历提供最高的精度。利用IGS精密星历和钟差参数能达到高精度定位和定轨，这使得GPS在维护和更新ITRF(international terrestrial reference frame)参考框架、固体地球的变形监测，监测地极移动，确定卫星轨道和卫星监测电离层和对流层的科学任务中起着非常重要的作用^[10]。

2 消电离层组合解算模糊度及坐标 2.1 消电离层组合计算模糊度

这里的模糊度为实数，可将它分解为宽巷模糊度和L₁模糊度来进行解算。设Φ₁、Φ₂分别为L₁、L₂的相位观测值，则宽巷观测量为:

$ \begin{array}{l} {\mathit{\Phi }_w} = {\mathit{\Phi }_1}-{\mathit{\Phi }_2} = (\frac{1}{{{\lambda _1}}}-\frac{1}{{{\lambda _2}}})\rho + \\ ({f_1}-{f_2})\Delta \delta - ({N_1} - {N_2}) - (\frac{{\Delta {I_1}}}{{{\lambda _1}}} - \frac{{\Delta {I_2}}}{{{\lambda _2}}}){\rm{ = }}\\ \frac{{{f_w}}}{c}\rho + {f_w}\Delta \delta - {N_w} - 40.3\frac{{{\rm{TEC}}}}{c}\left( {\frac{1}{{{f_1}}} - \frac{1}{{{f_2}}}} \right) \end{array} $

(7)

式中，频率f_w=f₁-f₂；宽巷模糊度N_w=N₁-N₂; $\Delta {I_1} = \frac{{40.3{\rm{TEC}}}}{{f_1^2}}$; $\Delta {I_2} = \frac{{40.3{\rm{TEC}}}}{{f_2^2}}$。

由宽巷观测方程很容易求解出宽巷模糊度N_w：

$ {N_W} = \frac{{{f_w}}}{c} \times \rho + {f_w}\Delta \delta-{\mathit{\Phi }_w}-40.3\frac{{{\rm{TEC}}}}{c}(\frac{1}{{{f_1}}}-\frac{1}{{{f_2}}}) $

(8)

一般情况下，经常使用最小二乘降相关分解法(least squares ambiguity decorrelation adjustment, LAMBDA)方法来解算模糊度。LAMBDA方法是由序贯条件最小二乘模糊度搜索方法发展而来的。其基本思想是将模糊度进行降相关变换，然后在变换后的模糊度空间进行条件搜索，这样可有效缩小模糊度的搜索空间，进而提高搜索效率^[11]。当求得变换后的整数模糊度最优解后，可通过模糊度逆变换反求出原始模糊度空间中的模糊度。因此，LAMBDA方法通常可分为整周模糊度正变换、条件搜索以及整周模糊度逆变换等3个步骤。由此可得L₁的模糊度为：

$ \begin{array}{l} {N_1} =-{\mathit{\Phi }_1} + ({\mathit{\Phi }_w} + {N_W})\frac{{{f_1}}}{{{f_w}}} + \\ \;\;\;\;\;\;\frac{{40.3{\rm{TEC}}}}{{{f_1}}} \times \frac{{{f_1} + {f_2}}}{{{f_1}{f_2}}} \end{array} $

(9)

2.2 解算坐标

考虑有一个基准站和一个流动站的情况，使用两颗卫星的观测方程便可组成一个双差观测方程。设基准站为A测站，流动站为B测站，同步观测j、k两颗卫星，利用双差载波相位观测方程式得：

$ \Delta \nabla \mathit{\Phi } _{AB}^{jk} = \frac{f}{c}\left[{\rho _B^k-\rho _A^k-\rho _B^j + \rho _A^j} \right] -\Delta \nabla N_{AB}^{jk} + {\varepsilon _\mathit{\Phi } } $

(10)

式中，△▽为双差算子；ε_Φ为观测噪声；

$ \left\{ \begin{array}{l} \rho _A^k = {\left[{{{({X^k}-{X_A})}^2} + {{({Y^k}-{Y_A})}^2} + {{({Z^k}-{Z_A})}^2}} \right]^{1/2}}\\ \rho _{B_0}^k = {\left[{{{({X^k}-{X_{B_0}})}^2} + {{({Y^k}-{Y_{B_0}})}^2} + {{({Z^k}-{Z_{B_0}})}^2}} \right]^{1/2}}\\ \rho _A^j = {\left[{{{({X^j}-{X_A})}^2} + {{({Y^j}-{Y_A})}^2} + {{({Z^j}-{Z_A})}^2}} \right]^{1/2}}\\ \rho _{B_0}^j = {\left[{{{({X^j}-{X_{B_0}})}^2} + {{({Y^j}-{Y_{B_0}})}^2} + {{({Z^j}-{Z_{B_0}})}^2}} \right]^{1/2}} \end{array} \right. $

(11)

式(10)中的ε_Φ包括A、B两个测站双差残差，由于无法估计且较小，可忽略不计，因此，式(10)中的未知量就只剩下流动站B的3个位置参数X_B, Y_B, Z_B。通过前面的改正模型已经将电离层延迟误差与对流层延迟误差加以改正，其改正数用A表示，站距离相关的误差残差和多路径误差的残差无法计算，将其忽略。则可得双差载波相位观测方程的线性化形式为:

$ \begin{array}{l} \;\;\;\;\;\;\;\Delta \nabla \mathit{\Phi }_{AB}^{jk} = \frac{f}{c}[(l_B^j-l_B^k)(m_B^j-m_B^k)\\ \times (n_B^j-n_B^k)]\left[{\begin{array}{*{20}{c}} {\partial {X_B}}\\ {\partial {Y_B}}\\ {\partial {Z_B}} \end{array}} \right] + \frac{f}{c}\left[{\rho _{B_0}^k-\rho _A^k-\rho _{B_0}^j + \rho _A^j} \right] -\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\Delta \nabla N_{AB}^{jk} + A \end{array} $

(12)

其中，

$ \left\{ \begin{array}{l} l_B^j-l_B^k = (\frac{{{X_{B_0}}-{X^j}}}{{\rho _{B_0}^j}}-\frac{{{X_{B_0}} - {X^k}}}{{\rho _{B_0}^k}})\\ m_B^j - m_B^k = (\frac{{{Y_{B_0}} - {Y^j}}}{{\rho _{B_0}^j}} - \frac{{{Y_{B_0}} - {Y^k}}}{{\rho _{B_0}^k}})\\ n_B^j - n_B^k = (\frac{{{Z_{B_0}} - {Z^j}}}{{\rho _{B_0}^j}} - \frac{{{Z_{B_0}} - {Z^k}}}{{\rho _{B_0}^k}}) \end{array} \right. $

(13)

通过观测值以及B站的近似坐标，可以计算出B站的改正值∂X、∂Y、∂Z，从而便可以计算出B点的精确坐标值。

2.3 数据分析 2.3.1 实验情况说明

在2014年5月18日和5月29日两天，在某校校区利用3台合众思壮E660、6台南方S86T共9台接收机进行观测，共布设了20个控制点，其中GPS6、GPSB、GPSC、GPSD、GPSE、GPSF、GPSG、GPSH、GPSI、GPSM等10个点(称为基站点)的观测时间都在4 h以上，另外的10个点GPS1、GPS4、GPS5、GPS7、GPS8、GPS9、GPSA、GPSJ、GPSK、GPSL(称为扩展点)与基站点的同步观测时间都在1 h以上，截止高度角设置为10°，采样间隔为30 s。为了保证高精度的高程测量，仪器高都在不同的时间和地点分别量取了3次，天线也全部调整为北方向。

对于基站点，其观测时间较长，可通过精密单点定位技术(precise point positioning, PPP)或者中长基线解算技术获取其在观测时的ITRF2008框架下的坐标；而对于扩展点，可通过短基线解算软件(HD2003)解算得到。

2.3.2 精密单点定位解算与中长基线处理结果

本次实验数据处理采用网络精密单点定位软件GAPS处理，新伯伦瑞克大学颁布的GAPS数据处理软件的数据处理速度比较快，且数据报告极为详细^[12-16]。数据处理时，其仪器基本配置如表 2所示。

随着国内外定位技术的完善和IGS产品精度的大幅度提升，PPP技术极大程度地提高了GPS的单点定位精度。因此，通过比较精密单点定位解算结果与中长基线解算结果，可以分析中长基线解算经过误差改正后的数据精度。利用AUSPOS 2.1网络解算长基线软件进行解算，所用的IGS参考站有AIRA、BADG、BJNM、HKFN、HKSC、HKSL、HKWS、IRKJ、TCMS、YONS等10个。利用精密单点定位与AUSPOS 2.1网络解算长基线软件解算结果见表 3。

实验中，对于GPS6、GPSB、GPSC、GPSD、GPSE、GPSF、GPSG、GPSH、GPSI、GPSM这10个点，部分点间的基线长度D如下：D_6-B=931.099 0 m，D_6-c=1 023.855 5 m，D_6-D=1 230.311 2 m，D_6-E=1 483.860 2 m, D_6-F=1 315.540 6 m, D_6-G=1 206.461 1 m, D_6-H=1 128.422 4 m, D_6-I=1 025.320 2 m, D_6-M=835.064 2 m。由表 3可看出，精密单点定位解算结果与中长基线处理结果并无显著差异。

3 结束语

在中长基线解算过程中，电离层、对流层传播误差、卫星星历的选择、模糊度的解算都是影响解算结果的重要因素。本文利用UNB3天顶延迟模型和VMF1映射函数对对流层天顶延迟进行改正，利用双频改正模型对电离层延迟进行改正，利用双频观测消电离层组合固定并计算模糊度，从而解算待定点位坐标。通过对对流层、电离层的延迟改正和模糊度的解算，并将解算结果与精度较高的精密单点定位解算结果进行比较得出，两种解算结果差异较小，说明部分误差已被消除或减弱，解算结果的精度有了明显提高。