CN102636800B - 一种基于整理预处理的gnss整周模糊度解相关方法 - Google Patents

Abstract

一种基于整理预处理的GNSS整周模糊度解相关方法，其特征是在基于LDLT分解的解相关处理中，先对模糊度协方差矩阵根据解相关后对角线元素大小进行整理预处理，再进行解相关处理。方法以迭代的方式进行。在每次迭代处理中以行数增大的顺序先用设计的公式计算LDLT分解中L和D矩阵中备选元素值，然后进行备选元素比较和筛选方法，确定主选元素，并记录主选元素序号，再利用设计的调整方法对协方差矩阵进行预处理整理，最后进行解相关运算；各行都处理后，得到一次迭代转换后的协方差矩阵。迭代处理一直进行到协方差矩阵解相关程度不再提高为止。本发明能达到比逆整数乔列斯基算法更好的解相关程度。

Description

一种基于整理预处理的GNSS整周模糊度解相关方法

技术领域

本发明涉及一种基于整理预处理的GNSS整周模糊度解相关方法，用于全球卫星导航系统。

背景技术

全球卫星导航系统(Global Navigation Satellite System，GNSS)飞速发展，一方面正以难以置信的速度改变着航天、航空、陆地、海洋等运动载体的导航方式和导航精度；另一方面也迅速地改变着对物体的测量方式和人们的生活方式。在系统设计水平(卫星钟差、卫星星历预报误差、观测频点)一定的前提下，基于载波相位的差分定位方式是提高GNSS的定位精度最好的解决办法，得到了广泛的应用。RTK(Real-time kinematic)和利用GNSS差分载波相位进行载体姿态测量是载波相位差分方式的两个典型应用方向。它们都涉及双差载波整周模糊度的动态确定。双差整周模糊度的解算非常困难，尤其当注意到其数值计算效率的时候。尽管国内外提出了很多方法，目前为止，还没有一种公认是最好的快速整周模糊度的解算方法，快速解算双差整周模糊度至今仍然是GNSS领域的一个研究热点和难题。Teunissen于1993年提出的LAMBDA(least-squares ambiguity decorrelation adjustment)算法是目前动态GNSS模糊度解算效果最好的方法之一。该算法的关键步骤之一是解相关处理(又称Z变换)，即通过变换矩阵对描述模糊度之间相关关系的协方差矩阵进行解相关处理，使模糊度协方差矩阵尽可能对角化。它直接影响到整周模糊度搜索的速度和成败。

国外提出以下很多典型的解相关算法，其中逆整数乔列斯基算法是利用模糊度协方差矩阵

的对称性和正定性对

进行LDL^T分解，对L进行取整、求逆变换后，对模糊度协方差矩阵

进行变换，使得变换后的相关矩阵Qu趋近于对角阵，用公式表示为

由于对L矩阵进行了取整，变换后得到的矩阵D不再是对角阵，因此变换需要迭代进行。

然而逆整数乔列斯基算法在解相关程度的性能并不理想，本发明的方法对其进行改进。本发明的方法基于逆整数乔列斯基算法对矩阵元素预先进行排序整理，排序预整理可以提高算法的解相关程度。这是由于LDL^T分解具有随着行数增加，对角线元素需减去更大的值才能进行分解的特点。所以按照对角线元素大小的预排序过程可以使LDL^T分解后D的对角线元素小于零的概率减小、转换后对角线元素之间差距减小，从而提高算法的解相关程度和预防病态Z变换的能力。本发明考虑以LDL^T分解后对角线的大小进行对角线元素排序，能达到比逆整数乔列斯基算法更好的解相关程度。

发明内容

本发明解决的技术问题是：克服现有技术的不足，提供一种基于整理预处理的GNSS整周模糊度解相关方法，可以提高解相关程度性能。

本发明技术解决方案的特点在于包括下列步骤：

1)假设

按照如下步骤①～③进行LDL^T分解在第i(i＝1，2，...，n)行的整理预处理过程。

①首先用设计的计算公式计算LDL^T分解中L和D矩阵中备选的元素值h_pj和s_p(p＝i，i+1，...，n，j＜i)；其中，L是单位下三角阵，D是对角阵。注意：若矩阵元素下标号小于零，则此元素不存在，不用计算。

②然后用备选元素比较和筛选方法，确定主选元素s_i，并记录主选元素序号M_i。

③再利用设计的调整方法对模糊度协方差矩阵进行处理。

Q_i＝S_i(i，M_i)Q_i-1S_i(i，M_i)^T。

其中，S_i(i，M_i)是第i次变换的初等置换矩阵。

2)步骤1)结束后，得到矩阵

和D＝diag(d₁，d₂，…，dn)。

3)将得到的矩阵L和D结合得到的主选元素序号综合对模糊度协方差矩阵进行转换，得到一次迭代处理后的模糊度协方差矩阵。

4)以上步骤1)～3)循环执行，直到L为单位阵或者解相关程度不再向好的方向进行时为止。得到最终转化的值即是解相关后的模糊度协方差矩阵。

所述步骤1)中LDL^T分解中L和D矩阵中备选元素值h_pj和s_p(p＝i，i+1，...，n，j＜i)的计算公式为：

$h_{\mathrm{pj}}=\frac{a_{\mathrm{pj}}-\underset{k=1}{\overset{j-1}{\mathrm{\Σ}}}{h}_{\mathrm{pk}}{d}_k{h}_{\mathrm{jk}}}{d_j},$

$s_p=a_{\mathrm{pp}}-\underset{k=1}{\overset{i-1}{\mathrm{\Σ}}}{l}_{\mathrm{pk}}{d}_k{l}_{\mathrm{pk}}$

所述步骤1)中确定主选元素的方法为：

$s_i=\underset{i\≤p\≤n}{\min }\left\{s_p\right\}$

令M_i的值为使s_p取得最小值的p的值。那么，

$l_{\mathrm{ij}}=h_{M_{i}j}$

$d_i=s_{M_i}$

所述步骤1)中设计的对模糊度协方差矩阵进行处理方法为Q_i＝S_i(i，M_i)Q_i-1S_i(i，M_i)^T。

其中，S_i(_i，M_i)是第i次变换的初等置换矩阵，S_i(i，M_i)是将单位矩阵I_n×n的第i行和第M_i行进行交换得到的。I_n×n是n×n单位矩阵。显然，i＝M_i时，

S_i(i，M_i)＝I_n×n。

所述步骤3)中对模糊度协方差矩阵进行转换，得到一次迭代处理后的模糊度协方差矩阵的具体步骤为：

先计算利用记录的M_i序号构造矩阵

S＝S_n(n，M_n)S_n-1(n-1，M_n-1)·····S₁(1，M₁)，

然后计算

$Q_{\hat{u}}=\left[L^{-1}\right]{\mathrm{SQ}}_0{S}^T{\left[L^{-1}\right]}^T.$

其中[]表示取括号内实数最近的整数值。

本发明的原理是：逆整数乔列斯基算法是利用模糊度协方差矩阵

的对称性和正定性对

进行LDL^T分解，对L进行取整、求逆变换后，对模糊度协方差矩阵

进行变换，使得变换后的相关矩阵Q_u趋近于对角阵。由于取整运算，变换的对角化能力被削弱，逆整数乔列斯基变换需要迭代进行。

LDL^T分解具有随着行数增加，对角线元素需减去更大的值才能进行分解的特点。所以在基于LDL^T分解的解相关处理中，考虑对角线元素的变化情况，实时对对角线上的元素进行预计算、比较，并通过置换矩阵对矩阵进行整理，向使得LDL^T分解后D的对角线元素出现小于零的概率减小、对角线元素之间差距减小的方向调整和处理矩阵，能够提高算法的解相关程度(包括对角线元素按照大小顺序排列)。

本发明与现有技术相比的优点在于：本发明考虑以LDL^T分解后对角线的大小进行对角线元素排序，找到转换后最小对角线元素，并将此对角线元素转换到对应位置，本发明能够达到比逆整数乔列斯基算法更好的减小对角线元素之间差距的效果，可以提高算法解相关程度。

附图说明

图1为本发明算法的解算流程图。

具体实施方式

在GNSS定位应用中，可以利用常规最小二乘方法得到整周模糊度a的浮点解

和模糊度浮点解的协方差矩阵

由于利用模糊度的整数特性可以提高精密导航信息的估计精度。所以希望确定整周模糊度的固定解

一×般固定解解算没有统一的方法，实现起来相对困难一些，通常利用搜索的方法求解模糊度整数解，模糊度的搜索空间可以定义为：

${(\hat{a}-a)}^T{Q}_{\hat{a}}^{-1}(\hat{a}-a)<{\mathrm{\&chi;}}^2.---\left(1\right)$

其中，χ是决定搜索空间大小的常数。

模糊度解相关处理的过程表述为其中。矩阵Z是转换矩阵，

表示转换后模糊度z的协方差矩阵，z＝Za；在Z变换之后，对式(1)的整数搜索过程就变换为对

的整数搜索。协方差矩阵的解相关程度越高，搜索效率越高。现在需要解决的关键问题是如何选择合适的矩阵Z使变换后的解相关程度尽量高。

下面介绍解相关处理的具体实施步骤。

模糊度浮点解的协方差矩阵

一般是对称正定矩阵，设整周模糊度未知数个数是n，那么是n×n维，

a_ij表示第i与第j个整周模糊度的自相关(i＝j)或互相关(i≠j)。

以行数由小到大的次序迭代进行对

的解相关处理。

1)令中间变量

按照如下步骤①～③进行LDL^T分解在第i(i＝1，2，...，n)行的整理预处理过程。

①首先用下面的公式计算LDL^T分解中L和D矩阵中备选的元素值h_pj和s_p(p＝i，i+1，...，n，j＜i)；

$h_{\mathrm{pj}}=\frac{a_{\mathrm{pj}}-\underset{k=1}{\overset{j-1}{\mathrm{\&Sigma;}}}{h}_{\mathrm{pk}}{d}_k{h}_{\mathrm{jk}}}{d_j},$

$s_p=a_{\mathrm{pp}}-\underset{k=1}{\overset{i-1}{\mathrm{\&Sigma;}}}{h}_{\mathrm{pk}}{d}_k{h}_{\mathrm{pk}}.$

其中，L是单位下三角阵，D是对角阵，h_pj是L矩阵中备选的元素值，下标表示第p行、第j列的元素，由于L是单位下三角阵，p≥j的时候h_pj可能有非零值，p＝j时h_pj＝1。s_p是D矩阵中对角线上第p个元素的备选值，d_i是对角线上第i个元素的确定值。

注意：若矩阵元素下标号小于零，则此元素不存在，不用计算。

②然后用备选元素比较和筛选方法，确定主选元素s_i，并记录主选元素序号M_i。

$s_i=\underset{i\&le;p\&le;n}{\min }\left\{s_p\right\}$

令M_i的值为使s_p取得最小值的p的值，那么，

$l_{\mathrm{ij}}=h_{M_{i}j}$

$d_i=s_{M_i}$

s_i、M_i中的下标i表示在第i(i＝1，2，...，n)行的处理。

③再利用设计的调整方法对模糊度协方差矩阵进行处理

Q_i＝S_i(i，M_i)Q_i-1S_i(i，M_i)^T。

其中，S_i(i，M_i)是第i次变换的初等置换矩阵，S_i(i，M_i)是将单位矩阵I_n×n的第i行和第M_i行进行交换得到的。I_n×n是n×n单位矩阵。显然，i＝M_i时，S_i(i，M_i)＝I_n×n。

2)步骤1)结束后，得到矩阵

和D＝diag(d₁，d₂，…，d_n)。

其中，diag(·)表示以括号中元素为对角线元素组成对角阵。

3)将得到的矩阵L和D结合得到的主选元素序号综合对模糊度协方差矩阵进行转换，得到一次迭代处理后的模糊度协方差矩阵。

先计算利用记录的M_i序号构造矩阵

S＝S_n(n，M_n)S_n-1(n-1，M_n-1)·····S₁(1，M₁)

然后计算

$Q_{\hat{u}}=\left[L^{-1}\right]{\mathrm{SQ}}_0{S}^T{\left[L^{-1}\right]}^T.$

其中，[]表示取括号内实数最近的整数值。

最后令

以方便进行循环处理。

4)以上步骤1)～3)循环执行，直到L为单位阵或者解相关程度不再向好的方向进行时为止。得到最终转化的

值即是解相关后的模糊度协方差矩阵。

解相关程度一般指转化后协方差矩阵的对角化程度，对角线元素之间差距越小，非角线元素越小，解相关程度越高，搜索效率就越高。值得一提的是对角线元素按照大小顺序排列有助于搜索效率的提高，把它也列入解相关程度的评价内容之一。

解相关程度的衡量常用计算协方差矩阵的条件数的方法。协方差矩阵条件数定义为矩阵最大与最小特征值之比，即c＝λ_max/λ_min。一般来讲，条件数越小，说明协方差矩阵对角化程度越高。对角线元素按照大小顺序排列有助于搜索效率的提高，把它也列入解相关程度的评价内容之一。

为了说明本发明比逆乔列斯基解相关算法在提高解相关程度方面的优势，列举试验说明。

用随机仿真试验产生大量的对称正定的随机矩阵。产生200个随机正定矩阵，分别用逆整数乔列斯基算法和对角线预排算算法进行解相关，然后统计200个原始随机正定矩阵条件数的平均值和用两种方法解相关后矩阵条件数的平均值、并统计解相关程度得到改善和没有得到改善的矩阵的个数，注意试验中用条件数不断减小作为循环迭代的条件。表1给出了逆整数乔列斯基算法和本发明提出的方法的条件数指标统计结果比较。注意本发明提出的方法能够取得使对角线元素按照大小顺序排列的效果，而逆整数乔列斯基算法不能达到使对角线元素按照大小顺序排列的效果。试验结果表明本发明提出的方法的解相关程度比逆整数乔列斯基算法高，取得解相关程度提高的矩阵概率比逆整数乔列斯基算法的大。

表1逆整数乔列斯基算法和本发明算法的结果比较

算法	原始平均值	Nolma	解相关后平均值	NoEqb
					逆整数乔列斯基算法	1039.1	148	232.7	52
本发明提出的方法	1039.1	192	47.4	8

a Nolm：条件数较小的例子(解相关程度提高)。

b NoEq：条件数没有减小的例子(解相关程度没有变化)。

本发明未详细阐述部分属于本领域公知技术。

Claims (1)

1.一种基于整理预处理的GNSS整周模糊度解相关方法，包括进行LDL^T分解在第i行的预处理整理过程，得到L和D矩阵；将得到的矩阵L和D结合得到的主选元素序号综合对模糊度协方差矩阵进行转换，得到一次迭代处理后的模糊度协方差矩阵，其特征在于具体步骤如下：

（1）按照如下步骤①～③进行LDL^T分解在第i行的预处理整理过程，其中i＝1,2,...,n：

①首先用下面的公式计算LDL^T分解中L和D矩阵中备选的元素值h_pj和s_p，p＝i,i+1,...,n，j＜i；

$h_{\mathrm{pj}}=\frac{a_{\mathrm{pj}}-\underset{k=1}{\overset{j-1}{\mathrm{\&Sigma;}}}{h}_{\mathrm{pk}}{d}_k{h}_{\mathrm{jk}}}{d_j},$

$s_p=a_{\mathrm{pp}}-\underset{k=1}{\overset{i-1}{\mathrm{\&Sigma;}}}{h}_{\mathrm{pk}}{d}_k{h}_{\mathrm{pk}}$

其中，L是单位下三角阵，D是对角阵，h_pj是L矩阵中备选的元素值，下标表示第p行、第j列的元素，由于L是单位下三角阵，p≥j的时候h_pj可能有非零值，p＝j时h_pj＝1；s_p是D矩阵中对角线上第p个元素的备选值；d_i是对角线上第i个元素的确定值；若矩阵元素下标号小于零，则此元素不存在，不用计算；

②然后用备选元素比较和筛选方法，确定主选元素s_i，并记录主选元素序号M_i，

$s_i=\underset{i\&le;p\&le;n}{\min \left\{s_p\right\}}$

令M_i的值为使s_p取得最小值的p的值，

$l_{\mathrm{ij}}=h_{M_{i}j}$

$d_i=s_{M_i}$

s_i、M_i中的下标i表示在第i行的处理；

③再利用设计的调整方法对模糊度协方差矩阵进行处理，

Q_i＝S_i(i,M_i)Q_i-1S_i(i,M_i)^T

其中，S_i(i,M_i)是第i次变换的初等置换矩阵，S_i(i,M_i)是将单位矩阵I_n×n的第i行和第M_i行进行交换得到的，I_n×n是n×n单位矩阵，i＝M_i时，S_i(i,M_i)＝I_n×n；

（2）步骤（1）结束后，得到矩阵

和D＝diag(d₁,d₂,…,d_n)，

其中，diag(·)表示以括号中元素为对角线元素组成对角阵；

（3）将得到的矩阵L和D结合得到的主选元素序号综合对模糊度协方差矩阵进行转换，得到一次迭代处理后的模糊度协方差矩阵，具体为：

先计算利用记录的M_i序号构造矩阵

S＝S_n(n,M_n)S_n-1(n-1,M_n-1)·····S₁(1,M₁)

然后计算原始模糊度协方差矩阵经解相关后得到的新的模糊度协方差矩阵

$Q_{\hat{u}}=\left[L^{-1}\right]S{Q}_0{S}^T{\left[L^{-1}\right]}^T,$

其中，[]表示取括号内实数最近的整数值；

最后令中间变量

以方便进行循环处理；

（4）以上步骤（1）～（3）循环执行，直到L为单位阵或者解相关程度不再向好的方向进行时为止，得到最终转化的

值即是解相关后的模糊度协方差矩阵。

Images