随着美国GPS的现代化、俄罗斯的GLONASS的完善、中国北斗卫星导航系统(BeiDou navigation satellite system, BDS)和欧盟伽利略（Galileo）导航系统的建设^[1-2]，全球导航卫星系统（global navigation satellite system, GNSS）以其全天候、高精度、实时连续性的优势成为室外大众导航的主要手段。传统的标准单点定位(single point positioning，SPP)算法以伪距作为主要观测值，由于观测噪声较大（分米甚至米级），且受多路径影响严重，定位精度仅10 m左右，定位结果不连续，用户导航体验差。

文献[3-6]针对如何利用高精度相位和多普勒观测值辅助伪距定位，克服SPP算法导航精度低、结果不平滑的缺陷进行了深入研究，并提出许多有效算法。其中Hatch滤波算法，即相位平滑伪距方法最为经典，该算法对于GNSS信号跟踪的连续性具有较高的要求，在船载、机载等观测条件良好的情况下，使用载波相位平滑伪距的方法可以获得精度1.5 m的定位结果。然而城市环境复杂，GNSS信号中断频繁，无法保证相位平滑伪距的连续性。文献[5]提出一种利用伪距和历元间载波相位变化量作为观测值，用多普勒观测估计载体的速度或加速度，并以此建立载体运动的状态方程，进而采用卡尔曼（Kalman）滤波估计载体状态的算法。但是该算法用常速度或常加速度状态模型描述载体的运动状态，在运动状态连续变化的情况下，状态方程的动态噪声不满足卡尔曼滤波要求动态噪声为高斯白噪声的假设条件，难以得到最优导航结果。另外，多普勒观测值作为一种瞬时观测值也常用于载体运动速度和加速度的估计^[6]。文献[7]利用多普勒观测值导出的载体速度和加速度来构建状态方程，描述载体的运动，可以得到精度更高、更平滑的载体运动轨迹，该算法有效避免了相位观测值周跳探测的问题，但是其状态方程的构建精度依赖多普勒观测值计算的载体速度和加速度的精度，在载体运动速度变化较大或者观测值采样间隔大时，状态方程无法准确描述载体的运动状态，该算法的效果将会下降。

基于以上分析，本文提出一种基于历元间载波相位差分的导航算法。该算法首先以高精度历元间相位变化量作为观测值，采用最小二乘方法得到精确载体位置变化量，并以此描述载体的运动状态。然后以伪距作为主要观测值，采用扩展Kalman滤波方法估计载体接收机的位置和钟差信息。最后在充分分析滤波特性的基础上，利用自主编制软件，对低成本单频u-blox接收机静态和动态条件下的实测GNSS数据进行实时处理，验证本文导航算法的效果。

1 历元间相位差分滤波算法 1.1 观测模型

对于单频GNSS接收机，伪距和相位观测方程为：

$\left\{ {\begin{array}{*{20}{l}} {\begin{array}{*{20}{l}} {P = \rho - c \cdot {\delta ^s} + c \cdot {\delta _r} + {d_{{\rm{eph}}}} + }\\ {\;\;\;\;\;\;\;{d_{{\rm{ion}}}} + {d_{{\rm{trop}}}} + {\varepsilon _P}} \end{array}}\\ {\begin{array}{*{20}{l}} {L = \rho - c \cdot {\delta ^s} + c \cdot {\delta _r} + {d_{{\rm{eph}}}} - }\\ {\;\;\;\;\;\;\;{d_{{\rm{ion}}}} + {d_{{\rm{trop}}}} - \lambda \cdot N + {\varepsilon _L}} \end{array}} \end{array}} \right.$

(1)

式中，P和L分别为伪距观测值和相位测值；ρ=$ \sqrt{\left(X^{s}-X_{r}\right)^{2}+\left(Y^{s}-Y_{r}\right)^{2}+\left(Z^{s}-Z_{r}\right)^{2}} $为卫星到接收机的几何距离，其中(X^s，Y^s，Z^s)为根据广播星历计算的卫星位置，(X_r，Y_r，Z_r)为接收机位置；c为真空中的光速；δ^s为根据广播星历计算的卫星钟差；δ_r为接收机钟差；d_eph为广播星历的轨道和卫星钟误差在信号传播方向上的投影；d_ion、d_trop分别为电离层和对流层延迟，分别使用克罗布歇(Klobuchar)模型和萨斯塔莫尼(Saastamoinen)模型进行改正；λ为载波波长；N为模糊度参数；ε_P为伪距观测多路径和测量噪声的影响；ε_L为相位观测多路径和测量噪声的影响。

在GNSS差分改正数可用的情况下，伪距观测方程可简化为：

$P = \rho - c \cdot {\delta ^s} + c \cdot {\delta _r} + {p_{{\rm{rc}}}} + {\varepsilon _P}$

(2)

式中，p_rc为基站播发的差分改正数。

已知广播星历轨道误差和钟误差在信号传播方向上的投影具有缓变特性^[8]，对于高频(大于1 Hz)观测值，历元间求差可以消去大部分的广播星历残余误差变化（小于1 mm/s），可认为$\Delta d_{\mathrm{eph}} \approx 0$，且同时$\Delta d_{\mathrm{ion}} \approx 0, \Delta d_{\mathrm{trop}} \approx 0 $。则在未发生周跳时，对历元间相位观测值求差，得到：

$\begin{array}{*{20}{c}} {{\rm{\Delta }}L = L\left( {{t_i}} \right) - L\left( {{t_{i - 1}}} \right) = {\rm{\Delta }}\rho - c \cdot {\rm{\Delta }}{\delta ^s} + }\\ {c \cdot {\rm{\Delta }}{\delta _r} + {\varepsilon _{{\rm{\Delta }}L}}} \end{array}$

(3)

式中，Δ表示差分算子。在使用高精度相位观测的同时，避免模糊度参数的解算也是历元间相位差分的一个特点。

假设t_i+1到t_i时刻接收机的位置和钟差变化量为$\Delta \boldsymbol{X}=\left(\Delta x, \Delta y, \Delta z, \Delta \delta_{r}\right) $，采用最小二乘方法求解历元间位置变化量得：

${\rm{\Delta }}\mathit{\boldsymbol{X}} = {({\mathit{\boldsymbol{H}}^{\rm{T}}}\mathit{\boldsymbol{PH}})^{ - 1}}{\mathit{\boldsymbol{H}}^{\rm{T}}}\mathit{\boldsymbol{Pl}} $

(4)

式中，H为系数矩阵；l为观测矩阵；P为观测值权阵。

1.2 滤波模型

利用历元间差分相位观测值求取高精度的位置变化量之后，可以构建基于历元间相位差分观测值的滤波导航算法的状态方程与观测方程：

$ \left\{ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{X}}_{{t_i}}} = {\mathit{\boldsymbol{X}}_{{t_{i - 1}}}} + {\rm{\Delta }}{\mathit{\boldsymbol{X}}_{{t_{i - 1}}, {t_i}}} + {\mathit{\boldsymbol{w}}_{{t_{i - 1}}, {t_i}}}}\\ {{\mathit{\boldsymbol{Z}}_{{t_i}}} = \mathit{\boldsymbol{B}}{\mathit{\boldsymbol{X}}_{{t_i}}} + {\mathit{\boldsymbol{v}}_{{t_i}}}} \end{array}} \right. $

(5)

式中，X为待估参数，在GPS观测中，X=$\left(\begin{array}{llll} {x} & {y} & {z} & {\delta_{r}^{\mathrm{GPS}}} \end{array}\right)$，在多系统观测中，只需增加对应的接收机钟差参数即可；w为系统噪声向量，主要表示历元间接收机位置和钟差变化量计算的不准确度，满足E(w) = 0，E (ww^T) = Q，其中E ()为求期望函数，Q为系统噪声协方差阵；ΔX直接采用历元间相位差分数据计算得到；Z为伪距观测值向量；B为系数矩阵；v为观测噪声向量，满足E(v) = 0，E (vv^T) = R，其中R为观测噪声协方差阵。

对于传统的常速度/常加速度模型，当历元间载体的运动状态发生突变时，运动方程的构建与真实的运动不符，预报位置误差较大，导航效果有限。而本文基于历元间相位差分（phase difference，PD）算法的最大特点在于状态方程的构建，采用t₁时刻位置加上历元间相位差分观测值计算的位置变化量预报t₂时刻的位置，与载体真实运动情况一致。同时，基于历元间位置变化量的滤波算法只需要前后历元有不少于4颗卫星的连续观测即可，即使部分卫星跟踪发生中断，算法性能不受影响。因此理论上该算法更适用于城市内受建筑物、树木遮挡影响的环境下，相位观测经常不连续的导航。

2 定位性能分析实验 2.1 静态实验

静态实验采用u-blox公司生产的NEO-M8T型号单频多系统GNSS接收机，于2017-11-04至2017-11-06，静态条件下采集了BDS/GPS单频观测数据，采样间隔1 s，该观测点的准确坐标由高精度的实时动态（real-time kinematic, RTK）测量给出。同时，在距离该接收机2 m处架设另一基站，用于计算并通过网络播发GNSS差分改正数。

图 1为使用PD算法解算的连续3 d的结果在东/北/天（east/north/up，E/N/U）方向上的误差曲线图。如果定位结果中包含历元间差分后残余的星历误差、电离层误差、对流层误差，则理论上随着时间的推移，E/N/U方向误差会逐渐增大。然而由图 1可知，连续3 d的定位误差并没有明显的增大趋势，并且各天的误差在E/N/U方向上有明显的重复性，这应该是受到广播星历误差的影响。图 2（a）为2017-11-06 00:00—02:00，3种算法的静态定位结果在E/N/U 3个方向的误差序列图，图 2（b）为该时间段使用GNSS差分改正数改正了广播星历轨道和钟差误差后的PD与SPP算法定位结果在E/N/U 3个方向的误差序列。可以发现，相对于SPP结果，u-blox接收机内嵌算法（UBX）的结果相对平滑，定位精度在E/N/U 3个方向上均有提高。但结果仍然明显受伪距观测噪声的影响，误差曲线波动较大。而使用本文PD滤波算法，定位结果进一步平滑，且E/N/U 3个方向上的定位精度进一步提高。PD算法以低噪声的历元间的相位变化量作为主要观测值，即使对于多系统观测，仅需4颗连续可见的卫星即可完成滤波时间更新。同时历元间相位差分大幅降低了卫星位置和钟差等误差对定位的影响，其定位结果在稳定性与精度方面均优于SPP和UBX结果。

对于行人等大众导航，普遍更关心水平方向结果的精度和可靠性。2017-11-06 00:00—02:00, PD、UBX和SPP算法结果在平面方向上的误差分布如图 3（a）所示。可以发现，在平面上，PD算法结果的稳定性最优，UBX结果次之，SPP最差。同时3种算法的结果误差存在一个E方向上约0.2 m、N方向上约0.5 m的系统性偏差。这说明使用历元间相位差分观测值虽然可以计算高精度的历元间相对位置，但定位的绝对位置仍依赖伪距观测值给出，其必然会受到广播星历轨道和钟差误差的影响。使用GNSS差分改正数之后，3种算法定位的平面误差如图 3（b）所示，可以发现广播星历引起的定位系统误差得到克服，PD算法平面精度达到0.16 m。3种算法静态定位结果的均方根误差(root mean square，RMS)统计如表 1所示。由表 1可知，使用GNSS差分改正数后，相较于UBX和SPP结果，PD算法结果在E/N/U 3个方向上的精度均有明显提高，系统误差得到有效削弱。在未使用GNSS差分改正数时，PD算法平面RMS和三维(3D) RMS分别为0.56 m、1.60 m；在使用GNSS差分改正数之后，分别可达到0.16 m和0.75 m。

2.2 动态实验

本文在开阔和遮挡环境下设计了两组动态实验，进行PD算法效果验证。实验过程中，u-blox接收机的天线由行人手持，距离地面约1.4 m。

2.2.1 开阔环境

开阔环境的实验在武汉大学信息学部足球场上进行，采用u-blox公司生产的NEO-M8T型号单频多系统GNSS接收机，于2017-11-04至2017-11-06绕直径18 m的大圆连续行走7圈，采样间隔1 s。观测场地通视条件良好，导航结果如图 4所示，图中S和T分别代表导航轨迹的起点和终点，参考轨迹由RTK解算得出。3种结果与参考坐标求差之后的误差曲线如图 5所示。3种导航算法结果的精度统计见表 2。

从图 4、图 5和表 2的结果可以发现：（1)从图 4中可以发现，3种算法结果与参考轨迹的整体运动吻合，其中PD算法的结果最平滑。在使用GNSS差分改正数数据之前，PD算法在E/N/U方向上的标准差(standard deviation，STD)分别为0.07 m、0.08 m、0.43 m，而使用GNSS差分改正数之后，对应E/N/U方向上的STD分别为0.08 m、0.10 m、0.15 m；对于u-blox接收机内嵌某种与速度相关的滤波算法，依赖瞬时多普勒观测值计算的速度，在直线部分的运动轨迹与参考估计吻合，在绕大圈运动时候，运动方向不断发生变化，瞬时速度并不能反应历元间真实的平均运动速度，因此该算法的效果略有下降，E/N/U方向对应的STD依次为0.32 m、0.31 m、0.58 m。由于伪距观测噪声的影响，SPP结果使用差分改正数前后的改进不大，整体结果仍存在明显跳动，STD约为0.55 m、0.70 m、1.10 m。（2)从图 5、表 2中可以发现，在未使用GNSS差分改正数时，PD结果误差在N/U方向上存在的系统偏差依次约为0.4 m、0.8 m；使用GNSS差分改正数改正广播星历轨道和钟差的系统误差后，在N/U方向上的误差曲线与0更为接近。而对于SPP结果，同样使用GNSS差分改正数。虽然校正了系统误差，但是由于伪距观测噪声的影响，其E/N/U方向的误差序列呈现随机波动的特性，用户导航体验会较差。（3)从表 2中可以发现，在使用PD算法有效消除伪距观噪声的影响，且使用GNSS差分改正数消除广播星历轨道和钟差的系统误差后，可以获得高精度的定位结果。使用GNSS差分改正数前，平面精度为0.41 m，垂向精度0.86 m；而使用GNSS差分改正数后，平面和垂向精度进一步提高达到0.21 m和0.46 m，提高比例约49%和46%。（4）在开阔环境下，PD算法可以为用户提供连续、平滑且高精度的导航体服务。

2.2.2 遮挡环境

为了分析基于历元间差分相位观测值的滤波算法在信号遮挡条件下的导航效果，本实验使用u-blox接收机在林荫道上进行测量，测试环境包括完全开阔、部分遮挡等各种场景，测试路线及沿路环境如图 6所示。

整个测试轨迹可以分为两个部分，第一部分在视野开阔的操场上，第二部分在树木茂密的林荫道上，该部分GNSS观测中断严重，即使使用事后RTK解算仍不能给出可靠的参考坐标。在实验路径上设置23个参考点（图 6中红点）作为导航性能评估的参考。这些点坐标由双频RTK解算或者全站仪传递得到。

在整个测试过程中，GNSS卫星的可用情况和位置精度因子（position dilution of precision，PDOP）如图 7所示。实验前半部分，由于观测条件良好，通常有至少3颗BDS卫星和7颗GPS卫星可用，PDOP值稳定在2.5左右。而在实验后半段，由于树木对信号的遮挡，可用GNSS卫星数在0~13波动。常规的导航算法（如Hatch滤波）要求GNSS信号连续跟踪，在此情况中并不适用。

使用PD、UBX和SPP 3种方法导航的轨迹如图 8所示。可以看出，PD和UBX的轨迹是平滑的，而SPP结果波动较大。图 8(b)中的放大图为行人经过一面约3 m高墙(白色点)时，PD和UBX的导航结果，该处对应的可见卫星数(仅4~5颗)和PDOP值如图 7中绿色框所示。PD算法导航轨迹与参考点吻合，而UBX结果虽然同样平滑，但是存在一个明显的系统偏差，说明在前后历元仅4颗卫星连续观测时，本文PD算法仍能提供连续、可靠的导航结果。

对算法结果使用以下准则评估精度：在参考点附近，计算各算法估计的位置与参考位置之间的三维距离，将最小距离作为导航精度的度量。由于伪距观测噪声的影响，SPP的结果出现随机波动，因此该评价方法并不适用于SPP结果。PD和UBX结果在E/N/U方向上的导航误差如图 9所示，可以明显看出，PD算法结果误差更小，更适合该遮挡环境。最终PD算法的平面统计精度约0.75 m，远优于UBX的1.33 m。

3 结语

本文提出了一种基于历元间相位差分观测值约束的导航算法，该算法合理利用相位观测值噪声小的特点，计算高精度历元间载体位置变化量，构建状态方程。

利用低成本单频u-blox接收机进行静态和动态实验，本文提出的PD算法与SPP算法、UBX算法的对比结果表明：（1）在静态条件下，无论有无GNSS差分改正数信息，PD算法的结果均为最优。未使用GNSS差分改正数时，PD算法单机导航的平面精度达到0.56 m，垂向精度能够达到1.50 m；使用GNSS差分改正数后，平面和垂向精度分别提高到0.16 m和0.73 m。（2）动态条件下的两组实验结果类似，使用PD算法的导航结果最为平滑，平面精度能达到1 m左右，且即使在遮挡环境下，前后历元仅有4颗卫星连续可用时，该算法仍能输出平滑、可靠的导航结果，为用户提供优良的导航体验。这说明本文提出的基于历元间相位差分观测值约束的导航算法更适用于城市环境下的GNSS导航。