CN110907975B - 一种基于序贯最小二乘的模糊度固定方法 - Google Patents

Abstract

本发明一种基于序贯最小二乘的模糊度固定方法，能够提高在动态环境下使用单频或双频导航卫星信号实时测出两个天线间模糊度的成功率，为后续的高精度基线测量、姿态确定等提供数据支持。在模糊度实数解估计过程中引入先验信息作为状态新息传递给滤波器，在减小设计矩阵病态性的同时加速收敛。本发明在序贯条件下应用基于测量误差建模理论和矩阵三角化分解技术的思想，实现了累加测量新息的同时不增加额外的计算负担，通过构建最小二乘形式的模糊度滤波器实现对模糊度实数解和协因数矩阵信息的计算，最后利用LAMBDA法对模糊度实数解信息进行处理完成模糊度整数解的搜索，进而解决了传统方法中导航卫星信号在动态环境下求解模糊度成功率不高的问题。

Description

一种基于序贯最小二乘的模糊度固定方法

技术领域

本发明为一种基于序贯最小二乘的模糊度固定方法，涉及基于卫星导航系统的载波差分定位技术。

背景技术

当前基于北斗/GPS/Galileo卫星导航系统进行高精度相对定位主要通过载波相位差分技术来实现，而差分技术的关键在于能否快速可靠的在航解算(On The Fly,OTF)载波相位的双差整周模糊度(以下简称为模糊度)。对于双差方程中的模糊度而言，当不考虑其整数性质时，求得的解常称为实数解(或浮点解)；当顾及其整数性质时，求得的解称为整数解(或固定解)，经典的最小二乘方法只能求出模糊度的实数解，而模糊度的整数解的求取方法大致可分为如下四类：一般的最小二乘搜索法、模糊度函数法、模糊度协方差法以及直接取整法。目前，国际公认最小二乘降相关平差法(LAMBDA)是理论最严密、模糊度搜索成功率最高整周模糊度解算方法。近年来随着小卫星、无人机等相关应用市场的兴起，低成本、低功耗、小体积的单频接收机越来越多的用于相对定位(包括测姿、测向)过程中，其中单频模糊度的求解是该应用过程中的核心话题。

采用载波相位观测量过程中，由于其波长和码精度的限制，使双差观测方程的法方程(以下简称为法方程)的病态性很严重，系数矩阵求逆会不稳定，又由于观测噪声不可避免，采用现有的在航解算模糊度方法会导致模糊度实数解与准确值之间的偏差较大，即使使用LAMBDA方法也很难搜索出正确的模糊度值。

发明内容

本发明的目的在于克服现有技术的不足，提供了一种基于序贯最小二乘的模糊度固定方法，该方法能够确保模糊度实数解及其协方差阵的精度和收敛速度，而且可以减小模糊度搜索空间，从而减少了LAMBDA方法搜索的难度，从而达到提高模糊度固定速度和可靠性的目的。

本发明的上述目的通过以下技术方案实现：

一种基于序贯最小二乘的模糊度固定方法，包括以下步骤：

1)地面接收机利用位于地面不同位置的天线A和主天线B捕获跟踪卫星信号；

2)在当前第n个历元时刻，根据天线A和天线B同时捕获跟踪到的P+1 颗卫星相对于天线B的方向向量，确定第n个历元时刻对应的载波设计矩阵H_B、双差伪距观测矩阵

和载波相位观测矩阵

具体为：

对应第j颗卫星与第1～p颗卫星的伪距双差测量值；

对应第j颗卫星与第1～p颗卫星的载波相位双差测量值；

对应第1～p颗卫星相对于天线B的方向向量，

为第j颗卫星相对于天线B的方向向量；

其中，

为第p颗卫星相对于主天线B的方向向量，p∈[1，P]；

为接收机测量得到的第j颗卫星和第k颗卫星到天线A和主天线B的伪距双差测量量，所述伪距双差测量量根据伪距测量量确定；

为测量得到的第j颗卫星和第k颗卫星到天线A和主天线B的载波双差测量量，所述载波双差测量量根据载波测量量确定，j∈[1,P]，k∈[1,P]，且j≠k；

为双差测量的线性化初值；

为天线A与卫星j的空间距离，

为天线A与卫星k的空间距离，

为天线B与卫星j的空间距离，

为天线B与卫星k的空间距离，

均由系统外部直接输入获得；

3)提取n-1历元时刻的双差模糊度的实数解

以及所述实数解的协因数矩阵

确定测量方程的矩阵形式，具体为：

Z＝A·X；

其中，X为待估参数；n＝1时，

为零；0_p×3为p×3阶的矩阵；E_p×p为p×p 阶的单位矩阵；

为1×3阶的天线B的位置改正向量；

获得测量方程的协因数矩阵∑Z，具体为：

构建测量方程的误差模型：

其中，

为伪距测量误差的方差，

为载波测量误差的方差，∑P和∑Φ均为p×p维矩阵，n＝1时，

为p×p阶的矩阵，矩阵所有元素设置为1×10⁹；

4)对矩阵∑Z使用乔里斯基分解方法，求解设计矩阵R，使矩阵R和矩阵∑Z两者之间的数学关系如下：

∑Z＝(R^-1)(R^-1)^T

5)对步骤3)中测量方程Z＝A·X的等号两边同时左乘矩阵R，整理可得满足最小二乘形式的测量方程：

其中，

6)对步骤5)中的设计矩阵

进行QR分解，使得

其中Q为标准正交矩阵，R为上三角矩阵，对步骤4中的测量方程等号两边左乘矩阵Q^T整理得到：

其中，

此时测量方程的矩阵形式可以表述为：

其中，

表示天线B位置修正量相关的变换向量，

表示双差模糊度变换向量，

表示矩阵三角化后的残差向量；

此时测量方程的设计矩阵为上三角形式的矩阵R，使用简单的数学迭代就可直接计算出待估参数

和

的实数解，模糊度

的协因数阵

可以表示为：

7)将模糊度实数解

和协因数阵

送入LAMBDA算法进行模糊度的整数搜索以及所述整数解的可靠性判据；

和

为先验信息被引入下一历元时刻的解算中。

本发明与现有技术相比有益效果为：

1)本发明在模糊度滤波计算过程中，将模糊度随时间不变的特性作为滤波器的先验信息，降低了滤波方程的奇异性，提高了载波相位模糊度实数解以及实数解协因数阵的精度和收敛速度，能够更快地完成载波相位模糊度的分解、搜索以及整数估计，实现快速的高精度相对定位。

2)本发明在模糊度实数解及其协因数阵求解过程中，应用基于测量误差建模理论和矩阵三角化分解技术的思想，实现了累加测量新息的同时不增加额外的计算负担，避免了传统方法中高阶矩阵求逆运算带来的庞大计算量。

附图说明

图1为GNSS接收机测量示意图；

图2为本发明中应用的双差观测值的示意图。

具体实施方式

本发明提供的一种基于序贯最小二乘的模糊度固定方法，可以应用于基于码分多址的卫星导航系统，如北斗B1、B2、B3频点信号，GPS L1、L2频点信号，以及Galileo E1、E5、E6频点信号，进行单频、双频、三频形式的模糊度在航解算。

其中，利用载波相位测量值进行相对定位时，如果采用LAMBDA方法进行整周模糊度解算，则需要满足以下几点要求：一、需要确保模糊度实数解的稳定性和计算效率，也就是要尽可能减小系数矩阵发生奇异的可能性，并降低计算矩阵的维度；二、需要确保模糊度实数解及其协因数阵的精度和收敛速度。本发明采用如下措施达到以上两点要求：首先在模糊度实数解求解过程中，引入了双差伪距测量信息、双差载波相位测量信息，在增加观测方程个数的基础上可以有效降低模糊度滤波器中系数矩阵发生奇异性的可能，有效增强滤波器抗发散的能力；另外本发明在每个历元根据模糊度的不变特性通过引入前一历元的模糊度实数解信息作为先验信息，使之作为一个“伪测量”辅助解算当前模糊度的实数解，可以有效提高模糊度实数解及其协因数阵的精度和收敛速度。针对高阶矩阵的数学计算方面引入基于测量误差建模理论和矩阵三角化分解技术的思想构建模糊度分解技术，规避了传统方法中模糊度实数解计算过程中涉及的高阶矩阵求逆运算，提高了滤波器的计算效率。

采用以上技术措施后，相对于现有技术来说，本发明计算得到的模糊度搜索空间和产生的计算负担会明显减小，进而减小了LAMBDA方法搜索的难度。

本发明的基于序贯最小二乘的模糊度搜索方法，其原理推导过程如下：

如图1所示，GNSS接收机测量的是导航卫星与接收天线在空间通视方向上的信号传输时延，本质上体现的是两者间相对运动状态在通视方向的投影。由于GNSS接收机导航解算结果存在一定程度的误差，其原始观测量的测量值与理论真值之间存在以下关系：

其中，理论真值

代表卫星S^j与天线A的真实几何距离；卫星j与天线A 之间测量的几何距离

S^j为卫星j的位置坐标，一般由广播星历计算得到；

为天线A的初始位置坐标，一般通过导航解算得到；卫星j相对于天线A的方向向量

为天线A初始位置的一阶修正值。

同理可得，主天线B的测量值与理论真值的关系如下所示：

如图2所示，卫星导航系统进行载波相位差分定位时，接收机利用设置在不同位置的两个天线接收卫星信号，分别为天线A和主天线B，对应的测量接收机分别为接收机A和接收机B。采用站际-星间形式的双差模型的双差真值

具体表示为

其中，

和

分别为卫星j和卫星k与天线A的真实几何距离；

和

分别为卫星j和卫星k与主天线B的真实几何距离；将四个线性表达式代入公式 (3)，可得：

由于差分的数学模型在缺少外界辅助信息的情况下只能得到高精度的相对位置，绝对位置修正的精度不高，因此本发明中不修正天线A的位置即认为天线A的位置没有误差，对于式(4)可令

此时式(3)可整理为：

其中，第j颗卫星和第k颗卫星对应的星地距离双差观测量

GNSS接收机的测量分为伪距和载波相位两种，两者的主要区别在于载波相位测量存在模糊度。

假设，在某一历元时刻，GNSS接收机观测P+1颗卫星时，选取其中一颗作为主星，可建立P个载波相位双差方程和P个伪距双差方程，即：

其中：伪距双差测量量

和

分别为接收机A测量的卫星j和卫星k伪距测量值，

和

分别为接收机B测量的卫星j和卫星k伪距测量值；载波相位双差测量量为

，

和

表示接收机A测量的卫星j和卫星k的载波相位测量值，

和

表示接收机B测量的卫星j和卫星k的载波相位测量值；

和

分别表示伪距和载波相位双差测量量的噪声；λ为载波相位的波长；

为双差模糊度，具有整数特性。

将式(5)代入式(6)，可得：

式(7)中，方程左边的为已知项，可以根据解调的电文和定位解算结果计算得到；方程右边的

和λ为已知系数项；

和

为待求项；

和

为双差量的测量噪声，在方程解算过程中，可以忽略不计。

假设，在某一历元时刻，天线A与B共同观测到p+1颗卫星，其中设置卫星编号j＝0作为主星；设置设置卫星编号k＝1～p，对应其他p颗导航卫星作为副星，并且设定

则式(7) 写为如下矩阵形式：

为了简化式(8)，设定如下：

是一个p×1的伪距双差测量的残差向量；

是一个p×1的载波相位双差测量的残差向量；

是一个p×3的设计矩阵；

忽略双差量的测量噪声

则式(8)可简化整理为：

其中，整周模糊度

E_P×P为p阶单位矩阵，0_p×p为p阶零矩阵。

由于在没有周跳的情况下，整周模糊度具有稳定不变的特性，因此滤波器的状态空间的历元更新可以引入整周模糊度的均方根信息，即：

其中，n表示历元个数，

为n-1历元时刻

的协因数矩阵，合并式(8)与(9)，可得：

其中，使得

公式(10)可以简化为：

Z＝A·X (11)

由于双差观测方程具有高度耦合特性，式(10)的设计矩阵A属于稀疏病态阵，对矩阵A的求逆操作需要额外增加大量的计算负担，且计算精度有限。本文采用基于测量误差建模理论和矩阵三角化分解技术避免上述计算负担的增加，提高计算精度。下文对测量误差的建模与矩阵三角化分解技术进行详细描述。

普遍认为，站际-星间双差形式的伪距、载波观测方程的协因数矩阵形式∑P 和∑Φ如下所示：

其中，

为伪距测量误差的方差，

为载波测量误差的方差。因此测量方程(10)的协因数阵∑Z的具体形式如下所示：

结合公式(11)、(13)，针对本发明的测量方程可以完整描述为：

Z＝A·X，∑Z (14)

由于∑P、∑Φ、

均为非对角矩阵，因此测量方程(14)的待估参数使用传统的最小二乘算法进行计算是较为困难且计算负担较大的。最小二乘算法只有在针对零均值高斯白噪声的测量问题才能取得最优解，显而易见公式 (14)不属于这种情况。通过公式(13)的推导可以得知，测量方程的协因数阵为正定对称矩阵，因此使用乔里斯基分解的方法对∑Z进行处理，即

∑Z＝(R^-1)(R^-1)^T (15)

然后公式(14)的等号两边同时左乘矩阵R，公式(14)变化为：

R·Z＝R·A·X，∑R·(R^-1)(R^-1)^T·R^T (16)

其中∑R·(R^-1)(R^-1)^T·R^T＝E，新的观测方程(16)的协因数矩阵为单位阵，满足对其采用最小二乘方法进行参数估计。

为了进一步减小计算量，对公式(16)的做进一步的整理可得：

其中，

E为单位矩阵。

进一步对设计矩阵

进行QR分解，使得

其中Q为标准正交矩阵， R为上三角矩阵，对方程(17)左乘矩阵Q^T可以得到：

其中，

此时公式(18)的矩阵形式可以表述为：

此时测量方程(19)的设计矩阵为上三角形式的矩阵R，使用简单的数学迭代就可直接计算出待估参数

和

的实数解，从而避免了高阶病态矩阵的求逆运算问题，此时根据最小二乘原理，模糊度

的协因数阵可以表示为：

随后模糊度实数解和协因数阵送入LAMBDA算法进行模糊度的整数搜索。

和

为先验信息被引入下一历元时刻的解算中。

根据以上的理论推导过程可知，本发明的提高单频模糊度在航计算成功率的方法，包括如下步骤：

1)地面接收机利用位于地面不同位置的天线A和主天线B捕获跟踪卫星信号；

2)在当前第n个历元时刻，根据天线A和天线B同时捕获跟踪到的P+1 颗卫星相对于天线B的方向向量，确定第n个历元时刻对应的载波设计矩阵H_B、双差伪距观测矩阵

和载波相位观测矩阵

具体为：

对应第j颗卫星与第1～p颗卫星的伪距双差测量值；

对应第j颗卫星与第1～p颗卫星的载波相位双差测量值；

对应第1～p颗卫星相对于天线B的方向向量，

为第j颗卫星相对于天线B的方向向量；

其中，

为第p颗卫星相对于主天线B的方向向量，p∈[1，P]；

为接收机测量得到的第j颗卫星和第k颗卫星到天线A和主天线B的伪距双差测量量，所述伪距双差测量量根据伪距测量量确定；

为测量得到的第j颗卫星和第k颗卫星到天线A和主天线B的载波双差测量量，所述载波双差测量量根据载波测量量确定，j∈[1,P]，k∈[1,P]，且j≠k；

为双差测量的线性化初值；

为天线A与卫星j的空间距离，

为天线A与卫星k的空间距离，

为天线B与卫星j的空间距离，

为天线B与卫星k的空间距离，

均由系统外部直接输入获得；

3)提取n-1历元时刻的双差模糊度的实数解

以及所述实数解的协因数矩阵

确定测量方程的矩阵形式，具体为：

Z＝A·X；

其中，X为待估参数；n＝1时，

为零；0_p×3为p×3阶的矩阵；E_p×p为p×p 阶的单位矩阵；

为1×3阶的天线B的位置改正向量；

获得测量方程的协因数矩阵∑Z，具体为：

构建测量方程的误差模型：

其中，

为伪距测量误差的方差，

为载波测量误差的方差，∑P和∑Φ均为p×p维矩阵，n＝1时，

为p×p阶的矩阵，矩阵所有元素设置为1×10⁹；

4)对矩阵∑Z使用乔里斯基分解方法，求解设计矩阵R，使矩阵R和矩阵∑Z两者之间的数学关系如下：

∑Z＝(R^-1)(R^-1)^T

5)对步骤3)中测量方程Z＝A·X的等号两边同时左乘矩阵R，整理得到满足最小二乘形式的测量方程：

其中，

6)对步骤5)中的设计矩阵

进行QR分解，使得

其中Q为标准正交矩阵，R为上三角矩阵，对步骤5)中的测量方程等号两边左乘矩阵Q^T整理得到：

其中，

此时测量方程的矩阵形式可以表述为：

其中，

表示天线B位置修正量相关的变换向量，

表示双差模糊度变换向量，

表示矩阵三角化后的残差向量；

此时测量方程的设计矩阵为上三角形式的矩阵R，使用简单的数学迭代就可直接计算出待估参数

和

的实数解，模糊度

的协因数阵

可以表示为：

7)将模糊度实数解

和协因数阵

送入LAMBDA算法进行模糊度的整数搜索以及所述整数解的可靠性判据；

和

为先验信息被引入下一历元时刻的解算中。

以上所述，仅为本发明一个具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。

本发明说明书中未作详细描述的内容属于本领域专业技术人员的公知技术。

Claims (1)

1.一种基于序贯最小二乘的模糊度固定方法，其特征在于，包括如下步骤：

1)地面接收机利用位于地面不同位置的天线A和主天线B捕获跟踪卫星信号；

2)在当前第n个历元时刻，根据天线A和天线B同时捕获跟踪到的P+1颗卫星相对于天线B的方向向量，确定第n个历元时刻对应的载波设计矩阵H_B、双差伪距观测矩阵

和载波相位观测矩阵

具体为：

对应第j颗卫星与第1～p颗卫星的伪距双差测量值；

对应第j颗卫星与第1～p颗卫星的载波相位双差测量值；

对应第1～p颗卫星相对于天线B的方向向量，

为第j颗卫星相对于天线B的方向向量；

其中，

为第p颗卫星相对于主天线B的方向向量，p∈[1，P]；

为接收机测量得到的第j颗卫星和第k颗卫星到天线A和主天线B的伪距双差测量量，所述伪距双差测量量根据伪距测量量确定；

为测量得到的第j颗卫星和第k颗卫星到天线A和主天线B的载波双差测量量，所述载波双差测量量根据载波测量量确定，j∈[1,P]，k∈[1,P]，且j≠k；

为双差测量的线性化初值；

为天线A与卫星j的空间距离，

为天线A与卫星k的空间距离，

为天线B与卫星j的空间距离，

为天线B与卫星k的空间距离，

均由系统外部直接输入获得；

3)提取n-1历元时刻的双差模糊度的实数解

以及所述实数解的协因数矩阵

确定测量方程的矩阵形式，具体为：

Z＝A·X；

其中，X为待估参数；n＝1时，

为零；0_p×3为p×3阶的矩阵；E_p×p为p×p阶的单位矩阵；

为1×3阶的天线B的位置改正向量；0_p×p为p阶零矩阵；

获得测量方程的协因数矩阵∑Z，具体为：

构建测量方程的误差模型：

其中，

为伪距测量误差的方差，

为载波测量误差的方差，∑P和∑Φ均为p×p维矩阵，n＝1时，

为p×p阶的矩阵，矩阵所有元素设置为1×10⁹；

4)对矩阵∑Z使用乔里斯基分解方法，求解设计矩阵R，使矩阵R和矩阵∑Z两者之间的数学关系如下：

∑Z＝(R^-1)(R^-1)^T

5)对步骤3)中测量方程Z＝A·X的等号两边同时左乘矩阵R，整理得到满足最小二乘形式的测量方程：

其中，

6)对步骤5)中的设计矩阵

进行QR分解，使得

其中Q为标准正交矩阵，R为上三角矩阵，对步骤5)中的测量方程等号两边左乘矩阵Q^T整理得到：

其中，

此时测量方程的矩阵形式表述为：

其中，

表示天线B位置修正量相关的变换向量，

表示双差模糊度变换向量，

表示矩阵三角化后的残差向量；

此时测量方程的设计矩阵为上三角形式的矩阵R，使用简单的数学迭代就可直接计算出待估参数

和

模糊度实数解

的协因数阵

表示为：

7)将模糊度实数解

和协因数阵

送入LAMBDA算法进行模糊度的整数解搜索获得整数解及整数解的可靠性判据；

和

为先验信息被引入下一历元时刻的解算中。

Images