CN113670337B - 一种用于gnss/ins组合导航卫星缓变故障检测方法 - Google Patents

Abstract

本发明提供了一种用于GNSS/INS组合导航卫星缓变故障检测方法；该方法采用标准t分布和IGG‑III方案设计自适应增益矩阵以缓解卡尔曼滤波故障跟踪的影响，从而缩短组合导航系统对缓变故障的检测延迟时间，其中自适应增益阵的构造形式简洁，抗差性能高；同时，结合传统外推法的故障检测状态，设计由外推法和残差卡方检验故障检测统计量构成的

统计量，在外推法检测到缓变故障的状态下利用样本分位数原理对

序列进行异常值检测，从而判断缓变故障结束时刻，及时避免虚警现象的发生，提高系统定位精度和可靠性；整个检测方法流程简单，思路清晰。

Description

一种用于GNSS/INS组合导航卫星缓变故障检测方法

技术领域

本发明属于导航定位技术领域，具体涉及一种用于GNSS/INS组合导航卫星缓变故障检测方法。

背景技术

在全球导航卫星系统(GNSS)/惯性导航系统(INS)组合导航系统中，卫星伪距极可能因时钟漂移、卫星轨道建模和电离干扰等因素而发生缓变故障，严重影响系统的可靠性和稳定。因此，研究如何及时、有效地检测缓变故障这一问题具有十分重要的现实意义。

目前用于组合导航卫星故障的检测主要有以下几种方法：残差卡方检测法(RCTM)，多解分离法，自主完好性检测外推法(AIME)和最优故障检测法等。其中，残差卡方检测法只对幅值较大的突变故障检测有效；多解分离法能够有效检测缓变故障，但是计算成本较高；而最优故障检测法工程实现复杂。自主完好性检测外推法的主要思想是利用卡尔曼滤波在外推过程中的新息序列构造故障检测统计量，然后基于假设检验理论判断卫星是否发生缓变故障。有报道指出：相比于其它的故障检测方法，外推法更适合用于组合导航系统的卫星缓变故障检测。

虽然外推法更适应于组合导航卫星缓变故障检测，但是其仍然存在以下两方面的不足：一是卡尔曼滤波器的故障跟踪作用会降低外推法对缓变故障检测的敏感性，从而导致缓变故障的检测延时；另外，由于传统外推法利用卡尔曼滤波在外推过程中的新息序列构造缓变故障检测统计量，因此当缓变故障结束后系统会发生一段时间的虚警现象，导致定位结果精度下降。

为了解决上述检测延迟问题，基于抗差扩展卡尔曼滤波REKF的外推-积累缓变故障检测法和基于标准化新息构造抗差增益阵的抗差卡尔曼滤波算法被相继提出。但是，上述抗差估计中的等效权函数大都是基于正态分布统计量构造的，并没有充分利用多余观测信息，同时等效权函数的临界值过于单一，因此在一定程度上限制了系统的抗差性能。另外，为了解决虚警问题，一种基于层次滤波器的改进型外推-残差卡方检验的联合故障检测算法被提出，其能够更快地检测到缓变故障并准确判定故障结束时间，但是算法实现过程复杂，计算量大。

术语解释：GNSS：全球导航卫星系统(Global Navigation Satellite Systems)；

INS：惯性导航系统(Inertial Navigation System)；

IMU：惯性测量单元(Inertial Measurement Unit)；

AIME：自主完好性监测外推法(Autonomous Integrity MonitoringExtrapolation)；

REKF：抗差扩展卡尔曼滤波(Robust Extended Kalman Filter)；

RCTM：残差卡方检测法(Residual Chi-square Test Method)。

发明内容

为解决上述问题，本发明公开了一种用于GNSS/INS组合导航卫星缓变故障检测方法，流程简单，思路清晰，系统定位精度和可靠性高。

为达到上述目的，本发明的技术方案如下：

一种用于GNSS/INS组合导航卫星缓变故障检测方法，包括如下步骤：

步骤1，对卫星接收机所接收的卫星数据和惯性测量单元IMU的测量数据进行抗差扩展卡尔曼滤波；

步骤2，针对GNSS/INS融合系统，提取新息和其协方差阵以进行卫星缓变故障检测。

步骤1包括：

步骤1-1，将惯性测量单元IMU输出的加速度数据以及角速度数据进行力学编排，得出用户在东、北、天三个方向上的位置，速度和姿态解；

步骤1-2，将惯性导航系统INS推导出的伪距和接收机接收到的伪距数据进行融合，经过抗差扩展卡尔曼滤波REKF，输出用户在东、北、天三个方向上的位置误差，速度误差和姿态误差；

步骤1-3，利用步骤1-2得到的误差状态修正步骤1中的状态解，同时将误差状态反馈到INS中；

步骤1-4，时间更新，回到步骤1-1。

步骤1-2包括：将GNSS/INS组合导航中的INS姿态、速度、位置在东、北、天三个方向上的误差和陀螺、加速度计的在相应方向上的误差以及接收机钟差等效的距离误差和与接收机频率误差等效的距离率误差作为状态量x，将惯性导航系统INS推导出的伪距和接收机接收到的伪距数据作差作为观测量z，抗差扩展卡尔曼滤波包括两个方程：状态方程x_k＝Φ_k/k-1x_k-1+w_k-1和量测方程z_k＝h(x_k)+v_k，其中，x_k是状态向量；Φ_k/k-1表示从历元k-1到k的状态转移矩阵；z_k为观测向量；h(·)表示状态x_k和观测z_k之间的非线性关系；w_k是历元k-1的过程噪声向量；v_k是历元k的观测噪声向量，其中w_k，v_k相互独立，且两者均服从零均值的高斯分布，协方差矩阵分别用Q_k和R_k表示。抗差扩展卡尔曼滤波步骤如下：

步骤1-2-1，初始化，初始化组合导航的状态x₀和误差协方差矩阵P₀；

步骤1-2-2，预测：

其中，符号“^”表示状态的估计；下标k/k-1表示从历元k-1转移到k；状态估计的误差协方差矩阵用P_k表示。

步骤1-2-3，更新：

其中，γ_k表示新息，是计算故障检测统计量的基本元素；

是抗差增益矩阵；/>

是等价权重矩阵；H_k表示非线性矢量函数h(·)的Jacobi矩阵。

另外，新息γ_k的协方差矩阵用符号S_k表示，

步骤1-2-3包括：在原始卡尔曼增益K_k的基础上，利用观测数据构造等价权重矩阵

抗差增益阵与原始增益阵的关系：/>

其中/>

为对角阵；选择IGG-III方案计算该矩阵；依照IGG-III方案的原理，设计了一种基于标准t分布的等价权矩阵计算方法；该抗差方法记为REKF-t；等价权矩阵的对角元素：

式中，临界值T₀,T₁分别取t(m-2)分布显著性水平为α₀,α₁的分位数，其中，α₀,α₁人为设定；t_i是观测值z_i对应新息的标准化t分布统计量，可通过当前历元k的新息向量γ_k和其对应的协方差矩阵S_k计算：

其中，

且m≠1，为标准t分布统计量的系数；/>

为标准化新息，

γ_i是第i个观测量对应的新息，σ_γi表示其对应的标准差，可由新息协方差矩阵S_k获得。在无故障情况下，/>

另外，由t_i的表达式可知，m≠1,2，即在历元k只要可视星数不低于3颗，就能有效使用REKF-t算法。

步骤2包括：

步骤2-1，设定误警率P_fa和外推法中滑窗长度L，并将统计值r_A/R置1；

步骤2-2，确定故障判决门限Th，基于外推法计算缓变故障检测统计值；

步骤2-3，判断卫星是否发生缓变故障并输出。若卫星发生缓变故障，则执行步骤2-4，否则返回步骤2-2；

步骤2-4，选择样本分位数，阈值ξ和窗口长度μ，计算系统的统计量r_A/R；

步骤2-5，对r_A/R序列进行异常值检测，判定缓变故障是否结束。若判断出异常值，即确定缓变故障结束时刻T_e，则执行步骤2-6，否则返回步骤2-2；

步骤2-6，将统计值r_A/R置1，基于外推法计算缓变故障检测统计值置0，使用残差卡方检验法作为缓变故障的检测方法；

步骤2-7，判断此时历元是否处于区间[T_e,T_e+L-1]中，若是，则返回步骤2-6，否则返回步骤2-2。

步骤2-2包括：确定故障判决门限Th，

式中，P_fa表示故障检测的误警概率，F(x|m)表示自由度为m的中心卡方分布累积分布函数，另外基于外推法的缓变故障检测统计值的计算：

其中，

表示平均新息γ_avg协方差阵的逆，L表示新息序列的滑窗长度。

步骤2-3包括：利用二元假设检验理论判断组合导航系统是否发生缓变故障，①原假设H₀：无故障，

s_avg<Th；②备择假设H₁：存在故障，/>

s_avg≥Th，其中,λ是自由度为m的卡方分布非中心参数，m表示可视星的数量。

步骤2-4包括：计算用于组合导航系统缓变故障结束时刻判定的基本元素，其定义为外推法和残差卡方检验法故障检验统计量的比率：

由于/>

因此r_A/R总是存在且有意义的。另外，窗口长度μ的选择应考虑所使用的样本分位数类型和外推法中的滑窗长度。过大的窗口可能会导致计算窗口内的样本分位数时对异常值发生漏检，窗口过小则丧失了统计样本分位数的意义。同时，阈值的选取也与外推法中的滑窗长度有关。

步骤2-5包括：

步骤2-5-1，将r_A/R序列转换为以固定时间窗口μ的多个时间序列片段X^(μ)；

步骤2-5-2，对每个时间序列片段进行样本分位数提取，组成的数据序列X^(med)；

步骤2-5-3，将窗口内的r_A/R数据与提取的样本分位数分别作差，得到序列△X^(μ)；

步骤2-5-4，将△X^(μ)中的差值逐个与设定阈值ξ比较，若超出阈值，则判定该数据为异常值，此时对应的历元即为缓变故障结束时刻T_e，反之则不是。

本发明的有益效果为：

本发明提供了一种用于GNSS/INS组合导航卫星缓变故障检测的方法。该方法采用标准t分布和IGG-III方案设计自适应增益矩阵以缓解卡尔曼滤波故障跟踪的影响，从而缩短组合导航系统对缓变故障的检测延迟时间。由于卡尔曼滤波更新过程中可以直接获得新息，因此基于新息计算等价权矩阵相比传统抗差阵的构造更为直接，简化了计算过程。另外，t分布统计量的设计充分利用多余观测信息并引入自由度指标适时调整临界值，因此系统的抗差性能得以提升。同时，结合传统外推法的故障检测状态，设计由外推法和残差卡方检验故障检测统计量构成的r_A/R统计量，在外推法检测到缓变故障的状态下利用样本分位数原理对r_A/R序列进行异常值检测，从而判断缓变故障结束时刻，及时避免虚警现象的发生，提高系统定位精度和可靠性。本发明设计的组合导航卫星缓变故障检测方法流程简单，思路清晰。

附图说明

图1是本发明为具体实施方式提供的REKF-t算法执行流程，

图2为提供的缓变故障检测方法流程图。

具体实施方式

下面结合附图和具体实施方式，进一步阐明本发明，应理解下述具体实施方式仅用于说明本发明而不用于限制本发明的范围。

在全球导航卫星系统(GNSS)/惯性导航系统(INS)组合导航系统中，卫星伪距极可能因时钟漂移、卫星轨道建模和电离干扰等因素而发生缓变故障，严重影响系统的可靠性和稳定。因此，研究如何及时、有效地检测缓变故障这一问题具有十分重要的现实意义。

自主完好性检测外推法利用卡尔曼滤波在外推过程中的新息序列构造故障检测统计量，然后基于假设检验理论判断卫星是否发生缓变故障。然而，卡尔曼滤波器的故障跟踪作用会降低外推法对缓变故障检测的敏感性，导致缓变故障的检测延时。同时，由于传统外推法的缓变故障检测统计量依赖于卡尔曼滤波在外推过程中的新息序列，因此当缓变故障结束后系统会发生一段时间的虚警现象，以致定位结果精度下降。为了解决经典缓变故障检测法—自主完好性监测外推法在GNSS/INS组合导航中检测故障延迟时间较长和不能准确判定故障结束时刻的问题，本发明设计了一种基于抗差估计和改进外推法的GNSS/INS组合导航卫星缓变故障检测方法，其中抗差估计的使用目的是提高组合导航系统对缓变故障检测的敏感性，改进外推法则是在传统外推法的基础上通过加入异常值检测的环节以达到对缓变故障结束时刻的判定。

如图2所示，本发明主要分为两个部分，第一部分为对卫星接收机所接收的卫星数据和惯性测量单元IMU的测量数据进行抗差扩展卡尔曼滤波该部分在图中用REKF-t估计模块表示。REKF-t算法的具体执行流程可参考图1。

第1步，将惯性测量单元IMU输出的加速度数据以及角速度数据进行力学编排，得出用户在东、北、天三个方向上的位置，速度和姿态解；

第2步，将惯性导航系统INS推导出的伪距和接收机接收到的伪距数据进行融合，经过抗差扩展卡尔曼滤波(REKF)，输出用户在东、北、天三个方向上的位置误差，速度误差和姿态误差。由于组合导航中卫星伪距发生缓变故障时，传统卡尔曼滤波的故障跟踪现象会导致作为缓变故障检测重要元素的新息无法如实地跟踪到故障幅值的变化情况，从而延长外推法对缓变故障的检测延迟时间。而REKF能够通过自适应增益阵减弱发生缓变故障的观测值对状态估计的影响，这有利于对缓变故障的及时检测。本发明利用当前历元下所有卫星观测值对应新息构造服从标准t分布统计量，然后基于IGG-III方案自适应卡尔曼增益阵。由于卡尔曼滤波更新过程中可以直接获得新息，因此基于新息计算等价权矩阵相比传统抗差阵的构造更为直接，简化了计算过程。另外，t分布统计量的设计充分利用多余观测信息并引入自由度指标适时调整临界值，系统的抗差性能得以提升。具体来说，新设计的REKF-t算法相较于传统的扩展卡尔曼滤波最大的区别就是卡尔曼增益的设计，REKF-t中的抗差增益阵是在原始卡尔曼滤波的基础上增加了等价权阵的部分，其中等价权阵的设计基于标准t分布和IGG-III方案。IGG-III方案的原则是数据的有效性决定了相应权值的大小，包括零权、降权和等权三种加权方式。数据的有效性主要是通过标准化处理，然后进行等效权重函数映射的方式衡量。由于新息包含了量测的全部信息，并且对异常的观测值极为敏感，因此利用同一历元下所有卫星观测值对应新息构造标准化统计量作为IGG-III函数输入提升系统的抗差性能。标准t分布统计量的构造主要是基于同一历元下各卫星观测值对应新息间相互独立且均服从高斯分布。将每颗可视星对应的标准t分布统计量作为IGG-III方案的输入，此时根据IGG-III方案的加权原则，可将数据分为3类，即完全信任类(权重为1)，可疑类(权重介于0和1之间)，完全不信任类(权重为0)。针对发生缓变故障的卫视数据，通过这种加权处理的方式可以保证卫星和惯导组合系统输出的误差状态解基本不受异常观测值干扰，从而保证作为缓变故障检测重要元素的新息能够如实地跟踪到故障幅值的变化情况，有利于对缓变故障的及时检测；

第3步，利用REKF-t得到的误差状态修正步骤1中的状态解，同时将误差状态反馈到INS中；

第4步，时间更新，回到第1步。

第二部分是针对GNSS/INS融合系统，提取新息和其协方差阵以进行卫星缓变故障检测。

第1步，设定误警率P_fa和外推法中滑窗长度L，并将统计值r_A/R置1。误警率和滑窗长度的设定需要结合组合系统的实际场景需求，例如在飞机的航行过程中，一般取误警率为10的-6量级；滑窗长度一般不超过3分钟。统计值r_A/R的状态则反映此时卫星系统是否发生缓变故障，值为1表示无缓变故障，非1则发生缓变故障。具体r_A/R非1时的值根据具体的公式计算，将在后面的步骤中说明。开始时，默认统计值r_A/R为1；

第2步，利用满足一定高度角范围内的可视星数(自由度)和卡方分布的密度函数，求取第1步中设定误警率的缓变故障判决门限Th。然后，利用以多个历元(历元k+1-L到当前历元k)新息序列的加权累加以提高对缓变故障的检测能力的外推法计算缓变故障检测统计值；

第3步，判断卫星是否发生缓变故障并输出。当第2步中的检测统计值大于门限，则说明卫星中发生缓变故障，否则卫星系统无故障。若卫星发生缓变故障，则执行下一步以判断缓变故障结束时刻，否则返回上一步；

第4步，选择样本分位数，阈值ξ和窗口长度μ，计算系统的统计量r_A/R。其中样本分位数的取值需要考虑数据的分布特征，阈值的选取与外推法中的滑窗长度有关，窗口长度的选择应考虑所使用的样本分位数类型和外推法中滑窗长度的大小。过大的窗口可能会导致计算窗口内的样本分位数时对异常值发生漏检，窗口过小则丧失了统计样本分位数的意义。另外，将外推法和残差卡方检测统计量的比率定义为新的统计量r_A/R以表征两者的近似线性关系，用作缓变故障结束判定的基本元素。

第5步，对r_A/R序列进行异常值检测，判定缓变故障是否结束。若判断出异常值，即确定缓变故障结束时刻T_e，则执行下一步，否则返回第2步。结合缓变故障的特点，即随着时间的推移，故障的幅值越来越大，其变化为一个随着时间累积的慢变增加过程。根据r_A/R的定义，在故障的缓慢变化过程中，r_A/R值应在一个固定的范围内波动。然而，当缓变故障结束时，由于外推法的累积效应，此时的故障检测统计量通常无法快速地降低到门限以下，不会发生显著变化。同时，由于残差卡方检验统计量的计算仅依赖于当前历元的观测信息，所以该统计值在此时会迅速减小到系统未发生故障时的范围。残差卡方统计值的变化导致r_A/R值发生较大抖动，破坏了r_A/R序列相对稳定的状态。由以上分析可知，r_A/R序列出现异常值的时刻就对应缓变故障的结束时刻。下面举例说明基于样本分位数判断缓变故障结束时刻的具体过程，取样本分位数50％，阈值

μ＝3s，缓变故障开始时刻记为T_s。首先，由第4步生成数据样本X；然后将r_A/R序列X转换为以固定时间窗口μ＝3s的多个时间序列片段X^(μ)：

接着，对每个时间序列片段进行样本分位数提取,即X^(μ)中每个片段的中位数。由样本分位数组成的数据序列X^(med)：

然后，将窗口内的r_A/R数据与提取的样本分位数分别作差，得到序列△X^(μ)：

最后，将△X^(μ)中的差值逐个与设定阈值ξ比较，若超出阈值，则判定该数据为异常值，此时对应的历元即为缓变故障结束时刻T_e，反之则不是；

第6步，将统计值r_A/R置1，基于外推法计算缓变故障检测统计值置0，使用残差卡方检验法作为缓变故障的检测方法；

第7步，判断此时历元是否处于区间[T_e,T_e+L-1]中，若是，则返回上一步，否则返回第2步。

本发明提供了一种用于GNSS/INS组合导航卫星缓变故障检测的方法，一方面能够快速觉察到卫星缓变故障的存在，另一方面可以准确判定缓变故障结束时刻。应当指出，针对本技术领域的普通技术人员来说，以本发明的发明原理为基础，依然可以作出进一步的优化和改进工作，这些优化和改进仍是本发明的保护范围。特别指出的是，具体实施方式中未进一步明确论述的技术部分都可以用现有的技术实现。

Claims (2)

1.一种用于GNSS/INS组合导航卫星缓变故障检测方法，其特征在于：包括如下步骤：

步骤1，对卫星接收机所接收的卫星数据和惯性测量单元IMU的测量数据进行抗差扩展卡尔曼滤波；

步骤1-1，将惯性测量单元IMU输出的加速度数据以及角速度数据进行力学编排，得出用户在东、北、天三个方向上的位置，速度和姿态解；

步骤1-2，将惯性导航系统INS推导出的伪距和接收机接收到的伪距数据进行融合，经过抗差扩展卡尔曼滤波，输出用户在东、北、天三个方向上的位置误差，速度误差和姿态误差；

将GNSS/INS组合导航中的INS姿态、速度、位置在东、北、天三个方向上的误差和陀螺、加速度计的在相应方向上的误差以及接收机钟差等效的距离误差和与接收机频率误差等效的距离率误差作为状态量x，将惯性导航系统INS推导出的伪距和接收机接收到的伪距数据作差作为观测量z，抗差扩展卡尔曼滤波包括两个方程：状态方程x_k＝Φ_k/k-1x_k-1+w_k-1和量测方程z_k＝h(x_k)+v_k，其中，x_k是状态向量；Φ_k/k-1表示从历元k-1到k的状态转移矩阵；z_k为观测向量；h(·)表示状态x_k和观测z_k之间的非线性关系；w_k是历元k-1的过程噪声向量；v_k是历元k的观测噪声向量，其中w_k，v_k相互独立，且两者均服从零均值的高斯分布，协方差矩阵分别用Q_k和R_k表示；抗差扩展卡尔曼滤波步骤如下：

步骤1-2-1，初始化，初始化组合导航的状态x₀和误差协方差矩阵P₀；

步骤1-2-2，预测：

其中，符号“^”表示状态的估计；下标k/k-1表示从历元k-1转移到k；状态估计的误差协方差矩阵用P_k表示；

步骤1-2-3，更新：

其中，γ_k表示新息，是计算故障检测统计量的基本元素；

是抗差增益矩阵；/>

是等价权重矩阵；H_k表示非线性矢量函数h(·)的Jacobi矩阵；

另外，新息γ_k的协方差矩阵用符号S_k表示，

步骤1-3，利用步骤1-2得到的误差状态修正步骤1中的状态解，同时将误差状态反馈到INS中；

步骤1-4，时间更新，回到步骤1-1；

步骤2，针对GNSS/INS融合系统，提取新息和其协方差阵以进行卫星缓变故障检测；

步骤2-1，设定误警率P_fa和外推法中滑窗长度L，并将统计值r_A/R置1；

步骤2-2，确定故障判决门限Th，基于外推法计算缓变故障检测统计值；

步骤2-2包括：确定故障判决门限Th，

式中，P_fa表示故障检测的误警概率，F(x|m)表示自由度为m的中心卡方分布累积分布函数，另外基于外推法的缓变故障检测统计值的计算：

其中，

表示平均新息γ_avg协方差阵的逆，L表示新息序列的滑窗长度；

步骤2-3，判断卫星是否发生缓变故障并输出：若卫星发生缓变故障，则执行步骤2-4，否则返回步骤2-2；

利用二元假设检验理论判断组合导航系统是否发生缓变故障，①原假设H₀：无故障，

s_avg<Th；②备择假设H₁：存在故障，/>

s_avg≥Th，其中,λ是自由度为m的卡方分布非中心参数，m表示可视星的数量；

步骤2-4，选择样本分位数，阈值ξ和窗口长度μ，计算系统的统计量r_A/R；

计算用于组合导航系统缓变故障结束时刻判定的基本元素，其定义为外推法和残差卡方检验法故障检验统计量的比率：

由于/>

因此r_A/R总是存在且有意义的；另外，窗口长度μ的选择应考虑所使用的样本分位数类型和外推法中的滑窗长度；过大的窗口会导致计算窗口内的样本分位数时对异常值发生漏检，窗口过小则丧失了统计样本分位数的意义；同时，阈值的选取也与外推法中的滑窗长度有关；

步骤2-5，对r_A/R序列进行异常值检测，判定缓变故障是否结束；若判断出异常值，即确定缓变故障结束时刻T_e，则执行步骤2-6，否则返回步骤2-2；

步骤2-5-1，将r_A/R序列转换为以固定时间窗口μ的多个时间序列片段X^(μ)；

步骤2-5-2，对每个时间序列片段进行样本分位数提取，组成的数据序列X^(med)；

步骤2-5-3，将窗口内的r_A/R数据与提取的样本分位数分别作差，得到序列△X^(μ)；

步骤2-5-4，将△X^(μ)中的差值逐个与设定阈值ξ比较，若超出阈值，则判定该数据为异常值，此时对应的历元即为缓变故障结束时刻T_e，反之则不是；

步骤2-6，将统计值r_A/R置1，基于外推法计算缓变故障检测统计值置0，使用残差卡方检验法作为缓变故障的检测方法；

步骤2-7，判断此时历元是否处于区间[T_e,T_e+L-1]中，若是，则返回步骤2-6，否则返回步骤2-2。

2.根据权利要求1所述的一种用于GNSS/INS组合导航卫星缓变故障检测方法，其特征在于，步骤1-2-3包括：在原始卡尔曼增益K_k的基础上，利用观测数据构造等价权重矩阵

抗差增益阵与原始增益阵的关系：/>

其中/>

为对角阵；选择等价权函数IGG-III计算该矩阵；依照IGG-III方案的原理，设计了一种基于标准t分布的等价权矩阵计算方法；抗差方法记为REKF-t；等价权矩阵的对角元素：/>

式中，临界值T₀,T₁分别取t(m-2)分布显著性水平为α₀,α₁的分位数，其中，α₀,α₁人为设定；t_i是观测值z_i对应新息的标准化t分布统计量，通过当前历元k的新息向量γ_k和其对应的协方差矩阵S_k计算：

其中，

且m≠1，为标准t分布统计量的系数；/>

为标准化新息，

γ_i是第i个观测量对应的新息，/>

表示其对应的标准差，由新息协方差矩阵S_k获得；在无故障情况下，/>

另外，由t_i的表达式可知，m≠1,2，即在历元k只要可视星数不低于3颗，就能有效使用REKF-t算法。/>

Images