CN111999750A - 针对杆臂不准的实时单站周跳探测改进方法 - Google Patents

针对杆臂不准的实时单站周跳探测改进方法 Download PDF

Info

Publication number
CN111999750A
CN111999750A CN202011175778.8A CN202011175778A CN111999750A CN 111999750 A CN111999750 A CN 111999750A CN 202011175778 A CN202011175778 A CN 202011175778A CN 111999750 A CN111999750 A CN 111999750A
Authority
CN
China
Prior art keywords
cycle slip
slip detection
epoch
satellite
test statistic
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.)
Granted
Application number
CN202011175778.8A
Other languages
English (en)
Other versions
CN111999750B (zh
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.)
National University of Defense Technology
Original Assignee
National University of Defense Technology
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 National University of Defense Technology filed Critical National University of Defense Technology
Priority to CN202011175778.8A priority Critical patent/CN111999750B/zh
Publication of CN111999750A publication Critical patent/CN111999750A/zh
Application granted granted Critical
Publication of CN111999750B publication Critical patent/CN111999750B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13—Receivers
    • G01S19/24—Acquisition or tracking or demodulation of signals transmitted by the system
    • G01S19/243—Demodulation of navigation message
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13—Receivers
    • G01S19/24—Acquisition or tracking or demodulation of signals transmitted by the system
    • G01S19/25—Acquisition or tracking or demodulation of signals transmitted by the system involving aiding data received from a cooperating element, e.g. assisted GPS
    • G01S19/256—Acquisition or tracking or demodulation of signals transmitted by the system involving aiding data received from a cooperating element, e.g. assisted GPS relating to timing, e.g. time of week, code phase, timing offset
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining 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/40—Correcting position, velocity or attitude
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO 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/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining 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/42—Determining position
    • G01S19/48—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system
    • G01S19/49—Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system whereby the further system is an inertial position system, e.g. loosely-coupled
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Evolutionary Biology (AREA)
  • Operations Research (AREA)
  • Probability & Statistics with Applications (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Algebra (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

针对基于伪距和多普勒紧组合的单站惯导辅助周跳探测方法在杆臂不准条件下的误警率高的缺陷，本发明提出了一种针对杆臂不准的实时单站周跳探测改进方法。首先利用ARIMA模型去除基础检验统计量趋势项，然后将GARCH模型用于估计带有波动聚集特性残差的条件方差。该方法不仅可以在较大杆臂误差条件下极大地降低误警率，而且可以准确地检测1周仿真周跳，同时又具有较低的计算量。

Description

针对杆臂不准的实时单站周跳探测改进方法
技术领域
本发明属于卫星定位中的周跳探测技术领域，特别涉及一种基于自回归综合滑动平均-广义自回归条件异方差 （Autoregressive Integrated Moving Average-Generalized Autoregressive Conditional Heteroscedastic， ARIMA-GARCH）模型的实时惯导辅助周跳探测改进方法。
背景技术
周跳是载波相位高精度定位中最常见的故障之一。周跳发生的原因通常是接收机跟踪环暂时失锁导致测量的载波相位突发整数跳变，并且该跳变在之后时间一直维持。周跳很容易使定位精度恶化，因此在将载波相位测量信息用于高精度定位之前必须探测和排除周跳故障。
杆臂是全球导航卫星系统（GNSS）天线与惯性测量器件（IMU）参考坐标系之间的相对位置。通常GNSS天线安装在车辆外部，IMU安装在车辆内部。因此，直接测量GNSS天线和IMU之间的距离常常十分困难，特别是大型车辆、船等设备。杆臂误差是GNSS/INS组合系统的主要误差源之一。
常见的周跳探测方法有码相比较法、电离层残差方法、多普勒积分法和载波相位历元间差分法（参见[1] Xu G (2007) GPS: Theory, Algorithms and Applications.Springer Berlin）。这些方法通常需要接收来自基准站播发的观测数据构造双差检验统计量，并且无法探测单站单频周跳。惯性导航系统通过提供短时间高精度定位结果可以辅助进行单频单接收机周跳探测。现有文献中有许多针对惯导辅助周跳探测的研究（参见[2]Karaim M (2013) Real-time Cycle-slip Detection and Correction for LandVehicle Navigation using Inertial Aiding. Queen’s University. [3] Du S, Gao Y(2012) Inertial Aided Cycle Slip Detection and Identification for IntegratedPPP GPS and INS. Sensors. 12:14344-14362 doi:10.3390/s121114344. [4] TakasuT, Yasuda A (2008) Cycle Slip Detection and Fixing by MEMS-IMU/GPSIntegration for Mobile Environment RTK-GPS. In: Proceedings of ION GNSS 2008,Institute of Navigation, Savannah, GA, USA, September 16-19 2008）。
然而，文献中还没有方法考虑杆臂误差对惯导辅助周跳探测影响的问题，也没有将ARIMA-GARCH模型应用到杆臂不准条件下惯导辅助探测周跳中的相关研究公布。
发明内容
针对基于伪距和多普勒紧组合的单站惯导辅助周跳探测方法在杆臂不准条件下的误警率高的缺陷，本发明提出了一种针对杆臂不准的实时单站周跳探测改进方法。
为实现上述技术目的，本发明采用的具体技术方案如下：
针对杆臂不准的实时单站周跳探测改进方法，包括以下步骤：
S1：设置参数更新周期
Figure 113798DEST_PATH_IMAGE001
和滑动窗口大小
Figure 57483DEST_PATH_IMAGE002
；
S2：计算当前历元的基础惯辅检验统计量及其方差，并保存在滑动窗口；
S3：当滑动窗口中保存的历史历元的基础惯辅检验统计量数量
Figure 554323DEST_PATH_IMAGE003
≥
Figure 408010DEST_PATH_IMAGE004
且当前建模时间
Figure 394421DEST_PATH_IMAGE005
＞
Figure 28664DEST_PATH_IMAGE006
，则计算ARIMA模型参数和GARCH模型参数，然后初始化建模时刻
Figure 1299DEST_PATH_IMAGE007
，并进入下一步；若
Figure 834126DEST_PATH_IMAGE008
≤
Figure 866804DEST_PATH_IMAGE009
，则直接进入下一步。
S4：基于滑动窗口中保存的历史惯辅检验统计量，利用ARIMA模型通过一阶时间差分对基础惯辅检验统计量去趋势项，并计算去趋势项后的残差及其方差；
S5：周跳探测；
判断前一历元去趋势项后的基础惯辅检验统计量的残差及其方差是否保存，若前一历元去趋势项后的基础惯辅检验统计量的残差及其方差没有保存，则进行第一步周跳探测，若检验不通过，则有周跳；若检验通过，则无周跳，并更新滑动窗口；周跳探测结束；
若前一历元去趋势项后的基础惯辅检验统计量的残差及其方差已保存，则基于GARCH模型计算条件方差，进入第二步周跳探测，若检验不通过，则有周跳；若检验通过，则无周跳，并更新滑动窗口，周跳探测结束。
作为优选方案，本发明中所述ARIMA模型为ARIMA (n, 1, 0)模型，ARIMA模型参数为ARIMA (n, 1, 0)模型参数。ARIMA (n, 1, 0)模型参数由最小二乘估计法、Yule-Walker方法或Burg方法计算得到。GARCH模型为GARCH (1,1)模型，GARCH模型为GARCH (1,1)模型参数。GARCH (1,1)模型参数采用极大似然法计算得到。
本发明中，所述基础惯辅检验统计量的计算方法是：
Figure 519502DEST_PATH_IMAGE010
其中
Figure 528785DEST_PATH_IMAGE011
和
Figure 216118DEST_PATH_IMAGE012
分别表示待检卫星和参考星编号，
Figure 685276DEST_PATH_IMAGE013
表示频点，
Figure 825271DEST_PATH_IMAGE014
表示波长，
Figure 733184DEST_PATH_IMAGE015
和
Figure 150390DEST_PATH_IMAGE016
表示相邻历元，
Figure 915084DEST_PATH_IMAGE017
表示星间历元间双差载波相位观测值，
Figure 745636DEST_PATH_IMAGE018
表示基于惯导提供的位置得到的双差星地距离，
Figure 863765DEST_PATH_IMAGE019
表示光速，
Figure 525690DEST_PATH_IMAGE020
表示双差卫星钟差，通过广播星历计算获得。
所述基础惯辅检验统计量方差为：
Figure 398969DEST_PATH_IMAGE021
其中
Figure 910894DEST_PATH_IMAGE022
表示历元
Figure 957347DEST_PATH_IMAGE023
和
Figure 677042DEST_PATH_IMAGE024
间位置增量协方差阵；
Figure 393325DEST_PATH_IMAGE025
和
Figure 260787DEST_PATH_IMAGE026
分别表示相邻历元
Figure 48614DEST_PATH_IMAGE027
和
Figure 560498DEST_PATH_IMAGE028
的载波相位测量误差方差。
Figure 572317DEST_PATH_IMAGE029
，
Figure 599178DEST_PATH_IMAGE030
和
Figure 862801DEST_PATH_IMAGE031
分别卫星
Figure 353825DEST_PATH_IMAGE032
与站点
Figure 910446DEST_PATH_IMAGE033
在
Figure 690183DEST_PATH_IMAGE034
和
Figure 882130DEST_PATH_IMAGE035
历元的视线矢量;
Figure 837448DEST_PATH_IMAGE036
和
Figure 659910DEST_PATH_IMAGE037
分别卫星
Figure 458102DEST_PATH_IMAGE038
与站点
Figure 63527DEST_PATH_IMAGE039
在
Figure 263564DEST_PATH_IMAGE040
和
Figure 929031DEST_PATH_IMAGE041
历元的视线矢量;
Figure 948940DEST_PATH_IMAGE042
表示
Figure 856591DEST_PATH_IMAGE043
到
Figure 911135DEST_PATH_IMAGE044
历元间状态转移矩阵的逆。
在本发明的S4中，利用ARIMA (n, 1, 0)模型计算去趋势项后的检验统计量的残差及其方差，方法如下：
ARIMA (n, 1, 0)表达式为
Figure 809821DEST_PATH_IMAGE045
其中
Figure 457971DEST_PATH_IMAGE046
和
Figure 795411DEST_PATH_IMAGE047
分别表示相邻历元的时间序列,
Figure 579827DEST_PATH_IMAGE048
表示第
Figure 914994DEST_PATH_IMAGE049
项ARIMA (n, 1, 0)模型参数，
Figure 175074DEST_PATH_IMAGE050
为白噪声。基础惯辅检验统计量的趋势项通过
Figure 925992DEST_PATH_IMAGE051
和
Figure 955128DEST_PATH_IMAGE052
的一阶时间差分消除。给定历史时间滑动窗口
Figure 195617DEST_PATH_IMAGE053
, 定义窗口内保存的卫星
Figure 316894DEST_PATH_IMAGE054
的基础检验统计量序列为
Figure 730558DEST_PATH_IMAGE055
式中
Figure 286304DEST_PATH_IMAGE056
表示
Figure 697694DEST_PATH_IMAGE057
历元卫星
Figure 542153DEST_PATH_IMAGE058
在滑动窗口中
Figure 25087DEST_PATH_IMAGE059
历元的基础检验统计量，可由在S2中计算获得。
对
Figure 638602DEST_PATH_IMAGE060
进行一阶差分，可以得到以下时间序列
Figure 486472DEST_PATH_IMAGE061
其中
Figure 942862DEST_PATH_IMAGE062
在确定ARIMA (n, 1, 0)模型参数
Figure 72229DEST_PATH_IMAGE063
之后，可获得当前历元（
Figure 930464DEST_PATH_IMAGE064
历元）与其前一历元（
Figure 621339DEST_PATH_IMAGE065
历元）的基础惯辅检验统计量的差分预报值为：
Figure 299445DEST_PATH_IMAGE066
其中，
Figure 733969DEST_PATH_IMAGE067
表示ARIMA (n, 1, 0)模型参数，
Figure 712289DEST_PATH_IMAGE068
表示
Figure 574066DEST_PATH_IMAGE069
历元基础惯辅检验统计量一阶差分值。
Figure 5047DEST_PATH_IMAGE070
为基础惯辅检验统计量一阶差分预报值，趋势项通过一阶差分去除。
从而，可得当前历元卫星
Figure 305579DEST_PATH_IMAGE071
的去趋势项后的基础惯辅检验统计量的残差为
Figure 512307DEST_PATH_IMAGE072
相应的，当前历元卫星
Figure 669619DEST_PATH_IMAGE073
的基础惯辅检验统计量的方差估计为：
Figure 525579DEST_PATH_IMAGE074
其中
Figure 301905DEST_PATH_IMAGE075
表示历史时间滑动窗口
Figure 254818DEST_PATH_IMAGE076
中的样本总数。
本发明的S5中，第一步周跳探测方法为：
第一步周跳探测中检验统计量为当前历元卫星
Figure 458397DEST_PATH_IMAGE077
的基础惯辅检验统计量的残差
Figure 801654DEST_PATH_IMAGE078
，即
Figure 506305DEST_PATH_IMAGE079
第一步周跳探测中的检验门限计算为：
Figure 923511DEST_PATH_IMAGE080
其中，
Figure 625888DEST_PATH_IMAGE081
表示设定的误警率；
Figure 518757DEST_PATH_IMAGE082
表示标准正态分布函数的逆函数。
第一步周跳探测中，周跳判断准则为：
Figure 135421DEST_PATH_IMAGE083
。
本发明的S5中，第二步周跳探测方法为：
基于GARCH模型计算条件方差，公式为：
Figure 609DEST_PATH_IMAGE084
假设当前历元是
Figure 342728DEST_PATH_IMAGE085
，式中
Figure 394998DEST_PATH_IMAGE086
表示前一历元卫星
Figure 441452DEST_PATH_IMAGE077
的去趋势项后的基础惯辅检验统计量的残差，
Figure 98829DEST_PATH_IMAGE087
表示前一历元卫星
Figure 939746DEST_PATH_IMAGE088
由GARCH模型计算得到的方差，需要通过存储获得。其中
Figure 479312DEST_PATH_IMAGE089
,
Figure 939243DEST_PATH_IMAGE090
,
Figure 575761DEST_PATH_IMAGE091
和
Figure 790841DEST_PATH_IMAGE092
均为GARCH (1,1)模型参数。
第二步周跳探测中检验统计量为当前历元卫星
Figure 253922DEST_PATH_IMAGE054
的基础惯辅检验统计量的残差
Figure 642178DEST_PATH_IMAGE093
即
Figure 8568DEST_PATH_IMAGE094
第二步周跳探测中的检验门限计算为：
Figure 394550DEST_PATH_IMAGE095
其中，
Figure 705446DEST_PATH_IMAGE096
表示设定的误警率；
Figure 772759DEST_PATH_IMAGE097
表示标准正态分布函数的逆函数。
第二步周跳探测中，周跳判断准则为：
Figure 118289DEST_PATH_IMAGE098
。
与现有技术相比，本发明具有以下优点：
通过分析杆臂不准条件下的检验统计量的误差特性，本发明采用了针对性的去趋势项和计算条件方差的处理，可以有效降低该条件下误警率，
同时可以探测1周大小的小周跳。本发明既能够提高观测信息的可用性，又能够保证可靠性。
本发明将ARIMA-GARCH模型应用于惯导辅助周跳探测，针对实时应用提出了定时参数更新的方法即通过设定设置参数更新周期
Figure 675173DEST_PATH_IMAGE099
，可以极大地降低计算量。
本发明可以探测和识别单站单频周跳，无需使用基准站播发的数据，特别是在杆臂难以准确测量的大型车辆、船等应用场景中具有特殊的优势。
附图说明
图1为实施例1 的流程图；
图2为本发明某次试验中不同方法北斗C03卫星的B1频点的误警历元个数比较图；
图3为本发明某次试验中不同方法北斗C08卫星的B1频点的误警历元个数比较图；
图4为本发明某次试验中不同方法门限与检验统计量的比较图；
图5 显示了本发明某次试验中添加
Figure 614310DEST_PATH_IMAGE100
周小周跳条件下，不同方法北斗C08卫星的B1频点的检验统计量与检验门限对比图；
图6 显示了本发明某次试验中添加
Figure 344368DEST_PATH_IMAGE100
周小周跳条件下，不同方法北斗C08卫星的B3频点的检验统计量与检验门限对比图。
具体实施方案
为使本发明实施例的目的、技术方案和优点更加清楚明白，下面将以附图及详细叙述清楚说明本发明所揭示内容的精神，任何所属技术领域技术人员在了解本发明内容的实施例后，当可由本发明内容所教示的技术，
加以改变及修饰，其并不脱离本发明内容的精神与范围。本发明的示意性实施例及其说明用于解释本发明，但并不作为对本发明的限定。
参照图1，本实施例提供一种针对杆臂不准的实时单站周跳探测改进方法，包括以下步骤：
步骤一：设置参数更新周期
Figure 747668DEST_PATH_IMAGE101
和滑动窗口大小
Figure 646092DEST_PATH_IMAGE102
。
步骤二：计算当前历元基础惯辅检验统计量及其方差。
步骤三：初始化准备用于建模的滑动窗口，判断滑动窗口保存的历史基础惯辅检验统计量数量
Figure 197159DEST_PATH_IMAGE103
是否小于
Figure 934170DEST_PATH_IMAGE104
，若
Figure 864080DEST_PATH_IMAGE105
，进入下一步；若
Figure 825083DEST_PATH_IMAGE106
＜
Figure 801129DEST_PATH_IMAGE107
，则将当前历元基础惯辅检验统计量保存进入滑动窗口后返回步骤二，计算下一历元基础惯辅检验统计量及其方差。
步骤四：判断当前建模时间
Figure 748357DEST_PATH_IMAGE108
是否大于参数更新周期
Figure 922986DEST_PATH_IMAGE109
，若
Figure 930256DEST_PATH_IMAGE110
＞
Figure 128020DEST_PATH_IMAGE111
，则分别利用最小二乘法和极大似然法计算ARIMA(1,1,0)模型参数和GARCH(1,1)模型参数，然后初始化建模时刻
Figure 3572DEST_PATH_IMAGE112
，并进入下一步；若
Figure 406609DEST_PATH_IMAGE113
≤
Figure 709414DEST_PATH_IMAGE114
，则直接进入下一步。
步骤五：判断ARIMA(1,1,0)模型参数和GARCH(1,1)模型参数是否求解成功，若求解成功，假设相邻历元时间间隔为
Figure 66577DEST_PATH_IMAGE115
，则建模时刻
Figure 214662DEST_PATH_IMAGE116
，则进入下一步；若求解不成功，则返回步骤二，计算下一历元基础惯辅检验统计量及其方差。
步骤六：基于滑动窗口中保存的历史惯辅检验统计量，利用ARIMA模型通过一阶时间差分对基础惯辅检验统计量去趋势项，并计算去趋势项后的残差及其方差。
步骤七：判断前一历元去趋势项后的基础惯辅检验统计量的残差及其方差是否保存；若前一历元去趋势项后的基础惯辅检验统计量的残差及其方差没有保存，则进行第一步周跳探测，若检验不通过，则有周跳，周跳探测结束。若检验通过，则无周跳，并更新滑动窗口，周跳探测结束。若前一历元去趋势项后的基础惯辅检验统计量的残差及其方差已保存，则进入下一步；
步骤八：基于GARCH(1,1)模型计算条件方差；
步骤九：进入第二步周跳探测，若检验不通过，则有周跳，周跳探测结束；若检验通过，则无周跳，并更新滑动窗口，周跳探测结束。
本发明中，基础惯辅检验统计量的计算方法为：
Figure 239250DEST_PATH_IMAGE117
其中
Figure 712957DEST_PATH_IMAGE118
和
Figure 885312DEST_PATH_IMAGE119
分别表示待检卫星和参考星编号，
Figure 978033DEST_PATH_IMAGE120
表示频点，
Figure 981761DEST_PATH_IMAGE121
表示波长，
Figure 270DEST_PATH_IMAGE015
和
Figure 659922DEST_PATH_IMAGE122
表示相邻历元，
Figure 415388DEST_PATH_IMAGE123
表示星间历元间双差载波相位观测值，
Figure 414568DEST_PATH_IMAGE124
表示基于惯导提供的位置得到的双差星地距离，
Figure 964498DEST_PATH_IMAGE125
表示光速，
Figure 111446DEST_PATH_IMAGE126
表示双差卫星钟差，可以通过广播星历计算获得。
Figure 545969DEST_PATH_IMAGE127
可以表示为
Figure 258710DEST_PATH_IMAGE128
Figure 120487DEST_PATH_IMAGE129
其中算子
Figure 285889DEST_PATH_IMAGE130
表示绝对值;
Figure 22639DEST_PATH_IMAGE131
、
Figure 855466DEST_PATH_IMAGE132
分别表示卫星
Figure 216040DEST_PATH_IMAGE133
在
Figure 744104DEST_PATH_IMAGE134
和
Figure 910643DEST_PATH_IMAGE135
历元的坐标;
Figure 473343DEST_PATH_IMAGE136
、
Figure 4818DEST_PATH_IMAGE137
分别表示卫星
Figure 410392DEST_PATH_IMAGE138
在
Figure 724830DEST_PATH_IMAGE139
和
Figure 266669DEST_PATH_IMAGE140
历元的坐标;
Figure 139685DEST_PATH_IMAGE141
和
Figure 970238DEST_PATH_IMAGE142
分别表示惯导提供的站点
Figure 213001DEST_PATH_IMAGE143
在
Figure 750292DEST_PATH_IMAGE144
和
Figure 623570DEST_PATH_IMAGE145
历元的坐标；
Figure 738157DEST_PATH_IMAGE146
表示
Figure 925556DEST_PATH_IMAGE147
和
Figure 176408DEST_PATH_IMAGE148
历元间位置增量。
假设忽略历元间视线矢量变化，可得基础惯辅检验统计量误差为
Figure 892692DEST_PATH_IMAGE149
其中
Figure 494574DEST_PATH_IMAGE150
式中
Figure 282402DEST_PATH_IMAGE151
表示星间历元间双差载波相位观测值误差，
Figure 27242DEST_PATH_IMAGE152
表示
Figure 304639DEST_PATH_IMAGE153
频点波长，
Figure 269184DEST_PATH_IMAGE154
表示
Figure 595123DEST_PATH_IMAGE155
和
Figure 820568DEST_PATH_IMAGE156
历元的模糊度差，即周跳；
Figure 409813DEST_PATH_IMAGE157
和
Figure 923971DEST_PATH_IMAGE158
分别卫星
Figure 850338DEST_PATH_IMAGE159
与站点
Figure 336814DEST_PATH_IMAGE160
在
Figure 893698DEST_PATH_IMAGE161
和
Figure 691889DEST_PATH_IMAGE162
历元的视线矢量;
Figure 795849DEST_PATH_IMAGE163
和
Figure 261466DEST_PATH_IMAGE164
分别卫星
Figure 661354DEST_PATH_IMAGE165
与站点
Figure 150104DEST_PATH_IMAGE166
在
Figure 683854DEST_PATH_IMAGE167
和
Figure 613764DEST_PATH_IMAGE168
历元的视线矢量;
Figure 574767DEST_PATH_IMAGE169
表示
Figure 550813DEST_PATH_IMAGE170
到
Figure 232461DEST_PATH_IMAGE171
历元间状态转移矩阵的逆。
电离层在平静状态下随时间变化平缓，可以通过历元间差分消除。根据公式(3),无故障条件下基础惯辅检验统计量方差为：
Figure 672670DEST_PATH_IMAGE172
其中
Figure 912896DEST_PATH_IMAGE173
表示历元
Figure 376238DEST_PATH_IMAGE174
和
Figure 251790DEST_PATH_IMAGE175
间位置增量协方差阵，
Figure 156292DEST_PATH_IMAGE176
和
Figure 662360DEST_PATH_IMAGE177
分别表示相邻历元
Figure 144157DEST_PATH_IMAGE178
和
Figure 698766DEST_PATH_IMAGE179
的载波相位测量误差方差。
由此可得基础检验门限为：
Figure 582409DEST_PATH_IMAGE180
式中
Figure 259378DEST_PATH_IMAGE181
表示设定的误警率，事先给定；
Figure 103837DEST_PATH_IMAGE182
表示标准正态分布函数的逆函数，定义为
Figure 321192DEST_PATH_IMAGE183
基础周跳探测的检验标准是检验统计量与检验门限做对比，若检验统计量大于则有周跳，反之，则无周跳。基础周跳探测在本文中对后续步骤不产生影响，仅作对比分析用途。
在步骤五和六中，考虑在实时应用中减少估计模型参数的在线计算量，本发明采用被广泛使用的ARIMA (n, 1, 0) 模型。ARIMA (n, 1, 0) 模型的表达式为：
Figure 722259DEST_PATH_IMAGE184
其中
Figure 570129DEST_PATH_IMAGE185
和
Figure 26518DEST_PATH_IMAGE186
分别表示相邻历元的时间序列，
Figure 922930DEST_PATH_IMAGE187
表示第
Figure 718847DEST_PATH_IMAGE188
项ARIMA (n, 1, 0)模型参数，
Figure 534357DEST_PATH_IMAGE189
为白噪声。
基础惯辅检验统计量的趋势项通过
Figure 618987DEST_PATH_IMAGE190
和
Figure 381407DEST_PATH_IMAGE191
的一阶时间差分消除。给定历史时间窗口
Figure 94148DEST_PATH_IMAGE192
,
Figure 955925DEST_PATH_IMAGE193
, 定义窗口内保存的卫星
Figure 386906DEST_PATH_IMAGE194
的基础检验统计量序列为：
Figure 687438DEST_PATH_IMAGE195
Figure 159745DEST_PATH_IMAGE196
式中
Figure 51478DEST_PATH_IMAGE197
表示
Figure 173017DEST_PATH_IMAGE198
历元卫星
Figure 949343DEST_PATH_IMAGE199
在滑动窗口中
Figure 371098DEST_PATH_IMAGE200
历元的基础检验统计量。
对
Figure 840256DEST_PATH_IMAGE201
进行一阶差分，可以得到以下时间序列
Figure 183513DEST_PATH_IMAGE202
其中
Figure 153743DEST_PATH_IMAGE203
现有很多方法可以用来求解
Figure 570949DEST_PATH_IMAGE204
中的模型参数
Figure 273326DEST_PATH_IMAGE205
，例如最小二乘估计，Yule-Walker方法和Burg方法。
在确定ARIMA (n, 1, 0)模型参数
Figure 166195DEST_PATH_IMAGE206
之后，可获得当前历元
Figure 782859DEST_PATH_IMAGE207
与其前一历元
Figure 382468DEST_PATH_IMAGE208
的基础惯辅检验统计量的差分预报值为：
Figure 318063DEST_PATH_IMAGE209
其中，
Figure 308015DEST_PATH_IMAGE210
表示
Figure 292152DEST_PATH_IMAGE211
历元基础惯辅检验统计量一阶差分值；
Figure 74163DEST_PATH_IMAGE212
为基础惯辅检验统计量一阶差分预报值，趋势项通过一阶差分去除。
从而，可得当前历元卫星
Figure 790446DEST_PATH_IMAGE213
的去趋势项后的基础惯辅检验统计量的残差为：
Figure 595591DEST_PATH_IMAGE214
相应的，当前历元卫星
Figure 445736DEST_PATH_IMAGE215
的基础惯辅检验统计量的方差估计为：
Figure 692040DEST_PATH_IMAGE216
Figure 907121DEST_PATH_IMAGE217
其中
Figure 996300DEST_PATH_IMAGE218
表示历史时间滑动窗口
Figure 758457DEST_PATH_IMAGE219
中的样本总数。
在本发明中去趋势项之后的残差用来作为实际的检验统计量。经过大量数据综合分析，本发明最终采用ARIMA (1, 1, 0)模型。
在步骤七中，第一步周跳探测中，检验统计量为式
Figure 187164DEST_PATH_IMAGE220
计算的当前历元的基础惯辅检验统计量的残差
Figure 635463DEST_PATH_IMAGE221
，即
Figure 87304DEST_PATH_IMAGE222
第一步周跳探测中，检验门限计算为：
Figure 216934DEST_PATH_IMAGE223
Figure 296886DEST_PATH_IMAGE224
由式
Figure 57032DEST_PATH_IMAGE225
计算获得，第一步周跳探测中的周跳判断准则为：
Figure 855223DEST_PATH_IMAGE226
在步骤八中，GARCH模型被广泛应用于估计波动聚集现象的时间序列的条件方差。经过大量数据分析，GARCH (1,1)模型可以满足对残差进行条件方差建模的要求，条件方差计算表达式为
Figure 788544DEST_PATH_IMAGE227
上式中假设当前历元是
Figure 863948DEST_PATH_IMAGE228
，式中
Figure 654049DEST_PATH_IMAGE229
表示前一历元卫星
Figure 142799DEST_PATH_IMAGE230
的去趋势项后的基础惯辅检验统计量的残差，
Figure 784871DEST_PATH_IMAGE231
表示前一历元卫星
Figure 104994DEST_PATH_IMAGE232
由GARCH模型计算得到的方差。
Figure 941363DEST_PATH_IMAGE233
。
Figure 448567DEST_PATH_IMAGE234
,
Figure 723691DEST_PATH_IMAGE235
和
Figure 39266DEST_PATH_IMAGE236
均为GARCH (1,1)模型参数，可由迭代的最大似然法计算获得。当
Figure 905591DEST_PATH_IMAGE237
和
Figure 41037DEST_PATH_IMAGE238
时，该模型也可以适用同方差情形。
在步骤九中，最后一步周跳探测中，检验统计量同样为式
Figure 854272DEST_PATH_IMAGE239
计算的当前历元的基础惯辅检验统计量的残差
Figure 148987DEST_PATH_IMAGE240
，即
Figure 560115DEST_PATH_IMAGE241
最后一步周跳探测中，检验门限计算为：
Figure 510753DEST_PATH_IMAGE242
式中
Figure 189996DEST_PATH_IMAGE243
由式
Figure 949005DEST_PATH_IMAGE244
计算获得。最后一步周跳探测中，周跳判断准则为：
Figure 625974DEST_PATH_IMAGE245
本发明提出了针对杆臂不准的实时单站周跳探测改进方法，首先利用ARIMA (1, 1,0) 模型去除基础检验统计量趋势项，然后将GARCH (1,1) 模型用于估计带有波动聚集特性残差的条件方差。针对实时应用，本发明进一步提出了一种模型参数定时更新的方法。该方法不仅可以在较大杆臂误差条件下极大地降低误警率，而且可以准确地检测1周仿真周跳，同时又具有较低的计算量。
为了测试本发明中的方法，采集场地测试实测数据作为仿真试验的基础数据。分别采用Sensonor STIM300 MEMS与ComNav OEM-K508面板采集惯导与GNSS原始观测数据，采样率分别为2hz (GNSS)和125hz (MEMS)。
首先基于事后分析，在场地试验实测数据中精确剔除周跳时刻历元数据，获得无周跳的“干净”数据。在进行试验之前，测量杆臂值并将杆臂误差补偿到0。为了考虑杆臂不准对周跳的影响，假设杆臂误差的三个方向分量都为0均值的高斯分布，且三个方向采用相同大小的方差
Figure 595067DEST_PATH_IMAGE246
。基于文献中对杆臂误差研究的结论, 杆臂误差方差取值为10 cm。
本实施例仅给出北斗系统（BDS）周跳探测结果，GPS系统周跳探测原理完全相同，位置增量由BDS /MEMS 紧组合算法提供，设定测量更新时间间隔为0.5s，使用的卫导观测信息为伪距和多普勒。
设计了两类数字试验分别分析改进方法的误警率和漏检率。第一个试验使用添加杆臂误差的“干净”试验数据来验证本发明可以大量减少误警率的能力。在第二个试验中，仿真的1周周跳等时间间隔添加到第一个实验使用的数据上，用来验证本发明对1周周跳的探测能力。本实施例中的误警率要求为
Figure 687788DEST_PATH_IMAGE247
。
1）设置参数更新周期
Figure 894778DEST_PATH_IMAGE248
和滑动窗口大小
Figure 539386DEST_PATH_IMAGE249
。
以某次试验为例，参数更新周期设置
Figure 871141DEST_PATH_IMAGE250
，滑动窗口大小设置为
Figure 829870DEST_PATH_IMAGE251
。
2）计算基础惯导辅助检验统计量及其方差。
以试验中北斗时
Figure 953684DEST_PATH_IMAGE252
时刻为例，参考星为C05号星，待检卫星C03星B1频点历元间星间双差载波相位、卫星钟差和惯导提供的几何距离分别为：
Figure 143095DEST_PATH_IMAGE253
其中北斗B1频点波长为
Figure 290042DEST_PATH_IMAGE254
，光速为
Figure 583620DEST_PATH_IMAGE255
。
从而，根据基础惯导辅助检验统计量计算公式
Figure 437307DEST_PATH_IMAGE256
，可得：
Figure 626980DEST_PATH_IMAGE257
Figure 57961DEST_PATH_IMAGE258
此时相邻历元双差载波相位观测值误差方差和惯导预报误差项分别为：
Figure 296175DEST_PATH_IMAGE259
根据式
Figure 66685DEST_PATH_IMAGE260
可得检验统计量直接理论方差为：
Figure 489576DEST_PATH_IMAGE261
Figure 283220DEST_PATH_IMAGE262
给定
Figure 121863DEST_PATH_IMAGE263
，由式
Figure 809196DEST_PATH_IMAGE264
可得检验门限为：
Figure 776890DEST_PATH_IMAGE265
若对比式
Figure 120147DEST_PATH_IMAGE266
中检验统计量数值与式
Figure 824798DEST_PATH_IMAGE267
中检验门限数值可知，此时检验不通过。由于“干净”数据中不存在周跳，若采用基础惯导辅助周跳探测方法，则会出现误警。
3）初始化准备样本滑动窗口。
为了在线使用该方法，基于历史检验统计量构造滑动窗口，该窗口内数据必须避免周跳的影响。当ARIMA滑动窗口未准备好时（保存的历史基础惯辅检验统计量数量
Figure 507583DEST_PATH_IMAGE268
小于
Figure 209960DEST_PATH_IMAGE269
），仅使用纯卫导测量信息的经典周跳探测方法如电离层残差方法作为补充。当GARCH滑动窗口未准备好时（前一历元残差及其方差未保存），基于ARIMA模型的第一步周跳探测方法用来避免周跳的负面影响。当所有的滑动窗口准备好之后（保存的历史基础惯辅检验统计量数量
Figure 837250DEST_PATH_IMAGE270
大于或等于
Figure 955379DEST_PATH_IMAGE269
，并且前一历元残差及其方差已保存），所有的模型参数均可以计算获得，最后一步周跳探测可以运行。一旦检验通过，基础检验统计量和ARIMA建模之后的残差用来更新滑动窗口（将当前历元基础惯辅检验统计量保存进入ARIMA的滑动窗口，并保存当前历元计算的趋势项后的残差及其方差）。如果检验不通过，则不更新滑动窗口。在这种策略下，滑动窗口在线更新不会受到周跳影响，从而可以为建模提供更准确的样本。
Figure 820567DEST_PATH_IMAGE271
确定ARIMA和GARCH模型参数。
利用初始滑动窗口历史数据，基于式
Figure 490582DEST_PATH_IMAGE272
中ARIMA (1,1,0)模型，利用最小二乘法可得ARIMA (1,1,0)模型参数为
Figure 746114DEST_PATH_IMAGE273
并且根据式
Figure 730251DEST_PATH_IMAGE274
可得残差方差为
Figure 246683DEST_PATH_IMAGE275
Figure 727080DEST_PATH_IMAGE276
和式
Figure 328963DEST_PATH_IMAGE277
中GARCH(1,1)模型参数计算方法，可得
Figure 116791DEST_PATH_IMAGE278
Figure 628674DEST_PATH_IMAGE279
基于ARIMA(1,1,0)计算去趋势项后的检验统计量残差及其方差。
根据式
Figure 843755DEST_PATH_IMAGE280
可得当前检验统计量残差为
Figure 667355DEST_PATH_IMAGE281
Figure 196556DEST_PATH_IMAGE282
对应的标准差由ARIMA确定模型参数时获得，根据式
Figure 625263DEST_PATH_IMAGE283
可得
Figure 807983DEST_PATH_IMAGE284
Figure 525403DEST_PATH_IMAGE279
第一步周跳探测。
当前时刻前一历元残差及其方差已准备好，因此不需要进行第一步周跳探测。但是为了体现最终一步周跳探测的优势，仍然根据式
Figure 451771DEST_PATH_IMAGE285
和
Figure 734985DEST_PATH_IMAGE286
计算其对应的检验统计量和检验门限为：
Figure 728086DEST_PATH_IMAGE287
根据式
Figure 791857DEST_PATH_IMAGE288
判断准则，可知第一步周跳探测可以避免误警，但是检验统计量残
差方差
Figure 459599DEST_PATH_IMAGE289
在ARIMA确定模型参数时获得，在周期内为常值，因此无法及时反映检验统计量实时的误差统计特性。
Figure 535002DEST_PATH_IMAGE290
基于GARCH(1,1)计算条件方差。
利用算好的模型参数根据式
Figure 590683DEST_PATH_IMAGE291
可得当前历元条件方差为
Figure 813854DEST_PATH_IMAGE292
Figure 222970DEST_PATH_IMAGE293
Figure 543093DEST_PATH_IMAGE294
最后一步周跳探测。
根据式
Figure 113882DEST_PATH_IMAGE295
和
Figure 89929DEST_PATH_IMAGE296
计算的检验统计量和检验门限为
Figure 427369DEST_PATH_IMAGE297
根据式
Figure 539682DEST_PATH_IMAGE298
判断准则，可知最后一步周跳探测可以避免误警。相比于第一步周跳探测，检验门限与检验统计量相差更小，可以更好地避免漏检。
Figure 45487DEST_PATH_IMAGE299
更新ARIMA和GARCH滑动窗口。
若检验通过，基础检验统计量和ARIMA建模之后的残差用来更新滑动窗口（将当前历元基础惯辅检验统计量保存进入ARIMA的滑动窗口，并保存当前历元计算的去趋势项后的残差及其方差）。若检验不通过，则不更新滑动窗口。
上述时刻算例直观地展示了本发明方法降低误警率的优势。两类仿真试验的统计结果可以进一步验证方法的良好性能。为了分析方便，传统基础方法和本发明中的第一步周跳探测方法分别定义为 ‘case 1’和‘case 2’。将本发明中模型参数定时变化的方法简称为 ‘case 3’。
表1 不同方法误警历元个数比较
Figure 305567DEST_PATH_IMAGE300
这里重点分析两颗代表性卫星的误警。C03是GEO卫星，C08是MEO卫星。杆臂误差设置了不同等级的标准差。历元表示周跳探测时间个数的单位。ARIMA滑动窗口大小为30个历元。Case 3 中滑动窗口大小为150个历元，模型参数在整个实验中保持不变。如表 1所示，可以发现
Figure 56485DEST_PATH_IMAGE301
值越大误警历元数越多，这反映了杆臂误差对周跳探测误警的影响。相比于case1和case 2，case 3对两颗卫星的两个频点都大幅减少了，其误警历元个数如图2和3所示，其中图2是不同方法北斗C03卫星的B1频点的误警历元个数比较图，图3是不同方法北斗C08卫星的B1频点的误警历元个数比较图。Case 2 中利用ARIMA (1, 1, 0)模型去趋势项初步避免了大量的误警。基于GARCH (1,1)模型估计自适应的条件方差，case 3中误警历元个数进一步减少。
表2 不同方法漏检历元个数比较
Figure 288884DEST_PATH_IMAGE302
为了分析所提方法探测小周跳的性能，
Figure 591689DEST_PATH_IMAGE303
周的仿真周跳随机产生并每隔3s独立地添加到不同卫星频点的载波相位测量值上。表2显示了场地试验中连续可见卫星的漏检历元个数。杆臂误差标准差设置为10 cm。添加过仿真周跳的历元个数为151。由表2可知，对于所有卫星的所有频点漏检历元个数为0。这意味着所有方法都可以探测1周周跳。图4显示了所提方法的性能，图中实心圆圈表示添加周跳的历元。 所有的1周小周跳都可以被改进方法识别。基于以上分析可知，本发明所提供的方法可以满足可靠性要求探测仿真的BDS所有卫星频点是1周小周跳。
图2 和图3显示了C03和C08号星B1频点不同方法误警历元个数比较。从图中可知，相比于case 1、case 2和case 3对两颗卫星的两个频点都大幅减少了误警历元个数如图2所示。Case 2中利用ARIMA (1, 1, 0)模型去趋势项初步避免了大量的误警。基于GARCH(1,1)模型估计自适应的条件方差，case 3中误警历元个数进一步减少。相比于case 1，case 3总是可以显著的减少误警历元个数，误警历元从几百下降到个位数。
图4显示了试验中典型时间段内无周跳条件下北斗C03卫星B1频点不同方法门限与检验统计量的比较，其中
Figure 152114DEST_PATH_IMAGE304
和
Figure 769041DEST_PATH_IMAGE305
分别表示case
Figure 652683DEST_PATH_IMAGE306
的检验统计量和门限。由图4可知，本发明方法case 3的检验门限可以更好地试验检验统计量误差特性变化，而case 2中检验统计量超过检验门限的次数较多。这直观地证明了本发明方法更加适应检验统计量误差特性变化，可以大幅降低误警能力。
图5和图6显示了仿真添加
Figure 500291DEST_PATH_IMAGE307
周小周跳条件下，不同方法北斗C08卫星的B1和B3频点的检验统计量与检验门限对比，其中图5显示添加
Figure 469384DEST_PATH_IMAGE308
周小周跳条件下，不同方法北斗C08卫星的B1频点的检验统计量与检验门限对比；图6显示了添加周小周跳条件下，不同方法北斗C08卫星的B3频点的检验统计量与检验门限对比。图中实心圆圈表示添加周跳的历元。由图可知，所有的1周小周跳都可以被本发明方法识别。
基于以上分析可知，本发明方法不仅可以在较大杆臂误差条件下极大地降低误警率，而且可以满足可靠性要求探测仿真的所有卫星频点的1周小周跳。
综上所述，虽然本发明已以较佳实施例揭露如上，然其并非用以限定本发明，任何本领域普通技术人员，在不脱离本发明的精神和范围内，当可作各种更动与润饰，因此本发明的保护范围当视权利要求书界定的范围为准。

Claims (10)

1.针对杆臂不准的实时单站周跳探测改进方法，其特征在于，包括以下步骤：
步骤一：设置参数更新周期
Figure 387325DEST_PATH_IMAGE001
和滑动窗口大小
Figure 835624DEST_PATH_IMAGE002
；
步骤二：计算当前历元的基础惯辅检验统计量及其方差，并保存在滑动窗口；
步骤三：当滑动窗口中保存的历史历元的基础惯辅检验统计量数量
Figure 349782DEST_PATH_IMAGE003
≥
Figure 151516DEST_PATH_IMAGE004
且当前建模时间
Figure 762626DEST_PATH_IMAGE005
＞
Figure 991613DEST_PATH_IMAGE006
，则计算ARIMA模型参数和GARCH模型参数，然后初始化建模时刻
Figure 993067DEST_PATH_IMAGE007
，并进入下一步；若
Figure 723126DEST_PATH_IMAGE008
≤
Figure 798529DEST_PATH_IMAGE009
，则直接进入下一步；
步骤四：基于滑动窗口中保存的历史惯辅检验统计量，利用ARIMA模型通过一阶时间差分对基础惯辅检验统计量去趋势项，并计算去趋势项后的残差及其方差；
步骤五：周跳探测；
判断前一历元去趋势项后的基础惯辅检验统计量的残差及其方差是否保存，若前一历元去趋势项后的基础惯辅检验统计量的残差及其方差没有保存，则进行第一步周跳探测，若检验不通过，则有周跳；若检验通过，则无周跳，并更新滑动窗口，周跳探测结束；
若前一历元去趋势项后的基础惯辅检验统计量的残差及其方差已保存，则基于GARCH模型计算条件方差，进入第二步周跳探测，若检验不通过，则有周跳；若检验通过，则无周跳，并更新滑动窗口，周跳探测结束。
2.根据权利要求1所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，所述ARIMA模型为ARIMA (n, 1, 0)模型，ARIMA模型参数为ARIMA (n, 1, 0)模型参数。
3.根据权利要求2所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，ARIMA (n, 1, 0)模型参数由最小二乘估计法、Yule-Walker方法或Burg方法计算得到。
4.根据权利要求1所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，GARCH模型为GARCH (1,1)模型，GARCH模型为GARCH (1,1)模型参数。
5.根据权利要求4所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，GARCH (1,1)模型参数采用极大似然法计算得到。
6.根据权利要求1所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，基础惯辅检验统计量的计算方法是：
Figure 588631DEST_PATH_IMAGE010
其中
Figure 77381DEST_PATH_IMAGE011
和
Figure 719453DEST_PATH_IMAGE012
分别表示待检卫星和参考星编号，
Figure 39576DEST_PATH_IMAGE013
表示频点，
Figure 875944DEST_PATH_IMAGE014
表示波长，
Figure 851991DEST_PATH_IMAGE015
和
Figure 923852DEST_PATH_IMAGE016
表示相邻历元，
Figure 973848DEST_PATH_IMAGE017
表示星间历元间双差载波相位观测值，
Figure 43435DEST_PATH_IMAGE018
表示基于惯导提供的位置得到的双差星地距离，
Figure 303515DEST_PATH_IMAGE019
表示光速，
Figure 54433DEST_PATH_IMAGE020
表示双差卫星钟差，通过广播星历计算获得。
7.根据权利要求4所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，基础惯辅检验统计量方差为：
Figure 83569DEST_PATH_IMAGE021
其中
Figure 589637DEST_PATH_IMAGE022
表示历元
Figure 445335DEST_PATH_IMAGE023
和
Figure 858999DEST_PATH_IMAGE024
间位置增量协方差阵；
Figure 883587DEST_PATH_IMAGE025
和
Figure 357293DEST_PATH_IMAGE026
分别表示相邻历元
Figure 201752DEST_PATH_IMAGE027
和
Figure 356790DEST_PATH_IMAGE028
的载波相位测量误差方差；
Figure 360518DEST_PATH_IMAGE029
，
Figure 146072DEST_PATH_IMAGE030
和
Figure 805723DEST_PATH_IMAGE031
分别卫星
Figure 561190DEST_PATH_IMAGE032
与站点
Figure 793326DEST_PATH_IMAGE033
在
Figure 608835DEST_PATH_IMAGE034
和
Figure 755782DEST_PATH_IMAGE035
历元的视线矢量；
Figure 190306DEST_PATH_IMAGE036
和
Figure 168626DEST_PATH_IMAGE037
分别卫星
Figure 92720DEST_PATH_IMAGE038
与站点
Figure 399067DEST_PATH_IMAGE039
在
Figure 761916DEST_PATH_IMAGE040
和
Figure 470109DEST_PATH_IMAGE041
历元的视线矢量，
Figure 565104DEST_PATH_IMAGE042
表示
Figure 483381DEST_PATH_IMAGE043
到
Figure 23822DEST_PATH_IMAGE044
历元间状态转移矩阵的逆。
8.根据权利要求1所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，第四步的实现方法是：
ARIMA (n, 1, 0)表达式为：
Figure 445576DEST_PATH_IMAGE045
其中
Figure 977051DEST_PATH_IMAGE046
和
Figure 257991DEST_PATH_IMAGE047
分别表示相邻历元的时间序列，
Figure 962642DEST_PATH_IMAGE048
表示第
Figure 707744DEST_PATH_IMAGE049
项ARIMA (n, 1, 0)模型参数，
Figure 82224DEST_PATH_IMAGE050
为白噪声；
基础惯辅检验统计量的趋势项通过
Figure 975094DEST_PATH_IMAGE051
和
Figure 421119DEST_PATH_IMAGE052
的一阶时间差分消除；给定历史时间滑动窗口
Figure 692831DEST_PATH_IMAGE053
, 定义窗口内保存的卫星
Figure 628426DEST_PATH_IMAGE054
的基础检验统计量序列为：
Figure 116914DEST_PATH_IMAGE055
式中
Figure 101051DEST_PATH_IMAGE056
表示
Figure 617483DEST_PATH_IMAGE057
历元卫星
Figure 599345DEST_PATH_IMAGE058
在滑动窗口中
Figure 404490DEST_PATH_IMAGE057
历元的基础检验统计量；
对
Figure 989055DEST_PATH_IMAGE059
进行一阶差分，得到以下时间序列：
Figure 500939DEST_PATH_IMAGE060
其中：
Figure 512758DEST_PATH_IMAGE061
在确定ARIMA (n, 1, 0)模型参数
Figure 539620DEST_PATH_IMAGE062
之后，可获得当前历元
Figure 803242DEST_PATH_IMAGE063
与其前一历元
Figure 294266DEST_PATH_IMAGE064
的基础惯辅检验统计量的差分预报值为：
Figure 850887DEST_PATH_IMAGE065
其中，
Figure 630624DEST_PATH_IMAGE066
表示
Figure 556992DEST_PATH_IMAGE067
历元基础惯辅检验统计量一阶差分值；
Figure 777889DEST_PATH_IMAGE068
为基础惯辅检验统计量一阶差分预报值，趋势项通过一阶差分去除；
从而，可得当前历元卫星
Figure 397089DEST_PATH_IMAGE069
的去趋势项后的基础惯辅检验统计量的残差为：
Figure 398543DEST_PATH_IMAGE070
相应的，当前历元卫星
Figure 269547DEST_PATH_IMAGE071
的基础惯辅检验统计量的方差估计为：
Figure 204005DEST_PATH_IMAGE072
其中
Figure 869473DEST_PATH_IMAGE073
表示历史时间滑动窗口
Figure 358223DEST_PATH_IMAGE074
中的样本总数。
9.根据权利要求8所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，步骤五中，第一步周跳探测方法为：
第一步周跳探测中检验统计量为当前历元的基础惯辅检验统计量的残差
Figure 891972DEST_PATH_IMAGE075
，即
Figure 320417DEST_PATH_IMAGE076
第一步周跳探测中的检验门限计算为：
Figure 281420DEST_PATH_IMAGE077
其中，
Figure 991887DEST_PATH_IMAGE078
表示设定的误警率；
Figure 204694DEST_PATH_IMAGE079
表示标准正态分布函数的逆函数；
第一步周跳探测中，周跳判断准则为：
Figure 379323DEST_PATH_IMAGE080
。
10.根据权利要求8所述的针对杆臂不准的实时单站周跳探测改进方法，其特征在于，步骤五中，第二步周跳探测方法为：
基于GARCH模型计算条件方差，公式为：
Figure 448910DEST_PATH_IMAGE081
上式中假设当前历元是
Figure 584357DEST_PATH_IMAGE082
，式中
Figure 459909DEST_PATH_IMAGE083
表示前一历元卫星
Figure 364411DEST_PATH_IMAGE084
的去趋势项后的基础惯辅检验统计量的残差，
Figure 401637DEST_PATH_IMAGE085
表示前一历元卫星
Figure 352275DEST_PATH_IMAGE086
由GARCH模型计算得到的方差，
Figure 405420DEST_PATH_IMAGE087
,
Figure 289062DEST_PATH_IMAGE088
,
Figure 966031DEST_PATH_IMAGE089
和
Figure 810490DEST_PATH_IMAGE090
均为GARCH (1,1)模型参数；
第二步周跳探测中检验统计量为当前历元的基础惯辅检验统计量的残差
Figure 27845DEST_PATH_IMAGE091
即
Figure 234836DEST_PATH_IMAGE092
第二步周跳探测中的检验门限计算为：
Figure 489231DEST_PATH_IMAGE093
其中，
Figure 211199DEST_PATH_IMAGE094
表示设定的误警率；
Figure 842031DEST_PATH_IMAGE095
表示标准正态分布函数的逆函数；
第二步周跳探测中，周跳判断准则为：
Figure 903528DEST_PATH_IMAGE096
。
CN202011175778.8A 2020-10-29 2020-10-29 针对杆臂不准的实时单站周跳探测改进方法 Active CN111999750B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011175778.8A CN111999750B (zh) 2020-10-29 2020-10-29 针对杆臂不准的实时单站周跳探测改进方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011175778.8A CN111999750B (zh) 2020-10-29 2020-10-29 针对杆臂不准的实时单站周跳探测改进方法

Publications (2)

Publication Number Publication Date
CN111999750A true CN111999750A (zh) 2020-11-27
CN111999750B CN111999750B (zh) 2021-01-08

Family

ID=73475232

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011175778.8A Active CN111999750B (zh) 2020-10-29 2020-10-29 针对杆臂不准的实时单站周跳探测改进方法

Country Status (1)

Country Link
CN (1) CN111999750B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112526573A (zh) * 2021-02-07 2021-03-19 腾讯科技（深圳）有限公司 对象定位方法和装置、存储介质及电子设备

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20100299066A1 (en) * 2009-05-20 2010-11-25 Gang Kevin Liu INS based GPS carrier phase cycle slip detection and repairing
CN104375157A (zh) * 2014-10-22 2015-02-25 南京航空航天大学 短基线下惯导辅助的北斗单频整周模糊度求解方法
CN106842236A (zh) * 2015-12-04 2017-06-13 航天恒星科技有限公司 Gnss接收机周跳探测与修复处理方法及装置
CN107505642A (zh) * 2017-10-23 2017-12-22 中国矿业大学 一种ins辅助的实时bds单频周跳探测方法
CN110441801A (zh) * 2019-09-03 2019-11-12 中国人民解放军国防科技大学 一种基于最优固定概率的三频周跳探测与修复方法
CN110441800A (zh) * 2019-09-03 2019-11-12 中国人民解放军国防科技大学 一种基于线性组合优选的四频周跳探测与修复方法
WO2019216261A1 (ja) * 2018-05-10 2019-11-14 ソニー株式会社 情報処理装置、情報処理方法およびプログラム
CN111781614A (zh) * 2020-08-25 2020-10-16 中国石油大学(华东) 一种周跳探测与修复方法、装置、设备及可存储介质

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20100299066A1 (en) * 2009-05-20 2010-11-25 Gang Kevin Liu INS based GPS carrier phase cycle slip detection and repairing
CN104375157A (zh) * 2014-10-22 2015-02-25 南京航空航天大学 短基线下惯导辅助的北斗单频整周模糊度求解方法
CN106842236A (zh) * 2015-12-04 2017-06-13 航天恒星科技有限公司 Gnss接收机周跳探测与修复处理方法及装置
CN107505642A (zh) * 2017-10-23 2017-12-22 中国矿业大学 一种ins辅助的实时bds单频周跳探测方法
WO2019216261A1 (ja) * 2018-05-10 2019-11-14 ソニー株式会社 情報処理装置、情報処理方法およびプログラム
CN110441801A (zh) * 2019-09-03 2019-11-12 中国人民解放军国防科技大学 一种基于最优固定概率的三频周跳探测与修复方法
CN110441800A (zh) * 2019-09-03 2019-11-12 中国人民解放军国防科技大学 一种基于线性组合优选的四频周跳探测与修复方法
CN111781614A (zh) * 2020-08-25 2020-10-16 中国石油大学(华东) 一种周跳探测与修复方法、装置、设备及可存储介质

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
HAIFENG HAN ET AL.: "Millimeter-Level Precision Deformation Monitoring Algorithms and Implementation of Beidou Satellite Navigation System", 《IOP CONF. SERIES: EARTH AND ENVIRONMENTAL SCIENCE376》 *
JINGZHOU XIN ET AL.: "Bridge Structure Deformation Prediction Based on GNSS Data Using Kalman-ARIMA-GARCH Model", 《SENSORS FOR TRANSPORTATION》 *
孙冲等: "单频GPS 相位观测值周跳探测方法研究", 《水利与建筑工程学报》 *
汤勇刚: "载波相位时间差分/捷联惯导组合导航方法研究", 《中国优秀硕士学位论文全文数据库·基础科学辑》 *
马朝忠等: "基于ARIMA模型的卫星钟差异常值探测的模型选择方法", 《武汉大学学报·信息科学版》 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112526573A (zh) * 2021-02-07 2021-03-19 腾讯科技（深圳）有限公司 对象定位方法和装置、存储介质及电子设备

Also Published As

Publication number Publication date
CN111999750B (zh) 2021-01-08

Similar Documents

Publication Publication Date Title
US8106823B2 (en) Method of operating a satellite navigation receiver
US7501981B2 (en) Methods and apparatus to detect and correct integrity failures in satellite positioning system receivers
CN105487088B (zh) 一种卫星导航系统中基于卡尔曼滤波的raim算法
CN101629997A (zh) 惯性辅助卫星导航完好性检测装置及检测方法
US8325086B2 (en) Methods and systems to diminish false-alarm rates in multi-hypothesis signal detection through combinatoric navigation
CN109085619B (zh) 多模gnss系统的定位方法及装置、存储介质、接收机
EP2037291A1 (en) Integrity monitoring method for GNSS navigation based on historical information
CN113138402B (zh) 基于rtk的模糊度固定方法及装置、存储介质
EP2685214B1 (en) Multiple truth reference system and method
CN115420284B (zh) 一种组合导航系统故障检测与识别方法
CN109633703B (zh) 一种应对遮挡场景的北斗导航无源定位方法
CN111623703A (zh) 一种基于新型卡尔曼滤波的北斗变形监测实时处理方法
CN117724125B (zh) 一种基于一致性的观测数据的质量控制方法及装置
CN111999750B (zh) 针对杆臂不准的实时单站周跳探测改进方法
CN103926596B (zh) 一种基于粒子滤波的稳健gnss抗欺骗方法
Wen et al. Mitigation of multiple outliers using consistency checking for GNSS standard point positioning in urban areas
CN116859415A (zh) 一种快速稳健高精度的多故障卫星识别及定位方法
EP3882665A1 (en) Implementing single differences within a solution separation framework
CN112926190B (zh) 一种基于vmd算法的多路径削弱方法与装置
CN114690210A (zh) 一种基于多普勒观测值的北斗卫星机动探测方法
CN111399019B (zh) 一种gnss中卫星定位方法和装置
CN105068092B (zh) 一种应用于星基增强系统机载接收机的周跳检测与修复方法
CN116736339B (zh) 一种面向禁航区管控的北斗自主监测及预警方法
US20240272303A1 (en) System and method for gnss spoofer detection with high integrity error bounding using multiple coasting solutions
US11714200B2 (en) Single-difference based pre-filter of measurements for use in solution separation framework

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
GR01 Patent grant
GR01 Patent grant