CN115540907A - 一种基于面向星间差分的gps/bds/ins紧组合导航的多故障检测与剔除方法 - Google Patents

一种基于面向星间差分的gps/bds/ins紧组合导航的多故障检测与剔除方法 Download PDF

Info

Publication number
CN115540907A
CN115540907A CN202211180980.9A CN202211180980A CN115540907A CN 115540907 A CN115540907 A CN 115540907A CN 202211180980 A CN202211180980 A CN 202211180980A CN 115540907 A CN115540907 A CN 115540907A
Authority
CN
China
Prior art keywords
satellite
fault detection
pseudo
fault
value
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202211180980.9A
Other languages
English (en)
Inventor
於周晴
许睿
胡哲熠
王景琪
刘建业
熊智
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing University of Aeronautics and Astronautics
Original Assignee
Nanjing University of Aeronautics and Astronautics
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing University of Aeronautics and Astronautics filed Critical Nanjing University of Aeronautics and Astronautics
Priority to CN202211180980.9A priority Critical patent/CN115540907A/zh
Publication of CN115540907A publication Critical patent/CN115540907A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C25/00Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C25/00Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass
    • G01C25/005Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass initial alignment, calibration or starting-up of inertial devices
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42Determining position
    • G01S19/45Determining position by combining measurements of signals from the satellite radio beacon positioning system with a supplementary measurement
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42Determining position
    • G01S19/45Determining position by combining measurements of signals from the satellite radio beacon positioning system with a supplementary measurement
    • G01S19/46Determining position by combining measurements of signals from the satellite radio beacon positioning system with a supplementary measurement the supplementary measurement being of a radio-wave signal type
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42Determining position
    • G01S19/45Determining position by combining measurements of signals from the satellite radio beacon positioning system with a supplementary measurement
    • G01S19/47Determining position by combining measurements of signals from the satellite radio beacon positioning system with a supplementary measurement the supplementary measurement being an inertial measurement, e.g. tightly coupled inertial
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Abstract

本发明公开了一种基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法,包括以下步骤:步骤1,初始化求得GPS与BDS的系统间偏差;步骤2,建立基于星际差分的卫星/惯性紧组合模型;步骤3,基于新息量对包括参考星在内的伪距观测量进行第一次故障检测与剔除;步骤4,在滤波后,基于残差值进行第二次故障检测与剔除,并输出最终定位结果。本发明的方法能够处理包括参考星在内的多卫星观测量故障,最大化利用正确的卫星观测量。

Description

一种基于面向星间差分的GPS/BDS/INS紧组合导航的多故障 检测与剔除方法
技术领域
本发明属于卫星定位导航技术领域,特别涉及一种基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法。
背景技术
卫星/惯性组合导航是当前室外导航的主要定位手段,其中紧组合导航技术是一种在卫星观测量端的组合导航技术,与松组合相比,其具有更高的定位精度并能够在卫星数量少于4颗时工作。随着各国卫星导航系统的建立与完善,在轨导航卫星的数量不断增加,更多的观测量能够应用于紧组合导航中,这将大大提升组合导航的定位精度,但是也增大了卫星多观测量同时发生故障的概率。故障往往是因为卫星受多路径效应、人为干扰、欺骗、卫星导航的软硬件随机故障等影响产生,尤其是在城市环境下,高楼耸立,遮挡、多路径效应等尤为明显,这些故障通常无法预知与建模。因此,实时的多观测量故障检测与剔除对于惯性/卫星紧组合导航系统有重要意义,可以保证导航定位的准确性与可靠性。
目前应用于组合导航系统中的卫星信息故障检测算法包括接收机自主完好性监测(RAIM),接收机自主完好性外推(AIME),这些方法对单故障敏感度很高但无法有效进行多故障检测。多假设解分离(MHSS)通过对所有的卫星进行组合运算将不同的故障分离到不同的子集中,根据子集解与全集解的差值构建检验统计量,是一种有效的多故障检测方法,被应用于高级接收机完好性监测(ARAIM),但是其在需要监测的观测量或故障数量较多时,易产生计算量爆炸的情况。同时,上述方法在面对星间差分卫星/惯性紧组合导航观测量故障时,无法有效分辨参考观测量故障与非参考观测量故障。
发明内容
为了解决现有技术中的问题,本发明提供一种基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法,该方法能够处理包括参考星在内的多卫星观测量故障,最大化利用正确的卫星观测量。
为实现上述目的,本发明采用的技术方案为:
一种基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法,包括以下步骤:
步骤1,初始化求得GPS与BDS的系统间偏差;
步骤2,建立基于星际差分的卫星/惯性紧组合模型;
所述步骤2中基于星际差分的卫星/惯性紧组合模型中量测方程的建立:
步骤21,对当前历元所得的卫星伪距观测值做卫星钟差、电离层误差、对流层误差的修正,并在BDS系统的伪距观测值中减去系统间偏差,最终所得量为之后用于解算的伪距观测值;
步骤22,选择卫星高度角最高的那颗卫星作为星际差分中的参考卫星L,分别计算第i个卫星的卫星伪距星际差分值
Figure BDA0003865490890000021
与惯导预测伪距星际差分值
Figure BDA0003865490890000022
Figure BDA0003865490890000023
Figure BDA0003865490890000024
其中:
Figure BDA0003865490890000025
Figure BDA0003865490890000026
分为第i个卫星和参考卫星L的伪距值,
Figure BDA0003865490890000027
Figure BDA0003865490890000028
分别为惯导解算的第i个卫星和参考卫星L的预测伪距,其值为对应卫星与接收机之间的几何距离,则假设卫星总数为N,量测量Zk和量测矩阵Hk/k-1分别为:
Figure BDA0003865490890000029
Figure BDA00038654908900000210
其中:
Figure BDA00038654908900000211
分别为历元k时的第1个,第i个,第N-1个卫星的惯导预测伪距星际差分值,
Figure BDA00038654908900000212
分别为历元k时的第1个,第i个,第N-1个卫星的伪距星际差分值;
Figure BDA00038654908900000213
分别对应着第1个、第i个,第N-1个卫星,其定义为各个卫星与参考卫星的单位观测矢量差值,即
Figure BDA00038654908900000214
Figure BDA00038654908900000215
Figure BDA00038654908900000216
分别为第i个卫星和参考卫星L的单位观测矢量,
Figure BDA00038654908900000217
为从导航系n至地球系e的旋转矩阵;
进一步的,步骤3,基于新息量对包括参考星在内的伪距观测量进行第一次故障检测与剔除,包括以下步骤:
步骤31,计算新息量γk,新息量γk定义为量测量与预测量测量的差值,即:
Figure BDA00038654908900000218
其中:Xk/k-1=Fk/k-1Xk-1为第k个历元时刻滤波器状态一步预测值,Fk/k-1为状态转移矩阵,Xk-1为上一历元状态量,由于在闭环滤波中,之前的误差都会被不断估计和修正,所以Xk-1为0,则Xk/k-1=0,即在闭环紧组合中新息量γk等于量测量Zk
步骤32,计算基于新息的第一步故障检测量Fγ,k,Fγ,k由N-1个元素组成,即
Figure BDA0003865490890000031
Figure BDA0003865490890000032
T为转置符号,
Figure BDA0003865490890000033
Figure BDA0003865490890000034
分别为第1个,第i个和第N-1个卫星对应的故障检测量,定义为:
Figure BDA0003865490890000035
其中:
Figure BDA0003865490890000036
为第i个卫星的新息量,等于γk中的第i个元素,
Figure BDA0003865490890000037
为第i个卫星的新息量的方差值,等于新息量协方差阵∑γ,k中位于第i行第i列的对角线元素,∑γ,k计算方法为:
Figure BDA0003865490890000038
其中:Pk/k-1为状态量一步预测协方差阵,Rk是量测噪声的协方差阵;
步骤33,分析Fγ,k的分布规律,筛选记录Fγ,k中符合标准正态分布的值,记录在集合A0中,即
Figure BDA0003865490890000039
并统计符合标准正态分布的元素值的个数,记为m0,即m0=length(A0);
步骤34,筛选记录Fγ,k中符合非0均值正态分布的值,记录在集合A1中,即
Figure BDA00038654908900000310
μ为正态分布的均值,μ≠0,统计符合非0均值正态分布的元素值的个数,记为m1,即m1=length(A1);由于μ值未知,所以只有当存在不低于三个元素皆满足非0均值分布时,才能够进行判断与记录,即m1≥3或m1=0;
步骤35,当m0≥2时,即存在两个以上的新息检测量满足标准正态分布,此时认为参考星伪距观测量无故障,转步骤36;当m0<2,认为参考星伪距观测量,转步骤37;
步骤36,所有不满足标准正态分布的检验统计量所对应的非参考星伪距观测量被检测为有故障,并被排除,转步骤310;
步骤37,参考星伪距观测量故障,需要更换参考星:当m1≥3时,有不低于3颗的非参考星伪距观测值正确,转步骤38;当m1<3时,无法进行无故障非参考星伪距观测值的判断,转步骤39;
步骤38,剔除原有的参考星,在集合A1中任意选择一个卫星作为参考星,A1中其余卫星作为非参考星,计算量测量Zk、量测矩阵Hk/k-1和量测噪声的协方差阵Rk,并转步骤310;
步骤39,剔除原有参考星,在剩余卫星中选择任意一个作为参考星,记录更换参考星的次数K,当K=N-2,转步骤311,否则转步骤22;
步骤310,进行扩展卡尔曼滤波,输出定位结果,基于新息的第一步故障检测方案结束;
步骤311,以惯导递推结果作为定位结果输出,基于新息的故障检测方案结束。
所述步骤1中,初始化阶段在系统正式开始工作前获得GPS与BDS的先验系统间偏差,历元k的系统间偏差求解方法为:
δtτ,k=δtuG,k-δtuB,k (1)
其中:δtuG,k和δtuB,k分别是历元k时的由伪距单点定位解算所得的GPS接收机钟差和BDS接收机钟差,δtτ,k为第k个历元的系统间偏差,对初始化阶段所得的系统间偏差求均值,即得先验系统间偏差δtτ
Figure BDA0003865490890000041
其中:j代表历元1至k中的某一历元,δtτ,j为第j个历元时的系统间偏差。
进一步的,步骤4,在滤波后,基于残差值进行第二次故障检测与剔除,并输出最终定位结果;
步骤41,计算量测量残差sk
sk=Zk-HkXk/k (10)
其中:Xk/k为状态量后验估计,Zk和Hk为经历过第一次故障检测后的量测量与量测方程,
步骤42,构造基于残差的检验统计量λk
Figure BDA0003865490890000042
其中:
Figure BDA0003865490890000043
是sk的转置矩阵,
Figure BDA0003865490890000044
为sk中第n行元素,
Figure BDA0003865490890000045
为残差协方差阵中的第n行第n列元素,κ为第一次故障检测与剔除后量测量与量矩阵的维数,若无故障发生,
Figure BDA0003865490890000046
E()为数学期望函数,χ2(κ)就是自由度为k的中心卡方分布,λk服从自由度为κ的中心卡方分布,即λk~χ2(κ);若有故障发生,
Figure BDA0003865490890000047
λk服从非中心卡方分布,即
Figure BDA0003865490890000048
Figure BDA0003865490890000049
是非中心化参数;
步骤43,进行基于残差的故障检测:
Figure BDA0003865490890000051
其中:Fd为实时计算的检测门限,根据预设的误警率Pf=1/1000来计算:
Figure BDA0003865490890000052
其中:P(λk<FD)表示λk<FD的概率,x是任意变量;
当基于残差的第二步故障检测方案无法通过时,认为基于新息的第一步故障检测方案是无效的,此时,定位将不再依赖卫星,而是直接输出惯导递推结果;若基于残差的第二步检验统计量通过,则基于新息的第一步故障检测统计量方案有效,输出组合导航滤波结果。
与现有技术相比,本发明具有以下有益效果:
本发明通过星间差分,GPS/BDS/INS紧组合不再需要估计接收机钟差,减少了需要估计的状态量,避免了接收机钟差建模不准确的问题,提高了定位精度;
本发明中两步故障检测法分别在滤波前后进行,最大程度上降低漏检可能性,保证定位的可靠性;
本发明中第一步故障检测方法根据检测量的分布提处一种新的故障检测策略,能够检测包括参考星观测量在内的多观测量故障,最大化利用正确的卫星观测量,从而保证定位的精度与稳定。
附图说明
图1是本发明的算法结构图;
图2是阶跃误差情况下的本发明故障检测与剔除方案与未加故障检测方案的定位误差对比图;
图3是斜坡误差情况下的本发明故障检测与剔除方案与未加故障检测方案的定位误差对比图。
具体实施方式
下面结合实施例对本发明作更进一步的说明。
实施例1
如图1所示,一种基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法,包括以下步骤:
步骤1,初始化求得GPS与BDS的系统间偏差;
步骤2,建立基于星际差分的卫星/惯性紧组合模型;
所述步骤2中基于星际差分的卫星/惯性紧组合模型中量测方程的建立:
步骤21,对当前历元所得的卫星伪距观测值做卫星钟差、电离层误差、对流层误差的修正,并在BDS系统的伪距观测值中减去系统间偏差,最终所得量为之后用于解算的伪距观测值;
步骤22,选择卫星高度角最高的那颗卫星作为星际差分中的参考卫星L,分别计算第i个卫星的卫星伪距星际差分值
Figure BDA0003865490890000061
与惯导预测伪距星际差分值
Figure BDA0003865490890000062
Figure BDA0003865490890000063
Figure BDA0003865490890000064
其中:
Figure BDA0003865490890000065
Figure BDA0003865490890000066
分为第i个卫星和参考卫星L的伪距值,
Figure BDA0003865490890000067
Figure BDA0003865490890000068
分别为惯导解算的第i个卫星和参考卫星L的预测伪距,其值为对应卫星与接收机之间的几何距离,则假设卫星总数为N,量测量Zk和量测矩阵Hk/k-1分别为:
Figure BDA0003865490890000069
Figure BDA00038654908900000610
其中:
Figure BDA00038654908900000611
分别为历元k时的第1个,第i个,第N-1个卫星的惯导预测伪距星际差分值,
Figure BDA00038654908900000612
分别为历元k时的第1个,第i个,第N-1个卫星的伪距星际差分值;
Figure BDA00038654908900000613
分别对应着第1个、第i个,第N-1个卫星,其定义为各个卫星与参考卫星的单位观测矢量差值,即
Figure BDA00038654908900000614
Figure BDA00038654908900000615
Figure BDA00038654908900000616
分别为第i个卫星和参考卫星L的单位观测矢量,
Figure BDA00038654908900000617
为从导航系n至地球系e的旋转矩阵;
步骤3,基于新息量对包括参考星在内的伪距观测量进行第一次故障检测与剔除,包括以下步骤:
步骤31,计算新息量γk,新息量γk定义为量测量与预测量测量的差值,即:
Figure BDA00038654908900000618
其中:Xk/k-1=Fk/k-1Xk-1为第k个历元时刻滤波器状态一步预测值,Fk/k-1为状态转移矩阵,Xk-1为上一历元状态量,由于在闭环滤波中,之前的误差都会被不断估计和修正,所以Xk-1为0,则Xk/k-1=0,即在闭环紧组合中新息量γk等于量测量Zk
步骤32,计算基于新息的第一步故障检测量Fγ,k,Fγ,k由N-1个元素组成,即
Figure BDA0003865490890000071
Figure BDA0003865490890000072
T为转置符号,
Figure BDA0003865490890000073
Figure BDA0003865490890000074
分别为第1个,第i个和第N-1个卫星对应的故障检测量,定义为:
Figure BDA0003865490890000075
其中:
Figure BDA0003865490890000076
为第i个卫星的新息量,等于γk中的第i个元素,
Figure BDA0003865490890000077
为第i个卫星的新息量的方差值,等于新息量协方差阵∑γ,k中位于第i行第i列的对角线元素,∑γ,k计算方法为:
Figure BDA0003865490890000078
其中:Pk/k-1为状态量一步预测协方差阵,Rk是量测噪声的协方差阵;
步骤33,分析Fγ,k的分布规律,筛选记录Fγ,k中符合标准正态分布的值,记录在集合A0中,即
Figure BDA0003865490890000079
并统计符合标准正态分布的元素值的个数,记为m0,即m0=length(A0);
步骤34,筛选记录Fγ,k中符合非0均值正态分布的值,记录在集合A1中,即
Figure BDA00038654908900000710
μ为正态分布的均值,μ≠0,统计符合非0均值正态分布的元素值的个数,记为m1,即m1=length(A1);由于μ值未知,所以只有当存在不低于三个元素皆满足非0均值分布时,才能够进行判断与记录,即m1≥3或m1=0;
步骤35,当m0≥2时,即存在两个以上的新息检测量满足标准正态分布,此时认为参考星伪距观测量无故障,转步骤36;当m0<2,认为参考星伪距观测量,转步骤37;
步骤36,所有不满足标准正态分布的检验统计量所对应的非参考星伪距观测量被检测为有故障,并被排除,转步骤310;
步骤37,参考星伪距观测量故障,需要更换参考星:当m1≥3时,有不低于3颗的非参考星伪距观测值正确,转步骤38;当m1<3时,无法进行无故障非参考星伪距观测值的判断,转步骤39;
步骤38,剔除原有的参考星,在集合A1中任意选择一个卫星作为参考星,A1中其余卫星作为非参考星,计算量测量Zk、量测矩阵Hk/k-1和量测噪声的协方差阵Rk,并转步骤310;
步骤39,剔除原有参考星,在剩余卫星中选择任意一个作为参考星,记录更换参考星的次数K,当K=N-2,转步骤311,否则转步骤22;
步骤310,进行扩展卡尔曼滤波,输出定位结果,基于新息的第一步故障检测方案结束;
步骤311,以惯导递推结果作为定位结果输出,基于新息的故障检测方案结束。
所述步骤1中,初始化阶段在系统正式开始工作前获得GPS与BDS的先验系统间偏差,历元k的系统间偏差求解方法为:
δtτ,k=δtuG,k-δtuB,k (1)
其中:δtuG,k和δtuB,k分别是历元k时的由伪距单点定位解算所得的GPS接收机钟差和BDS接收机钟差,δtτ,k为第k个历元的系统间偏差,对初始化阶段所得的系统间偏差求均值,即得先验系统间偏差δtτ
Figure BDA0003865490890000081
其中:j代表历元1至k中的某一历元,δtτ,j为第j个历元时的系统间偏差。
如图2所示,步骤4,在滤波后,基于残差值进行第二次故障检测与剔除,并输出最终定位结果;
步骤41,计算量测量残差sk
sk=Zk-HkXk/k (10)
其中:Xk/k为状态量后验估计,Zk和Hk为经历过第一次故障检测后的量测量与量测方程,
步骤42,构造基于残差的检验统计量λk
Figure BDA0003865490890000082
其中:
Figure BDA0003865490890000083
是sk的转置矩阵,
Figure BDA0003865490890000084
为sk中第n行元素,
Figure BDA0003865490890000085
为残差协方差阵中的第n行第n列元素,κ为第一次故障检测与剔除后量测量与量矩阵的维数,若无故障发生,
Figure BDA0003865490890000086
E()为数学期望函数,χ2(κ)就是自由度为k的中心卡方分布,λk服从自由度为κ的中心卡方分布,即λk~χ2(κ);若有故障发生,
Figure BDA0003865490890000087
λk服从非中心卡方分布,即
Figure BDA0003865490890000088
Figure BDA0003865490890000089
是非中心化参数;
步骤43,进行基于残差的故障检测:
Figure BDA0003865490890000091
其中:Fd为实时计算的检测门限,根据预设的误警率Pf=1/1000来计算:
Figure BDA0003865490890000092
其中:P(λk<FD)表示λk<FD的概率,x是任意变量;
当基于残差的第二步故障检测方案无法通过时,认为基于新息的第一步故障检测方案是无效的,此时,定位将不再依赖卫星,而是直接输出惯导递推结果;若基于残差的第二步检验统计量通过,则基于新息的第一步故障检测统计量方案有效,输出组合导航滤波结果。
针对故障检测部分:通过新息和残差分别构建检验统计量,首先根据基于新息的检验统计量的分布规律识别正确的非参考星观测量,并检测参考星是否有故障,若有故障则在正确的非参考星中挑选新的参考星进行替代并重新开始;若无参考星故障则剔除故障的非参考星观测量并进行滤波。在滤波完成后进行第二步故障检测,根据基于残差的检验统计量检测所使用卫星观测量,若通过,则证明第一步基于新息检测的方案有效,输出滤波定位结果;若不通过,则说明第一步基于新息检测的方案无效,直接输出惯导递推结果。
实施例2
仿真说明
故障仿真考虑如下因素:故障类型(阶跃故障、斜坡故障)、发生故障的卫星观测量个数,参考星是否有故障,故障的大小,故障持续时间。本实施例2中分别设置了不同类型的阶跃故障和斜坡故障仿真,具体仿真情况分别如表1和表2所示。
阶跃故障
阶跃故障类型详细情况如表1所示,各故障情况主要在于故障大小、故障数量、是否有参考星观测量故障。其中G23卫星是本次实验的参考星,即情况2、4含有参考星故障。
表1.阶跃故障类型仿真情况
Figure BDA0003865490890000093
参考图2,不同情况的阶跃故障被添加在观测伪距上,本发明所得的定位误差与未加故障检测方案的定位误差如图2所示,本发明所述方法能够有效识别并剔除包括参考星在内的多卫星观测量阶跃故障,可识别并剔除低至5m的阶跃故障,将定位误差维持在正常标准,保证定位的完好性与可靠性。
斜坡故障
斜坡故障类型详细情况如表2所示,各故障情况主要在于故障斜率、故障持续时间、故障数量、是否有参考星观测量故障。其中G23卫星是本次实验的参考星,即情况2、4含有参考星故障。
表2.斜坡故障类型仿真情况
Figure BDA0003865490890000101
参考图3,不同形式的阶跃故障被添加在观测伪距上,本发明所得的定位误差与未加故障检测方案的定位误差如附图3所示。未使用故障检测和剔除时,斜坡故障带来的定位误差是巨大的且在故障消失后无法及时恢复原有的定位精度。本发明所述方法能够有效识别并剔除包括参考星在内的多卫星观测量斜坡故障,定位精度明显提高,有效保证定位的完好性与可靠性。
以上所述仅是本发明的优选实施方式,应当指出:对于本技术领域的普通技术人员来说,在不脱离本发明原理的前提下,还可以做出若干改进和润饰,这些改进和润饰也应视为本发明的保护范围。

Claims (3)

1.一种基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法,其特征在于,包括以下步骤:
步骤1,初始化求得GPS与BDS的系统间偏差;
步骤2,建立基于星际差分的卫星/惯性紧组合模型;
步骤3,基于新息量对包括参考星在内的伪距观测量进行第一次故障检测与剔除;
步骤4,在滤波后,基于残差值进行第二次故障检测与剔除,并输出最终定位结果;
所述步骤2中基于星际差分的卫星/惯性紧组合模型中量测方程的建立:
步骤21,对当前历元所得的卫星伪距观测值做卫星钟差、电离层误差、对流层误差的修正,并在BDS系统的伪距观测值中减去系统间偏差,最终所得量为之后用于解算的伪距观测值;
步骤22,选择卫星高度角最高的那颗卫星作为星际差分中的参考卫星L,分别计算第i个卫星的卫星伪距星际差分值
Figure FDA0003865490880000011
与惯导预测伪距星际差分值
Figure FDA0003865490880000012
Figure FDA0003865490880000013
Figure FDA0003865490880000014
其中:
Figure FDA0003865490880000015
Figure FDA0003865490880000016
分为第i个卫星和参考卫星L的伪距值,
Figure FDA0003865490880000017
Figure FDA0003865490880000018
分别为惯导解算的第i个卫星和参考卫星L的预测伪距,其值为对应卫星与接收机之间的几何距离,则假设卫星总数为N,量测量Zk和量测矩阵Hk/k-1分别为:
Figure FDA0003865490880000019
Figure FDA00038654908800000110
其中:
Figure FDA00038654908800000111
分别为历元k时的第1个,第i个,第N-1个卫星的惯导预测伪距星际差分值,
Figure FDA00038654908800000112
分别为历元k时的第1个,第i个,第N-1个卫星的伪距星际差分值;
Figure FDA00038654908800000113
分别对应着第1个、第i个,第N-1个卫星,其定义为各个卫星与参考卫星的单位观测矢量差值,即
Figure FDA00038654908800000114
Figure FDA00038654908800000115
Figure FDA00038654908800000116
分别为第i个卫星和参考卫星L的单位观测矢量,
Figure FDA00038654908800000117
为从导航系n至地球系e的旋转矩阵;
所述步骤3包括以下步骤:
步骤31,计算新息量Υk,新息量Υk定义为量测量与预测量测量的差值,即:
Figure FDA0003865490880000021
其中:xk/k-1=Fk/k-1xk-1为第k个历元时刻滤波器状态一步预测值,Fk/k-1为状态转移矩阵,xk-1为上一历元状态量,由于在闭环滤波中,之前的误差都会被不断估计和修正,所以xk-1为0,则xk/k-1=0,即在闭环紧组合中新息量Υk等于量测量zk
步骤32,计算基于新息的第一步故障检测量FΥ,k,FΥ,k由N-1个元素组成,即
Figure FDA0003865490880000022
Figure FDA0003865490880000023
T为转置符号,
Figure FDA0003865490880000024
Figure FDA0003865490880000025
分别为第1个,第i个和第N-1个卫星对应的故障检测量,定义为:
Figure FDA0003865490880000026
其中:
Figure FDA0003865490880000027
为第i个卫星的新息量,等于Υk中的第i个元素,
Figure FDA0003865490880000028
为第i个卫星的新息量的方差值,等于新息量协方差阵∑Υ,k中位于第i行第i列的对角线元素,∑Υ,k计算方法为:
Figure FDA0003865490880000029
其中:Pk/k-1为状态量一步预测协方差阵,Rk是量测噪声的协方差阵;
步骤33,分析FΥ,k的分布规律,筛选记录FΥ,k中符合标准正态分布的值,记录在集合A0中,即
Figure FDA00038654908800000210
并统计符合标准正态分布的元素值的个数,记为m0,即m0=length(A0);
步骤34,筛选记录FΥ,k中符合非0均值正态分布的值,记录在集合A1中,即
Figure FDA00038654908800000211
μ为正态分布的均值,μ≠0,统计符合非0均值正态分布的元素值的个数,记为m1,即m1=length(A1);由于μ值未知,所以只有当存在不低于三个元素皆满足非0均值分布时,才能够进行判断与记录,即m1≥3或m1=0;
步骤35,当m0≥2时,即存在两个以上的新息检测量满足标准正态分布,此时认为参考星伪距观测量无故障,转步骤36;当m0<2,认为参考星伪距观测量,转步骤37;
步骤36,所有不满足标准正态分布的检验统计量所对应的非参考星伪距观测量被检测为有故障,并被排除,转步骤310;
步骤37,参考星伪距观测量故障,需要更换参考星:当m1≥3时,有不低于3颗的非参考星伪距观测值正确,转步骤38;当m1<3时,无法进行无故障非参考星伪距观测值的判断,转步骤39;
步骤38,剔除原有的参考星,在集合A1中任意选择一个卫星作为参考星,A1中其余卫星作为非参考星,计算量测量Zk、量测矩阵Hk/k-1和量测噪声的协方差阵Rk,并转步骤310;
步骤39,剔除原有参考星,在剩余卫星中选择任意一个作为参考星,记录更换参考星的次数K,当K=N-2,转步骤311,否则转步骤22;
步骤310,进行扩展卡尔曼滤波,输出定位结果,基于新息的第一步故障检测方案结束;
步骤311,以惯导递推结果作为定位结果输出,基于新息的故障检测方案结束。
2.根据权利要求1所述的基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法,其特征在于,
所述步骤1中,初始化阶段在系统正式开始工作前获得GPS与BDS的先验系统间偏差,历元k的系统间偏差求解方法为:
δtτ,k=δtuG,k-δtuB,k (1)
其中:δtuG,k和δtuB,k分别是历元k时的由伪距单点定位解算所得的GPS接收机钟差和BDS接收机钟差,δtτ,k为第k个历元的系统间偏差,对初始化阶段所得的系统间偏差求均值,即得先验系统间偏差δtτ
Figure FDA0003865490880000031
其中:j代表历元1至k中的某一历元,δtτ,j为第j个历元时的系统间偏差。
3.根据权利要求2所述的基于面向星间差分的GPS/BDS/INS紧组合导航的多故障检测与剔除方法,其特征在于,
所述步骤4:在滤波后,基于残差值进行第二次故障检测与剔除,并输出最终定位结果,
步骤41,计算量测量残差sk
sk=Zk-HkXk/k (10)
其中:Xk/k为状态量后验估计,Zk和Hk为经历过第一次故障检测后的量测量与量测方程,
步骤42,构造基于残差的检验统计量λk
Figure FDA0003865490880000041
其中:
Figure FDA0003865490880000042
是sk的转置矩阵,
Figure FDA0003865490880000043
为sk中第n行元素,
Figure FDA0003865490880000044
为残差协方差阵中的第n行第n列元素,κ为第一次故障检测与剔除后量测量与量矩阵的维数,若无故障发生,
Figure FDA0003865490880000045
E(·)为数学期望函数,χ2(κ)是自由度为κ的中心卡方分布,λk服从自由度为κ的中心卡方分布,即λk~χ2(κ);若有故障发生,
Figure FDA0003865490880000046
λk服从非中心卡方分布,即
Figure FDA0003865490880000047
Figure FDA0003865490880000048
是非中心化参数;
步骤43,进行基于残差的故障检测:
Figure FDA0003865490880000049
其中:Fd为实时计算的检测门限,根据预设的误警率Pf=1/1000来计算:
Figure FDA00038654908800000410
其中:P(λk<FD)表示λk<FD的概率,x是任意变量;
当基于残差的第二步故障检测方案无法通过时,认为基于新息的第一步故障检测方案是无效的,此时,定位将不再依赖卫星,而是直接输出惯导递推结果;若基于残差的第二步检验统计量通过,则基于新息的第一步故障检测统计量方案有效,输出组合导航滤波结果。
CN202211180980.9A 2022-09-27 2022-09-27 一种基于面向星间差分的gps/bds/ins紧组合导航的多故障检测与剔除方法 Pending CN115540907A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211180980.9A CN115540907A (zh) 2022-09-27 2022-09-27 一种基于面向星间差分的gps/bds/ins紧组合导航的多故障检测与剔除方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211180980.9A CN115540907A (zh) 2022-09-27 2022-09-27 一种基于面向星间差分的gps/bds/ins紧组合导航的多故障检测与剔除方法

Publications (1)

Publication Number Publication Date
CN115540907A true CN115540907A (zh) 2022-12-30

Family

ID=84729016

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211180980.9A Pending CN115540907A (zh) 2022-09-27 2022-09-27 一种基于面向星间差分的gps/bds/ins紧组合导航的多故障检测与剔除方法

Country Status (1)

Country Link
CN (1) CN115540907A (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117111101A (zh) * 2023-06-26 2023-11-24 北京航空航天大学 消除双层空基导航增强自组网杠杆效应的故障检测方法

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117111101A (zh) * 2023-06-26 2023-11-24 北京航空航天大学 消除双层空基导航增强自组网杠杆效应的故障检测方法
CN117111101B (zh) * 2023-06-26 2024-03-22 北京航空航天大学 消除双层空基导航增强自组网杠杆效应的故障检测方法

Similar Documents

Publication Publication Date Title
US7720642B2 (en) Measurement fault detection
US8106823B2 (en) Method of operating a satellite navigation receiver
RU2559842C2 (ru) Способ и устройство для обнаружения и исключения множественных отказов спутников системы гнсс
CN109359270B (zh) 北斗地基增强系统完好性风险监测的阈值模型建立方法
US9146320B2 (en) Method for detecting and excluding multiple failures in a satellite
CN109100748B (zh) 一种基于低轨星座的导航完好性监测系统及方法
JP5572877B2 (ja) 異常な擬似距離測定値から無線ナビゲーション受信機ユーザを保護するための方法
US11409002B2 (en) Method for operating a plurality of GNSS receivers for detecting satellite signal deformation
CN111060133B (zh) 一种用于城市复杂环境的组合导航完好性监测方法
CN110068840B (zh) 一种基于伪距测量特征值提取的araim故障检测方法
CN111965668A (zh) 一种面向卫星多故障的raim方法
CN115420284B (zh) 一种组合导航系统故障检测与识别方法
CN116075747A (zh) 用于卫星定位的系统和方法
CN110879407A (zh) 一种基于完好性风险模型的卫星导航观测量新息检测方法
CN106932793A (zh) 一种北斗三频信号的实时周跳探测与修复方法
CN112130177A (zh) 一种基于稳定分布的地基增强系统完好性监测方法
CN115540907A (zh) 一种基于面向星间差分的gps/bds/ins紧组合导航的多故障检测与剔除方法
CN113009520A (zh) 一种卫星导航矢量跟踪环的完好性检测方法
CN109946722B (zh) 一种多系统多频段定位方法及系统
CN113848570A (zh) 基于非高斯分布的自适应实时多径消除及抗差定位方法
Patino-Studencka et al. Approach for detection and identification of multiple faults in satellite navigation
Almagbile et al. Analysis of outlier separability in integrated GPS/INS systems
CN113777629B (zh) 一种地基增强系统星钟故障组合监测方法
JP3469939B2 (ja) 衛星測位装置及びその衛星の故障回復判定方法
Tanil et al. Kalman filter partial innovation sequence monitor

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination