On-the-fly Calculation of the Whole Cycle Ambiguity --Data Processing Method on GNSS Kinematic Carrier Phase Measurements (5)
--Data Processing Method on GNSS Kinematic Carrier Phase Measurements (5)
Liu Jiyu, Chen Xiaoming
(School of Geodesy and Geomatics, Wuhan University, Wuhan, 430079)
Abstract: The correct calculation of the whole cycle ambiguity is the key problem of high-accuracy kinematic positioning using GNSS carrier phase measurements. Taking the GPS carrier phase measurement kinematic positioning as an example, this paper discusses On-the-fly calculation of the whole cycle ambiguity.
Keywords: GNSS kinematic carrier phase measurements; Whole cycle ambiguity; On-the-fly calculation
整周模糊度的正确求解,是用GNSS载波相位
测量进行高精度动态定位的关键问题。一旦正确解算出整周模糊度,载波相位观测值即可转换为高精度的站星距离测量值,进而实现厘米级精度的动态
定位。本文以GPS载波相位测量动态定位为例,论
述整周模糊度的在航解算。
在连续跟踪而无周跳的情况下,一次卫星通过的载波相位观测值,均含有相同的初始整周模糊度
NJ。为了正确解算出这个整周模糊度,在传统的GPS
载波相位测量动态定位中,往往采用在动态定位之前进行一段时间的静态测量或者在已知基线上进行
短时间的静态测量,称之为静态初始化测量。1985年,美国学者Remondi提出了一种天线交换的方法 来解求初始整周模糊度;这种方法既无需在已知基线上开始动态测量,又无需进行较长时间的静态
测量,只需在基准站附近(3~5m)选择一个临时测站,分别架设GPS信号接收机,而采集几个时元的
静态测量数据;然后在基准站和临时测站之间进行天线交换,再采集几个时元的静态测量数据;最后
再将GPS信号接收天线交换回到原测站上去。这种
天线交换法,实施简单,且软件设计又不复杂,而广泛用于传统的高精度动态定位。
无论是静态初始化测量,还是天线交换法,都是利用动态定位实施之前的测量数据来确定整周模糊度;即在动态定位之前,正确解算出整周模糊度,并依据无周跳时整周模糊度的不变特性,方能
用后续的GPS载波相位观测值,进行高精度的动态
定位解算。上述天线交换等方法,显然难以适用海上的动态定位以及飞机的精密进场着陆等许多动态应用场合。此外,即使可以在动态定位之前采用上述方法确定整周模糊度,但是,因动态定位环境的
复杂性,失锁和周跳实在难以避免。一旦在GPS动
态定位过程中出现失锁和周跳,以致某一时间段内
连续跟踪的GPS卫星数少于4颗,即使后来又重新锁定4颗以上的GPS卫星,重新锁定后的卫星整周模糊度仍无法解求,而导致动态定位精度的降低。
为了实现连续而可靠的高精度GPS动态定位,
近年来国内外学者开始着手研究在动态定位过程中正确解算整周模糊度的方法,称之为整周模糊度的
在航解算(Ambiguity Resolution On the Fly),简称OTF法或AROF法,并取得了丰硕的研究成果,导出了多种OTF解算方法。许多文献也介绍了一些OTF方法;但这些OTF方法可综合划分为如下四类:双频P码伪距法;模糊度函数法;最小二乘搜索
法;模糊度协方差法。
1 双频P码伪距法
最简单的整周模糊度求解方法是直接利用伪距观测量来确定载波相位的整周模糊度。依前述, L1-P/L2-P的双频伪距观测值构成的宽巷和窄巷伪距观测值分别为
宽巷伪距: Pd = ( P 1 − P2 )λd = ρ * − I12 (1) λ λ 1 2 P P2窄巷伪距: Pa = ( 1 + )λa = ρ * + I12 (2) λ λ 1 2
式中,下标d, a, 1, 2分别表示宽巷、窄巷、L1、L2;Λ为波长;P为伪距;I12为电离层延迟;Ρ*为包
括所有与频率无关误差的距离项。
L1/L2双频载波相位测量值构成的宽巷和窄巷
观测值分别为
宽巷载波相位: (Φd + Nd )λd = ρ * + I12 (3)
窄巷载波相位: (Φa + Na )λa =ρ* − I12 (4)
式中,Φ表示载波相位测量值;N表示整周模糊度,其他各项意义同前。
比较式(2)和式(3),可以发现两式右端相同,
将两式相减并取多时元观测值的平均,且取为最接近的整数,可以得到宽巷载波相位测量的整周模糊度: Pa Nd [ ] = −Φ d 平均,取整λ d ( 5 )
由于宽巷载波相位观测量的波长为86cm,而现代的P(Y)码接收机伪距测量精度至少具有亚米级,
因此只需要少数时元的观测值平滑,即可确定Nd。此外,由于式(2)和式(3)中电离层延迟大小,符号
均相同,即使进行长时间平滑也不会引入因电离延迟所引起的误差。
当基准站和流动站之间距离较近时,可以认为差分电离层延迟近似为零,此时可以直接用Φd+nd
来解求窄巷载波相位测量整周模糊度,同样采用多时元观测值的平均,且取最接近整数的形式: λ Na = [ d (Φd + Nd ) −Φ ] a 平均,取整λ a
(6)
式(6)在基准站和流动站之间距离较远或电离
层较活跃,差分电离层延迟较显著时不适用。此时Na可由式(1)和式(4)相减,并取多时元观测值的
平均,且取最接近整数而得到 Pd Na [ ] = −Φ a 平均,取整λ a
(7)
对于一些可以减少多路径效应的新型GPS信号接收机,利用10~15分钟的观测值由式(7)可以
确定窄巷载波相位测量整周模糊度;若已经确定宽巷整周模糊度,则根据宽巷、窄巷整周模糊度的同奇、同偶的特性,确定窄巷整周模糊度的有效波长将扩大一倍,反之亦然。根据这一关系,德国学者Wubbena博士提出了所谓的超宽巷技术(extra Widelaning)。
在确定Nd, NA后,可以求得L1/L2观测值Φ1/Φ2
的整周模糊度: 1 N1 Nd N = + 2 a 1 N2 Na Nd = − 2
(8)
双频P码伪距法的优点在于只需利卫星的观测数据且与其他卫星无关,计算十分简单,当采用式(6)和式(7)分别计算Nd和Na时,整周模糊度的求解与测站站间距离无关,可用于数千公里相对定位的整周模糊度求解。其基本要求是,需要高精度的双频伪距观测量;当21世纪初叶GPS现代化后,这对广大GPS用户是不难达到的。
2 模糊度函数法
模糊度函数最早由美国学者Couselman教授于1981年提出,其后又由美国Remondi博士于1984年引入静态定位的数据处理。Remondi博士(1991)和Mader(1992)博士最早将这种方法用于在航模糊度解算。
模糊度函数定义为: k j−1 l { ( ) ∑ ∑∑exp ⎡Δ∇Φ ( Z0) A X ,Y , Z i2 jkl X 0,Y0, = π ⎣ obs k =1 j=1 l=1 jkl ( X ,Y ,Z ) ⎦} ⎤ − Δ∇Φ (9) cal
式中,K表示时元数;J为每个时元同步观测的GPS
卫星数; l是观测所用频率个数; ∇Δφobs jkl ( X 0Y0 Z0 )是动态测站正确位置为(X0Y0Z0)处的载波相位测量双
差观测值; ∇ΔΦ jkl ( XYZ )是由某一个检测点(XYZ) cal
计算出来的载波相位测量双差计算值。对于单个频
率、单个双差和单个观测时元,在(∇ΔΦ −∇ΔΦ ) 为obs cal 0或整周数时,模糊度函数具有最大值1。该最大值出现在检测点X,Y,Z=XO,YO,ZO (假设卫星轨道误差及其他传播延迟误差、量测误差等均为零)以及
(∇ΔΦ −∇ΔΦ )为整周数的检测点处。显然,对于单obs cal时元和单频率的单个双差观测值对应的模糊度函数最大值为1的检测点不是惟一的。为了克服这一现象,对于每个检测点都应利用来自所有观测的卫星,不同的观测频率以及不同的观测时元的双差载波相位观测值组合起来计算该点的模糊度函数值。正如式(9)所示,当有足够多的、来源不同的观测量组合起来,则除了在正确动态站位置(X0,Y0,Z0) 处,模糊度函数值仍始终保持最大值外,在其他的检测点上出现最大值的现象将逐个消失。
图1表示模糊度函数法的解算原理;它主要依如下三步确定整周模糊度:
2.1 估计整周模糊度的初始值
现代的GPS信号接收机,都能在作载波相位测量的同时,用伪噪声码进行伪距测量。特别是2005年开始发射的新型GPS卫星,都将增设第三个民用导航定位信号;以致GPS卫星的第一、第二和第三用
导航定位信号的载波频率分别为
f L1=1575.42MHZ,
f L2=1227.60MHZ,
f L5=1176.45MHZ。
当用调制在上列三个载波上的一个伪噪声码测
得三个伪距(PL1、PL2和PL5)时,则可按下列算式求得载波相位测量整周模糊度的初始值(N0 i):
式中,N01、N02和N03分别表示第一、二、三载波(L1,L2,L5)的整周模糊度的初始值。当考虑到f 1=1575.42Hmz,f 2=1227.60HMZ和F5=1176.45HMZ时,算得