CN110058281B - 动态定位方法及装置 - Google Patents

动态定位方法及装置 Download PDF

Info

Publication number
CN110058281B
CN110058281B CN201910356541.0A CN201910356541A CN110058281B CN 110058281 B CN110058281 B CN 110058281B CN 201910356541 A CN201910356541 A CN 201910356541A CN 110058281 B CN110058281 B CN 110058281B
Authority
CN
China
Prior art keywords
epoch
variance
satellite
state vector
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.)
Active
Application number
CN201910356541.0A
Other languages
English (en)
Other versions
CN110058281A (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.)
Hunan Goke Microelectronics Co Ltd
Original Assignee
Hunan Goke Microelectronics Co Ltd
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 Hunan Goke Microelectronics Co Ltd filed Critical Hunan Goke Microelectronics Co Ltd
Priority to CN201910356541.0A priority Critical patent/CN110058281B/zh
Publication of CN110058281A publication Critical patent/CN110058281A/zh
Application granted granted Critical
Publication of CN110058281B publication Critical patent/CN110058281B/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/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/43—Determining position using carrier phase measurements, e.g. kinematic positioning; using long or short baseline interferometry

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

本发明实施例提出一种动态定位方法及装置，涉及定位技术领域。该方法包括：计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，该卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星；依据信噪比方差值以及预先建立的卡尔曼滤波模型计算移动站在当前历元的位置。由于低仰角卫星更易受多路径的影响，所以在计算信噪比方差值前先剔除高度角低于预设角度的卫星，有效降低多路径效应的影响；而将信噪比方差值参与到卡尔曼滤波的计算过程，能够更好地对多路径误差进行滤波和消除，从而提高RTK定位结果的准确性。

Description

动态定位方法及装置
技术领域
本发明涉及定位技术领域，具体而言，涉及一种动态定位方法及装置。
背景技术
RTK(Real Time Kinematic，实时动态)定位技术，是实时处理两个测量站载波相位观测量的差分方法，将基准站采集的载波相位发给用户接收机，进行计算求得精确的站间距离。多路径误差是指在卫星信号传输到定位装置过程中，除了直射信号外还有其他反射信号被定位装置捕获跟踪，导致的定位误差。在RTK定位中，多路径误差是影响定位结果的主要误差之一，特别是在“城市峡谷”等复杂环境下，高楼林立，多路径频繁发生，且在动态环境下变化复杂，并且无法通过差分来消除或削弱，严重影响RTK定位结果。
目前，多路径影响可以在天线端、信号处理等方面采取一些措施进行处理。在天线端采用扼流圈天线对多路径信号有很好的屏蔽效果，但扼流圈天线一般体积较大，在很多场合中应用受限；在信号处理中，采用窄相关技术只能检测一部分的多径信号，多路径抑制效果较差。
发明内容
本发明实施例的目的在于提供一种动态定位方法及装置，以解决复杂环境下，RTK定位受多路径的影响而导致的定位误差问题。
为了实现上述目的，本发明实施例采用的技术方案如下：
第一方面，本发明提出一种动态定位方法，所述方法包括：
计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，所述卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星；
依据所述信噪比方差值以及预先建立的卡尔曼滤波模型计算所述移动站在当前历元的位置。
第二方面，本发明还提出一种动态定位装置，所述装置包括：
信噪比方差计算模块，用于计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，所述卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星；
位置计算模块，用于依据所述信噪比方差值以及预先建立的卡尔曼滤波模型计算所述移动站在当前历元的位置。
相对现有技术，本发明实施例具有以下有益效果：
本发明实施例提供的动态定位方法及装置，通过计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，所述卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星；依据所述信噪比方差值以及预先建立的卡尔曼滤波模型计算所述移动站在当前历元的位置。由于低仰角卫星更易受多路径的影响，所以在计算信噪比方差值前先剔除高度角低于预设角度的卫星，有效降低多路径效应的影响；而将信噪比方差值参与到卡尔曼滤波的计算过程，能够更好地对多路径误差进行滤波和消除，从而提高RTK定位结果的准确性。
本发明的其他特征和优点将在随后的说明书阐述，并且，部分地从说明书中变得显而易见，或者通过实施本发明实施例了解。本发明的目的和其他优点可通过在所写的说明书、权利要求书、以及附图中所特别指出的结构来实现和获得。
附图说明
为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，应当理解，以下附图仅示出了本发明的某些实施例，因此不应被看作是对范围的限定，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。
图1示出了本发明实施例所提供的动态定位方法及装置的应用环境示意图。
图2示出了图1所示的移动站的一种结构框图。
图3示出了本发明实施例所提供的一种动态定位方法的流程示意图。
图4示出了图3中步骤S102的具体流程示意图。
图5示出了图4中步骤S1022的具体流程示意图。
图6示出了本发明实施例所提供的另一种动态定位方法的流程示意图。
图7示出了本发明实施例所提供的动态定位装置的功能模块示意图。
图标：100-移动站；200-基准站；400-动态定位装置；110-存储器；120-处理器；130-通信接口；410-信噪比方差计算模块；420-位置计算模块；430-验后残差处理模块。
具体实施方式
下面将结合本发明实施例中附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。通常在此处附图中描述和示出的本发明实施例的组件可以以各种不同的配置来布置和设计。因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明的实施例，本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。
应注意到：相似的标号和字母在下面的附图中表示类似项，因此，一旦某一项在一个附图中被定义，则在随后的附图中不需要对其进行进一步定义和解释。同时，在本发明的描述中，术语“第一”、“第二”等仅用于区分描述，而不能理解为指示或暗示相对重要性。
请参照图1，为本发明实施例所提供的动态定位方法及装置的一种应用环境示意图。移动站100和基准站200通信连接，以实现移动站100与基准站200之间的数据通信或交互，移动站100和基准站200均可以观测到多个卫星(例如，卫星A、卫星B、卫星C、卫星D等)。基准站200通常置于经度、纬度、高度均已知的地方(即基准站200的位置已知)，基准站200采集卫星数据，并实时将基准站200的卫星观测值传送给移动站100；相应地，移动站100采集卫星数据，并根据获取的移动站100的卫星观测值和接收的基准站200的卫星观测值解算移动站100的实时位置。
在本实施例中，该移动站100可以是，但不限于，无人机、汽车、智能终端(例如，手机、平板电脑等)等设备。
请参照图2，为图1所示的移动站100的一种结构框图。该移动站100可以包括存储器110、处理器120和通信接口130，该存储器110、处理器120和通信接口130，各元件相互之间直接或间接地电性连接，以实现数据的传输或交互。例如，这些元件相互之间可通过一条或多条通讯总线或信号线实现电性连接。存储器110可用于存储软件程序及模块，如本发明实施例所提供的动态定位方法及装置对应的程序指令/模块，处理器120通过执行存储在存储器110内的软件程序及模块，从而执行各种功能应用以及数据处理。该通信接口130可用于与其他节点设备进行信令或数据的通信。
其中，存储器110可以是但不限于，随机存取存储器(Random Access Memory，RAM)，只读存储器(Read Only Memory，ROM)，可编程只读存储器(Programmable Read-OnlyMemory，PROM)，可擦除只读存储器(Erasable Programmable Read-Only Memory，EPROM)，电可擦除只读存储器(Electric Erasable Programmable Read-Only Memory，EEPROM)等。
处理器120可能是一种集成电路芯片，具有信号处理能力。该处理器120可以是通用处理器，包括中央处理器(Central Processing Unit，CPU)、网络处理器(NetworkProcessor，NP)等；还可以是数字信号处理器(DSP))、专用集成电路(ASIC)、现成可编程门阵列(FPGA)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件。
可以理解，图2所示的结构仅为示意，移动站100还可包括比图2中所示更多或者更少的组件，或者具有与图2所示不同的配置。图2中所示的各组件可以采用硬件、软件或其组合实现。
请参照图3，为本发明实施例所提供的动态定位方法的流程示意图。需要说明的是，本发明实施例所述的动态定位方法并不以图3以及以下所述的具体顺序为限制，应当理解，在其它实施例中，本发明所述的动态定位方法其中部分步骤的顺序可以根据实际需要相互交换，或者其中的部分步骤也可以省略或删除。本发明实施例所述的动态定位方法可应用于上述的移动站100中，下面将对图3所示的具体流程进行详细阐述。
步骤S101，计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，该卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星。
可以理解，移动站100在计算信噪比方差值之前，需要选取该基准站200和移动站100的共视卫星(能够同时观测到的卫星)，并将共视卫星中高度角低于预设角度的卫星剔除，从而确定出最终参与信噪比方差值计算的卫星。其中，该预设角度可以设置为15°，由于低仰角卫星更易受多路径的影响，故通过剔除高度角较低的卫星，可以降低多路径效应的影响。
当卫星信号中产生多路径效应时，该卫星的信噪比会出现相应的波动，根据这一特性，本申请设计一个固定时间长度的滑动窗口(即上述的预设滑动窗口)，通过计算每颗卫星在预设滑动窗口内的信噪比方差值，将该信噪比方差值作为检验量，来统计每颗卫星在预设滑动窗口内的信噪比波动情况。
在本实施例中，信噪比方差值可以通过下式(1-1)进行计算：
Figure GDA0002822260940000061
其中，N为该预设滑动窗口内的信噪比观测值的总数，xi为该预设滑动窗口内的信噪比观测值，u为该预设滑动窗口内的信噪比观测值的平均值。
步骤S102，依据信噪比方差值以及预先建立的卡尔曼滤波模型计算移动站在当前历元的位置。
如图4所述，在本实施例中，该步骤S102具体包括如下步骤：
步骤S1021，依据信噪比方差值计算每颗卫星的信噪比方差系数。
例如，假设在某一历元Ti时刻，观测到k颗高度角在15°以上的卫星，则会得到k个信噪比方差值[σ1 σ2 … σk]，找出其中的最小值σmin，将[σ1 σ2 … σk]中的全体值分别除以σmin并平方得到下式(1-2)：
Figure GDA0002822260940000071
在本实施例中，l1 l2...lk表示每颗卫星对应的信噪比方差系数，用于在卡尔曼测量更新时与对应卫星的测量值方差相乘。
步骤S1022，依据信噪比方差系数以及卡尔曼滤波模型计算移动站在当前历元的位置。
在动态定位中，多路径效应在不同的场景中对定位结果的影响也不一样，当周围环境快速变化，多路径表现为随机特性，其对定位结果的影响更类似于白噪声，故本申请采用卡尔曼滤波能更好地对多路径误差进行滤波和消除。卡尔曼滤波的系统方程为：
Figure GDA0002822260940000072
式(1-3)中，m和m-1分别表示第m个历元和第m-1个历元，
Figure GDA0002822260940000073
表示系统的状态向量，包括位置、速度、双差载波模糊度等待估参数；
Figure GDA0002822260940000074
表示测量信息，包括伪距观测向量ymp和载波相位观测向量
Figure GDA0002822260940000075
Figure GDA0002822260940000076
表示系统的状态转移矩阵，ΔT表示历元时间间隔，I表示单位矩阵，可以理解，I3表示3阶单位矩阵，In表示n阶单位矩阵；
Figure GDA0002822260940000081
表示系统的噪声系数矩阵；wm-1表示系统的过程噪声向量，vm表示观测噪声，系统的过程噪声wm-1和观测噪声vm是均值为零、互不相关的白噪声，其可以表述为：
E[wm]＝E[vm]＝0
E[wmwj]＝Qwδmj
E[vmvj]＝Rmδmj
cov[wm,vj]＝0
Figure GDA0002822260940000082
其中，E表示求均值，cov表示求协方差，Qw表示过程噪声协方差，
Figure GDA0002822260940000083
表示观测噪声方差，Rp表示伪距观测值的测量方差阵，
Figure GDA0002822260940000084
表示载波相位观测值的测量方差阵，
Figure GDA0002822260940000085
d为伪距观测值与载波相位观测值的测量精度比例，
Figure GDA0002822260940000086
其中，
Figure GDA0002822260940000087
Figure GDA0002822260940000088
Figure GDA0002822260940000089
为基准站200对应的参考卫星的非差观测值的方差，
Figure GDA00028222609400000810
为移动站100对应的参考卫星的非差观测值的方差，
Figure GDA0002822260940000091
为基准站200对应的第i颗卫星的非差观测值的方差，
Figure GDA0002822260940000092
为移动站100对应的第i颗卫星的非差观测值的方差，lr为根据式(1-2)得到的参考卫星的信噪比方差系数，li为根据式(1-2)得到的第i颗卫星的信噪比方差系数。
值得注意的是，在本实施例中，非差观测值的方差可利用公式σ＝a2+b2/(sin(el))2来计算，其中，a、b为固定系数，el为对应卫星的高度角。
式(1-3)中，
Figure GDA0002822260940000093
表示测量系数矩阵且该测量系数矩阵是根据获取的伪距观测值和载波相位观测值计算得到，λ为常数。
具体地，
Figure GDA0002822260940000094
Figure GDA0002822260940000097
表示双差伪距观测值，
Figure GDA0002822260940000095
Figure GDA0002822260940000096
表示双差载波相位观测值。
如图5所示，在本实施例中，卡尔曼滤波模型包括状态预测方程和测量更新方程，该步骤S1022具体包括如下子步骤：
子步骤S11，依据状态预测方程计算当前历元的状态向量预测值和当前历元的状态向量方差的预测值。
在本实施例中，该状态预测方程包括第一预测方程
Figure GDA0002822260940000101
和第二预测方程
Figure GDA0002822260940000102
移动站100根据第一预测方程
Figure GDA0002822260940000103
计算m历元的状态向量预测值，并将m历元的状态向量预测值作为当前历元的状态向量预测值，其中，
Figure GDA0002822260940000104
表示m历元的状态向量预测值，
Figure GDA0002822260940000105
表示m-1历元的状态向量滤波值，
Figure GDA0002822260940000106
表示状态转移矩阵，I表示单位矩阵，ΔT表示历元时间间隔；移动站100根据第二预测方程
Figure GDA0002822260940000107
计算m历元的状态向量方差的预测值，并将m历元的状态向量方差的预测值作为当前历元的状态向量方差的预测值；其中，
Figure GDA0002822260940000108
表示m历元的状态向量方差的预测值，
Figure GDA0002822260940000109
表示m-1历元的状态向量方差的滤波值，
Figure GDA00028222609400001010
表示噪声系数矩阵，Qm-1表示m-1历元的过程噪声，
Figure GDA00028222609400001011
表示状态转移矩阵Φm-1的转置矩阵，
Figure GDA00028222609400001012
表示噪声系数矩阵Γm-1的转置矩阵。
子步骤S12，依据信噪比方差系数、当前历元的状态向量预测值、当前历元的状态向量方差的预测值以及测量更新方程计算移动站在当前历元的位置。
在本实施例中，该测量更新方程包括增益方程
Figure GDA0002822260940000111
状态向量更新方程
Figure GDA0002822260940000112
以及方差更新方程
Figure GDA0002822260940000113
移动站100根据m历元的状态向量方差的预测值
Figure GDA0002822260940000114
及增益方程
Figure GDA0002822260940000115
计算卡尔曼增益矩阵；其中，km表示卡尔曼增益矩阵，Hm表示测量系数矩阵且测量系数矩阵是根据获取的伪距观测值和载波相位观测值计算得到，
Figure GDA0002822260940000116
表示m历元的状态向量方差的预测值
Figure GDA0002822260940000117
的逆矩阵，
Figure GDA0002822260940000118
表示观测噪声方差矩阵Rm的逆矩阵，
Figure GDA0002822260940000119
表示测量系数矩阵Hm的转置矩阵，Rm表示观测噪声方差且观测噪声方差是根据信噪比方差系数计算得到(Rm为包含信噪比方差系数的观测噪声方差矩阵)。
移动站100根据m历元的状态向量预测值
Figure GDA00028222609400001110
卡尔曼增益矩阵km及状态向量更新方程
Figure GDA00028222609400001111
计算m历元的状态向量滤波值，并将m历元的状态向量滤波值作为移动站100在当前历元的位置；其中，
Figure GDA00028222609400001112
表示m历元的状态向量滤波值，vm表示观测噪声。
移动站100根据m历元的状态向量方差的预测值
Figure GDA00028222609400001113
卡尔曼增益矩阵km、测量系数矩阵Hm及方差更新方程
Figure GDA00028222609400001114
计算m历元的状态向量方差的滤波值；其中，
Figure GDA00028222609400001118
表示m历元的状态向量方差的滤波值，
Figure GDA00028222609400001115
表示m历元的状态向量方差的滤波值的逆矩阵，
Figure GDA00028222609400001116
表示m历元的状态向量方差的预测值
Figure GDA00028222609400001117
的逆矩阵。
进一步地，在本实施例中，当滤波器没有发散时，当某颗卫星出现多路径效应，其对应的验后残差也会出现相应的波动，因此，本申请计算出移动站100在当前历元的位置后，对每一颗卫星的验后残差进行统计监测，在下一历元的卡尔曼滤波算法中进行相关的处理，进一步消除了多路径效应的影响。下面，对验后残差的检测及处理过程进行详细说明。
请参照图6，在步骤S102之后，该动态定位方法还包括：
步骤S103，根据移动站在当前历元的位置、每颗卫星在当前历元的位置以及每颗卫星在当前历元的观测值计算每颗卫星在当前历元的验后残差。
可选地，该步骤S103包括：根据移动站在当前历元的位置和每颗卫星在当前历元的位置计算移动站与每颗卫星之间的距离，并计算距离与每颗卫星的观测值之间的差值，从而得到每颗卫星在当前历元的验后残差。
步骤S104，当任一颗目标卫星在当前历元的验后残差大于预设值时，则对目标卫星打上不可用标记，其中，不可用标记用于指示目标卫星的观测值不参与计算移动站在下一历元的位置。
在本实施例中，观测值包括伪距观测值和载波相位观测值，伪距观测值对应的预设值与载波相位观测值对应的预设值不同。例如，伪距观测值对应的预设值可以为3m，载波相位观测值对应的预设值可以为3dm，移动站100在进行载波相位观测值对应的验后残差处理时，采用3dm作为上述的预设值，在进行伪距观测值对应的验后残差处理时，采用3m作为上述的预设值。
步骤S105，在得到移动站在下一历元的位置后，根据移动站在下一历元的位置、目标卫星在下一历元的位置以及目标卫星在下一历元的观测值计算目标卫星在下一历元的验后残差，若目标卫星在下一历元的验后残差仍大于预设值，则保持目标卫星的不可用标记不变；若目标卫星在下一历元的验后残差小于或等于预设值，则将目标卫星的不可用标记删除。
例如，移动站100对每一颗卫星的验后残差进行统计监测，当某一历元时刻某颗卫星对应的伪距验后残差大于预设值3m，则将该卫星打上伪距不可用的标记，在下一历元解算位置时，该卫星的伪距观测值不进行卡尔曼滤波的测量更新计算，等到下一历元的位置解算完成后，重新计算该被打上伪距不可用的标记的卫星的伪距验后残差，如果仍大于预设值3m，则保持该伪距不可用的标记不变，如果小于或等于预设值3m，则将该卫星的伪距不可用的标记去除；载波相位验后残差的检测及处理过程与伪距验后残差相同，只是采用的预设值由3m改为3dm，此处不再赘述。
需要说明的是，在本实施例中，为了提高移动站100的位置计算精度，可以对获取的载波相位观测值进行周跳检测以及模糊度固定，通过周跳检测，可使载波相位观测值保持连续性，而通过Lambda算法对载波相位观测值的整周模糊度进行固定，从而有效提高了定位结果的计算精确性。
请参照图6，为本发明实施例所提供的动态定位装置400的功能模块示意图。需要说明的是，本实施例所提供的动态定位装置400，其基本原理及产生的技术效果与前述方法实施例相同，为简要描述，本实施例中未提及部分，可参考前述方法实施例中的相应内容。动态定位装置400包括信噪比方差计算模块410、位置计算模块420以及验后残差处理模块430。
该信噪比方差计算模块410用于计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，卫星为基准站200和移动站100同时观测到的且高度角不低于预设角度的卫星。
可以理解，该信噪比方差计算模块410可以执行上述步骤S101。
该位置计算模块420用于依据信噪比方差值以及预先建立的卡尔曼滤波模型计算移动站100在当前历元的位置。
可选地，该位置计算模块420具体用于依据信噪比方差值计算每颗卫星的信噪比方差系数，并依据信噪比方差系数以及卡尔曼滤波模型计算移动站100在当前历元的位置。
其中，该卡尔曼滤波模型包括状态预测方程和测量更新方程，该位置计算模块420用于依据状态预测方程计算当前历元的状态向量预测值和当前历元的状态向量方差的预测值，依据信噪比方差系数、当前历元的状态向量预测值、当前历元的状态向量方差的预测值以及测量更新方程计算移动站100在当前历元的位置。
可选地，该状态预测方程包括第一预测方程
Figure GDA0002822260940000141
和第二预测方程
Figure GDA0002822260940000142
该位置计算模块420具体用于根据第一预测方程
Figure GDA0002822260940000143
计算m历元的状态向量预测值，并将m历元的状态向量预测值作为当前历元的状态向量预测值；其中，
Figure GDA0002822260940000144
表示m历元的状态向量预测值，
Figure GDA0002822260940000145
表示m-1历元的状态向量滤波值，
Figure GDA0002822260940000146
表示状态转移矩阵，I表示单位矩阵，ΔT表示历元时间间隔。
该位置计算模块420具体用于根据第二预测方程
Figure GDA0002822260940000147
计算m历元的状态向量方差的预测值，并将m历元的状态向量方差的预测值作为当前历元的状态向量方差的预测值；其中，
Figure GDA0002822260940000148
表示m历元的状态向量方差的预测值，
Figure GDA0002822260940000149
表示m-1历元的状态向量方差的滤波值，
Figure GDA0002822260940000151
表示噪声系数矩阵，Qm-1表示m-1历元的过程噪声。
可选地，该测量更新方程包括增益方程
Figure GDA0002822260940000152
状态向量更新方程
Figure GDA0002822260940000153
以及方差更新方程
Figure GDA0002822260940000154
该位置计算模块420具体用于根据m历元的状态向量方差的预测值及增益方程
Figure GDA0002822260940000155
计算卡尔曼增益矩阵；其中，km表示卡尔曼增益矩阵，Hm表示测量系数矩阵且测量系数矩阵是根据获取的伪距观测值和载波相位观测值计算得到，Rm表示观测噪声方差且观测噪声方差是根据信噪比方差系数计算得到。
该位置计算模块420具体用于根据m历元的状态向量预测值、卡尔曼增益矩阵及状态向量更新方程
Figure GDA0002822260940000156
计算m历元的状态向量滤波值，并将m历元的状态向量滤波值作为移动站100在当前历元的位置；其中，
Figure GDA0002822260940000157
表示m历元的状态向量滤波值，vm表示观测噪声。
该位置计算模块420还用于根据m历元的状态向量方差的预测值、卡尔曼增益矩阵、测量系数矩阵及方差更新方程
Figure GDA0002822260940000158
计算m历元的状态向量方差的滤波值；其中，
Figure GDA0002822260940000159
表示m历元的状态向量方差的滤波值。
可以理解，该位置计算模块420可以执行上述步骤S102。
该验后残差处理模块430用于根据移动站100在当前历元的位置、每颗卫星在当前历元的位置以及每颗卫星在当前历元的观测值计算每颗卫星在当前历元的验后残差，当任一颗目标卫星在当前历元的验后残差大于预设值时，则对目标卫星打上不可用标记，其中，不可用标记用于指示目标卫星的观测值不参与计算移动站100在下一历元的位置。
可选地，观测值包括伪距观测值和载波相位观测值，其中，伪距观测值对应的预设值与载波相位观测值对应的预设值不同。
该验后残差处理模块430具体用于根据移动站100在当前历元的位置和每颗卫星在当前历元的位置计算移动站100与每颗卫星之间的距离，并计算距离与每颗卫星的观测值之间的差值，从而得到每颗卫星在当前历元的验后残差。
可选地，该验后残差处理模块430还用于在得到移动站100在下一历元的位置后，根据移动站100在下一历元的位置、目标卫星在下一历元的位置以及目标卫星在下一历元的观测值计算目标卫星在下一历元的验后残差，若目标卫星在下一历元的验后残差仍大于预设值，则保持目标卫星的不可用标记不变；若目标卫星在下一历元的验后残差小于或等于预设值，则将目标卫星的不可用标记删除。
可以理解，该验后残差处理模块430可以执行上述步骤S103、S104以及S105。
综上，本发明实施例提供的动态定位方法及装置，通过计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星；依据信噪比方差值以及预先建立的卡尔曼滤波模型计算移动站在当前历元的位置。在解算出移动站的当前历元的位置后，根据移动站在当前历元的位置、每颗卫星在当前历元的位置以及每颗卫星在当前历元的观测值计算每颗卫星在当前历元的验后残差，当任一颗目标卫星在当前历元的验后残差大于预设值时，则对目标卫星打上不可用标记，其中，不可用标记用于指示目标卫星的观测值不参与计算移动站在下一历元的位置，在得到移动站在下一历元的位置后，根据移动站在下一历元的位置、目标卫星在下一历元的位置以及目标卫星在下一历元的观测值计算目标卫星在下一历元的验后残差，若目标卫星在下一历元的验后残差仍大于预设值，则保持目标卫星的不可用标记不变；若目标卫星在下一历元的验后残差小于或等于预设值，则将目标卫星的不可用标记删除。由于低仰角卫星更易受多路径的影响，所以在计算信噪比方差值前先剔除高度角低于预设角度的卫星，有效降低多路径效应的影响；而将信噪比方差值参与到卡尔曼滤波的计算过程，能够更好地对多路径误差进行滤波和消除，从而提高RTK定位结果的准确性；通过对验后残差的检测和处理，对于验后残差超过预设值的卫星，为其打上标记，被标记的卫星在卡尔曼滤波测量更新中进行降权处理。如此，本申请针对多路径影响，在信噪比方差值计算、卡尔曼滤波以及验后残差检测三个方面进行有效的检测和滤波，有效解决了RTK定位受多路径的影响而导致的定位误差问题，实现了有效的多路径抑制效果。
需要说明的是，在本文中，诸如“第一”和“第二”等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。应注意到：相似的标号和字母在下面的附图中表示类似项，因此，一旦某一项在一个附图中被定义，则在随后的附图中不需要对其进行进一步定义和解释。

Claims (8)

1.一种动态定位方法，其特征在于，所述方法包括：
计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，所述卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星；
获取所有卫星对应的信噪比方差值中的最小值，并将每颗卫星的信噪比方差值与所述最小值之间的比值进行平方，得到每颗卫星的信噪比方差系数；
依据所述信噪比方差系数以及预先建立的卡尔曼滤波模型计算所述移动站在当前历元的位置；其中，所述卡尔曼滤波模型包括状态预测方程和测量更新方程，所述依据所述信噪比方差系数以及所述卡尔曼滤波模型计算所述移动站在当前历元的位置的步骤包括：依据所述状态预测方程计算当前历元的状态向量预测值和当前历元的状态向量方差的预测值；依据所述信噪比方差系数、所述当前历元的状态向量预测值、所述当前历元的状态向量方差的预测值以及所述测量更新方程计算所述移动站在当前历元的位置；所述信噪比方差系数用于计算所述卡尔曼滤波模型中的观测噪声方差。
2.如权利要求1所述的动态定位方法，其特征在于，所述状态预测方程包括第一预测方程
Figure FDA0003112711970000011
和第二预测方程
Figure FDA0003112711970000012
所述依据所述状态预测方程计算当前历元的状态向量预测值和当前历元的状态向量方差的预测值的步骤包括：
根据所述第一预测方程
Figure FDA0003112711970000013
计算m历元的状态向量预测值，并将所述m历元的状态向量预测值作为所述当前历元的状态向量预测值；其中，
Figure FDA0003112711970000014
表示所述m历元的状态向量预测值，
Figure FDA0003112711970000015
表示m-1历元的状态向量滤波值，
Figure FDA0003112711970000021
表示状态转移矩阵，I表示单位矩阵，ΔT表示历元时间间隔；
根据所述第二预测方程
Figure FDA0003112711970000022
计算m历元的状态向量方差的预测值，并将所述m历元的状态向量方差的预测值作为所述当前历元的状态向量方差的预测值；其中，
Figure FDA0003112711970000023
表示所述m历元的状态向量方差的预测值，
Figure FDA0003112711970000024
表示m-1历元的状态向量方差的滤波值，
Figure FDA0003112711970000025
表示噪声系数矩阵，Qm-1表示m-1历元的过程噪声。
3.如权利要求2所述的动态定位方法，其特征在于，所述测量更新方程包括增益方程
Figure FDA0003112711970000026
及状态向量更新方程
Figure FDA0003112711970000027
所述依据所述信噪比方差系数、所述当前历元的状态向量预测值、所述当前历元的状态向量方差的预测值以及所述测量更新方程计算所述移动站在当前历元的位置的步骤包括：
根据所述m历元的状态向量方差的预测值及所述增益方程
Figure FDA0003112711970000028
计算卡尔曼增益矩阵；其中，km表示所述卡尔曼增益矩阵，Hm表示测量系数矩阵且所述测量系数矩阵是根据获取的伪距观测值和载波相位观测值计算得到，Rm表示观测噪声方差，
Figure FDA0003112711970000029
表示所述观测噪声方差的逆矩阵；
根据所述m历元的状态向量预测值、所述卡尔曼增益矩阵及所述状态向量更新方程
Figure FDA0003112711970000031
计算m历元的状态向量滤波值，并将所述m历元的状态向量滤波值作为所述移动站在当前历元的位置；其中，
Figure FDA0003112711970000032
表示所述m历元的状态向量滤波值，vm表示观测噪声。
4.如权利要求3所述的动态定位方法，其特征在于，所述测量更新方程还包括方差更新方程
Figure FDA0003112711970000033
所述依据所述信噪比方差系数、所述当前历元的状态向量预测值、所述当前历元的状态向量方差的预测值以及所述测量更新方程计算所述移动站在当前历元的位置的步骤之后，所述方法还包括：
根据所述m历元的状态向量方差的预测值、所述卡尔曼增益矩阵、所述测量系数矩阵及所述方差更新方程
Figure FDA0003112711970000034
计算m历元的状态向量方差的滤波值；其中，
Figure FDA0003112711970000035
表示所述m历元的状态向量方差的滤波值。
5.如权利要求1所述的动态定位方法，其特征在于，所述依据所述信噪比方差值以及预先建立的卡尔曼滤波模型计算所述移动站在当前历元的位置的步骤之后，所述方法还包括：
根据所述移动站在当前历元的位置、每颗卫星在当前历元的位置以及所述每颗卫星在当前历元的观测值计算所述每颗卫星在当前历元的验后残差；
当任一颗目标卫星在当前历元的验后残差大于预设值时，则对所述目标卫星打上不可用标记，其中，所述不可用标记用于指示所述目标卫星的观测值不参与计算所述移动站在下一历元的位置。
6.如权利要求5所述的动态定位方法，其特征在于，所述当任一颗目标卫星在当前历元的验后残差大于预设值时，则对所述目标卫星打上不可用标记的步骤之后，所述方法还包括：
在得到所述移动站在下一历元的位置后，根据所述移动站在下一历元的位置、所述目标卫星在下一历元的位置以及所述目标卫星在下一历元的观测值计算所述目标卫星在下一历元的验后残差，若所述目标卫星在下一历元的验后残差仍大于所述预设值，则保持所述目标卫星的不可用标记不变；若所述目标卫星在下一历元的验后残差小于或等于所述预设值，则将所述目标卫星的不可用标记删除。
7.如权利要求5所述的动态定位方法，其特征在于，所述根据所述移动站在当前历元的位置、每颗卫星在当前历元的位置以及所述每颗卫星在当前历元的观测值计算所述每颗卫星在当前历元的验后残差的步骤包括：
根据所述移动站在当前历元的位置和每颗卫星在当前历元的位置计算所述移动站与所述每颗卫星之间的距离，并计算所述距离与所述每颗卫星的观测值之间的差值，从而得到所述每颗卫星在当前历元的验后残差。
8.一种动态定位装置，其特征在于，所述装置包括：
信噪比方差计算模块，用于计算每颗卫星在预设滑动窗口内的信噪比方差值；其中，所述卫星为基准站和移动站同时观测到的且高度角不低于预设角度的卫星；
位置计算模块，用于获取所有卫星对应的信噪比方差值中的最小值，并将每颗卫星的信噪比方差值与所述最小值之间的比值进行平方，得到每颗卫星的信噪比方差系数；依据所述信噪比方差系数以及预先建立的卡尔曼滤波模型计算所述移动站在当前历元的位置；其中，所述卡尔曼滤波模型包括状态预测方程和测量更新方程，所述位置计算模块用于依据所述状态预测方程计算当前历元的状态向量预测值和当前历元的状态向量方差的预测值；依据所述信噪比方差系数、所述当前历元的状态向量预测值、所述当前历元的状态向量方差的预测值以及所述测量更新方程计算所述移动站在当前历元的位置；所述信噪比方差系数用于计算所述卡尔曼滤波模型中的观测噪声方差。
CN201910356541.0A 2019-04-29 2019-04-29 动态定位方法及装置 Active CN110058281B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910356541.0A CN110058281B (zh) 2019-04-29 2019-04-29 动态定位方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910356541.0A CN110058281B (zh) 2019-04-29 2019-04-29 动态定位方法及装置

Publications (2)

Publication Number Publication Date
CN110058281A CN110058281A (zh) 2019-07-26
CN110058281B true CN110058281B (zh) 2021-09-17

Family

ID=67321703

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910356541.0A Active CN110058281B (zh) 2019-04-29 2019-04-29 动态定位方法及装置

Country Status (1)

Country Link
CN (1) CN110058281B (zh)

Families Citing this family (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111045052B (zh) * 2019-10-14 2022-08-09 广东星舆科技有限公司 一种智能终端伪距差分定位及质量控制方法
CN110927757B (zh) * 2019-12-26 2022-03-04 广东星舆科技有限公司 卫星观测数据的质控方法、装置及定位装置
CN112904390B (zh) * 2021-01-14 2022-03-25 腾讯科技（深圳）有限公司 定位方法、装置、计算机设备和存储介质
CN113670953B (zh) * 2021-08-20 2022-07-22 北京航空航天大学 一种基于加权窗求直射与反射信号snr干涉峰值的方法
CN114035216A (zh) * 2021-11-10 2022-02-11 北京百度网讯科技有限公司 定位方法、装置、设备以及存储介质
CN115598686A (zh) * 2022-10-19 2023-01-13 中国民航大学（Cn） 一种基于伪卫星信号的测向方法

Family Cites Families (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CL2009000010A1 (es) * 2008-01-08 2010-05-07 Ezymine Pty Ltd Metodo para determinar la posicion global de una pala minera electrica.
CN101655546B (zh) * 2009-09-22 2012-02-01 西北工业大学 卫星导航直接信号参数与多径信号数目和参数的估计方法
RU2476904C1 (ru) * 2011-08-30 2013-02-27 Открытое акционерное общество "Конструкторское бюро приборостроения" Способ сопровождения объекта и способ формирования сигнала управления положением луча приемно-передающего антенного устройства системы сопровождения объекта
CN104931981A (zh) * 2015-05-11 2015-09-23 中国科学院光电研究院 一种基于信噪比波动信息的gnss抗多径选星方法
CN107534530B (zh) * 2015-09-25 2020-07-17 诸暨市尚诺五金经营部 信干噪比的计算方法、装置和接收机
CN106231553B (zh) * 2016-08-03 2019-03-29 西北工业大学 基于无线声传感器网络的多节点信息融合声源定位方法
CN106526634B (zh) * 2016-10-19 2019-02-26 闽江学院 一种基于自调整卡尔曼滤波的多普勒辅助载波相位平滑伪距方法
CN107748377B (zh) * 2017-12-06 2021-09-17 湖南国科微电子股份有限公司 基于gnss的差分定位方法及其定位系统
CN108614284B (zh) * 2018-04-24 2022-04-05 北京邮电大学 一种定位信号处理方法、装置及设备
CN108845340A (zh) * 2018-06-01 2018-11-20 浙江亚特电器有限公司 基于gnss-rtk的定位方法
CN109521443B (zh) * 2018-12-29 2021-04-23 广东电网有限责任公司 一种探测星历异常的方法

Also Published As

Publication number Publication date
CN110058281A (zh) 2019-07-26

Similar Documents

Publication Publication Date Title
CN110058281B (zh) 动态定位方法及装置
CN109521443B (zh) 一种探测星历异常的方法
Zangeneh-Nejad et al. Cycle slip detection and repair of undifferenced single-frequency GPS carrier phase observations
JP2003185728A (ja) キャリア位相相対測位装置
CN107367744B (zh) 基于自适应测量噪声方差估计的星载gps定轨方法
CN113138402B (zh) 基于rtk的模糊度固定方法及装置、存储介质
AU2004273434A1 (en) Method for receiver autonomous integrity monitoring and fault detection and elimination
CN112394383B (zh) 一种卫星与5g基站组合定位方法及装置
CN108614284B (zh) 一种定位信号处理方法、装置及设备
CN111077550A (zh) 一种应用于智能终端rtd定位的粗差探测方法及系统
US7299058B2 (en) Position calculation method and storage medium storing a program therefor
CN114488230B (zh) 一种多普勒定位方法、装置、电子设备和存储介质
CN114966762A (zh) 一种gnss/ins紧组合导航系统的故障检测方法
CN114662059A (zh) 一种海上卫星大地坐标的高程拟合方法与装置
CN116840876A (zh) 适用于相位失真的卫星导航天线的双天线定向方法及装置
CN111522032A (zh) 一种北斗三号系统用户完好性处理的优化方法及优化装置
CN111665529A (zh) 伪距观测数据的处理方法、装置及可读介质
CN117724125B (zh) 一种基于一致性的观测数据的质量控制方法及装置
Elsayed et al. Bounding of correlated double-differenced GNSS observation errors using NRTK for precise positioning of autonomous vehicles
CN111679307B (zh) 一种卫星定位信号解算方法及装置
CN110174691B (zh) 一种定位设备、方法及计算机可读存储介质
US7868822B2 (en) Method and apparatus for determining an integrity risk in a satellite reference system
CN114355393A (zh) 一种基于低成本接收机的三天线姿态估计方法
CN115327584B (zh) 综合多指标的gnss观测值定权方法、终端及介质
US11047992B2 (en) Positioning method and positioning terminal

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