37 baidu 差分GPS载波相位测量整周模糊度的快速求解

胡国辉孟浩袁信摘要：对Cholesky分解整周模糊度的求解进行了改进，在求解整周模糊度的过程中，首先采用LAMBDA法对整周模糊度进行整数线性变换再作Cholesky分解，然后利用最优剪枝法（best cut）对整周模糊度进行搜索，实验结果表明该方法具有快速搜索整周模糊度的能力，可以满足采用GPS载波相位测量确定姿态以及GPS载波相位测量与INS组合的实时性。

关键词：导航整周模糊度载波相位Cholesky分解中图分类号：V241.5FAST CARRIER PHASE AMBIGUITY RESOLUTION FORDIFFERENCE GPSHu Guohui1, Meng Hao2, Yuan Xin11(Department of Automatic Control, Naijing University of Aeronautics & Astronautics,Nanjing,210016)2(Department of Automatic Control, Harbin EngineeringUniversity,Harbin,150001)Abstract The paper presents a new development method for Cholesky ambiguity search method. The method makes use of an ambiguity reparametrization, Cholesky decomposition and best cut. Experiment results show that the method can achieve fast search ability, and satisfy real time attitude determination and GPS/INS integration with GPS carrier phase measurement.Key words navigation, ambiguity, carrier phase, Cholesky factorization单纯采用Cholesky分解整周模糊度的求解［1］往往搜索次数较多，采用LAMBDA法［2］对整周模糊度进行整数线性变换再作Cholesky分解，使变换后的整周模糊度方差更小，有效的提高了搜索速度，实验结果表明该算法能快速确定整周模糊度，能满足采用载波相位的姿态确定以及与惯导组合着陆的实时性要求。

1 整周模糊度的求解对于双天线GPS载波相位测量，系统的状态方程和观测方程为（1）式中：I i×j ， 0i×j 分别为i×j 单位矩阵、零矩阵；X x ＝（x y z)为基线向量和基线速度；X N ＝（N 1 N 2 … N n ）为双差模糊度；y k 为载波相位的双差测量值；为指向卫星的方向余弦矩阵。

对于式（1）的状态方程和观测方程可以采用卡尔曼滤波器进行估计。

k ＋1／k ＝Ak ／k（2）P k ＋1／k ＝AP k ／k A T ＋Q k （3）k ＋1／k ＋1＝k ＋1／k ＋K k ＋1［y k ＋1－H k ＋1k ＋1／k］ （4）K k ＋1＝P k ＋1／k H Tk ＋1［H k ＋1P k ＋1／k H Tk ＋1＋R k ＋1］－1 （5） P k ＋1／k ＋1＝（I ＋K k ＋1H k ＋1）P k ＋1／k （I － K k ＋1H k ＋1）T ＋K k ＋1R k ＋1K Tk ＋1 （6）P 0／0，X 0／0初值的选取取决于对基线向量和初始模糊度范围的了解，当卡尔曼滤波器的滤波达到一定精度后，可进行整周模糊度的搜索。

对于式（6）中的P k ／k 可以写成其中：P x ，P N 分别为基线向量与整周模糊度误差的协方差矩阵；P xN ，P Nx 分别为它们之间的误差协方差矩阵。

文献［2］证明了采用LAMBDA 的下三角高斯变换可以将置信椭球变成置信球，变换后参数方差更小。

这里为了使变换更为简单，同时采用整数上三角和下三角高斯变换，对P N 阵的变换步骤为：（1）上三角变换 P N ＝U 1D U1U T1。

其中：U U1为对角阵; U 1为上三角高斯变换阵;（2）计算整数变换阵 Z U1＝［round （U 1）］－1，round 表示取整; （3）计算整数变换后的方差阵 P Z^NU1\}＝Z U1P N Z TU1; （4）上三角变换 PNU1＝L 1D L1L T1。

其中：D L1为对角阵；L 1为上三角高斯变换阵;（5）计算整数变换阵 Z L1＝［round （L 1）］－1; （6）计算整数变换后的方差阵 P Z^NL1＝ Z L1PNU 1Z TL1;（7） 将P NL1代替P N ，重复上述过程，直到整数变换为单位阵。

最后，整数变换阵为Z ＝Z Lk －1Z Uk －1…Z L1Z U1 （7）对于Kalman 滤波器中给出的双差整周模糊度X^N 经过式（7）的整数线性变换得N ＝ZN（8）相应的X^N 状态协方差阵变为PN＝ZP N Z T （9）整周模糊度Cholesky 分解应使 J ＝（N －N ）T P －1N （N －N ）＝（N －N′）T P Z^N －1（N－N′） （10）最小，其中N 为待确定的整周模糊度，N′＝ZN ,取N＝（z 1 z 2 … z n ）；N ′＝（N ′1 N ′2 \: N ′n ）为了使搜索能满足树状规律，可进行Cholesky 分解。

PN －1＝CCT（11）其中：C＝将式（8）带入式（7）可得J＝（N －N′）T－P－1N（ZN－N′)＝（N －N′）TCCT（N－N′）＝［CT （N－N′）］T［CT（N－N′）］＝fTf＝f21＋f22＋…＋f2n（12）式中：f＝CT （N－N′）［f1f2… fn］（13）由于Z^N的误差方差阵为N′的取值范围为N′j ∈［round（zj－3）,round（zj+3）］（j＝1，…，n），如果round（zj －3）＝round（zj+3）,则取Nj ∈［round（zj－3）-1,round（zj+3）+1］（j＝1，…，n）。

计算量由fn →f1递增，由于f2i是大于或等于零，在计算中式（12）和式（13）可以交替计算，当某一检测模糊度f部分分量平方和已大于某一限值（可选为前面已经计算的J次小值），则该模糊度组将被拒绝而不用计算f的剩余分量及其平方和，从而大大减小计算量。

最优剪枝法的搜索为一历元的模糊度搜索顺序应从搜索范围的中间开始向两边进行，这是因为越靠近搜索中间的模糊度组为正确模糊度的可能性越大，其J数值一般较搜索空间的边缘要小，从而在上述搜索计算过程中可以进一步减小计算量。

搜索过程中只要保留已搜索过的模糊度、模糊度组中J为最小及次小的模糊度数值及相应的J数值。

对最后搜索的模糊度需进行ratio检验。

式中：VT R－1kV为模糊度参数估计为实数时的双差载波相位残差平方和，Jk为式（12）中的计算值，当ratio值大于某一限值时（一般取为常值），则认为残差平方和和Jk最小的模糊度组是正确模糊度组。

但仅用单历元的ratio检验来确定的模糊度有时不可靠，尤其是卫星图形较差或观测噪声较大时。

然而随观测历元的增加，通过ratio检验的不正确模糊度组会逐渐被正确模糊度组所取代，ratio检验一段连续检验时间连续大于限值，且Ωk最小的模糊度组为同一模糊度组，则双差整周模糊度可以固定作为已知值。

2 试验结果与分析试验采用2台Motorola VP oncore 接收机做双天线方位测量，基线长度为1.5m，数据是在南京航空航天大学15号楼楼顶某已知点采集的，数据采集的时间是1998年4月9日晚，观测期间共观测到7颗卫星(prn 6,10,13,17,24,26,27),数据采集的频率为1Hz, 载波相位标准差设为1cm, C/A码伪距的标准差设为3cm，ratio检验门限值选为4，OVT 检验选为10个历元，图1(a)为模糊度固定为整数前滤波的基线向量(x y z)和模糊度固定后仅利用该历元观测求解的基线向量，图1(b)为双差观测滤波残差。

图1 基线向量和残差(a))基线向量；(b)残差数据处理在奔腾586采用MATLAB编程。

从Motorola VP oncore的数据处理来看，由于采用LAMBDA法对整周模糊度进行整数线性变换再作Cholesky分解，然后利用最优剪枝法（best cut）对整周模糊度进行搜索，模糊度的求解时间由原来的几十秒缩短到0.50s，可以满足采用GPS 载波相位测量确定姿态以及GPS载波相位测量与INS组合的实时性。

模糊度固定后, 基线向量的精度由分米级上升到厘米级，由于天线接收的多路径误差较大。

根据：由试验结果可以看出，1.5m基线采用该历元观测求解方法的方位精度可达二十几毫弧度。

国防预研基金资助课题作者简介：胡国辉，男，副教授，1966年7月生，1989年7月在哈尔滨船舶工程学院获工学学士学位，1992年3月在该校获工学硕士学位，1995年3月在哈尔滨工业大学获工学博士学位，1996年12月在南京航空航天大学博士后出站留校工作，研究方向为：惯性导航组合导航、多传感器信息融合、GPS理论与应用。

袁信，教授，博士研究生导师，研究领域：导航、制导与控制。

孟浩男，讲师，研究方向：GPS应用、自动驾驶仪、数据融合和自适应控制。

作者单位：胡国辉袁信南京航空航天大学自动控制系, 南京, 210016孟浩哈尔滨工程大学自动控制系，哈尔滨，150001参考文献1 Hatch R, Euler H J. Comparison of several AROF kinematic techniques. Proceeding of the 8th Internatioal Technical Meeting of the Satellite Division of the Institute of Navigation. Palm Springs, California:1995. 363～3702 Teunissen P J G. A new method for fast carrier phase ambiguity estimation. Proceeding of IEEE Position Location and Navigation Symposium. Las Vegas: 1994. 562～573。