CN111025336B - 一种多星非高斯载波多普勒异常检测及剔除方法 - Google Patents

Abstract

一种实时多星非高斯载波多普勒异常检测及剔除(Fault Detection and Elimination，FDE)方法，属于低轨GNSS导航增强精密测速技术领域，通过最近邻低轨用户钟速约束实现非高斯分布下载波多普勒高完好性监测和故障剔除，可实现多颗卫星切换引起的载波多普勒异常跳变的检测和剔除，从而克服了多颗卫星切换引入的用户钟速序列的异常跳变问题，提高了采用GNSS对低轨卫星时钟驯服时用户钟速的估计精度和完好性，具有广阔的应用前景。

Description

一种多星非高斯载波多普勒异常检测及剔除方法

技术领域

本发明涉及一种多星非高斯载波多普勒异常检测及剔除方法，属于低轨GNSS导航增强精密测速技术领域。

背景技术

随着低轨GNSS导航增强卫星系统的发展，基于载波多普勒的时钟驯服技术在低轨卫星频率基准建立和维持方面具有巨大优势和潜力，高稳晶振经过驯服后的长期和短期稳定度性能可以媲美替代星载原子钟，具有重要的应用价值和市场竞争力。该技术主要通过GNSS系统高精度载波相位测量值计算出高精度载波多普勒(精度可达到1～3mm/s量级)，进而以此高精度的载波多普勒可估计出接收机自身的钟速，再根据三阶锁频环对接收机的本振频率进行闭环调整，使其锁定到其标称频率上，实现卫星频率基准的高精度建立和维持。整个过程对载波多普勒的完好性提出了较高的要求，但是，可见卫星的切换的因素对毫米每秒量级的载波多普勒测量值会造成毫米量级的抖动，其对载波多普勒测量的完好性造成了严重挑战。由于载波相位测量值的误差特性是非高斯分布的，其差分后的载波多普勒的误差特性也是非高斯分布，传统的基于卡方分布的接收机自主完好性监测RAIM(ReceiverAutonomous Integrity Monitoring)方法不适用。由于卫星切换导致的载波多普勒跳变的完好性监测，目前尚无好的解决方法。

针对载波多普勒高完好性监测存在的以上问题，世界上已经有学者对载波相位RAIM算法做了相关理论研究，其具有代表性的是韩国Ho Yun的基于高斯混合模型(Gaussian Mixture Model)的高斯和滤波器(Gaussian Sum Filter)载波相位RAIM方法。该方法针对载波相位非高斯分布的问题，采用多个高斯分布叠加后的分布来逼近载波相位的非高斯误差分布，并给出了贝叶斯框架下的载波相位完好性监测方法。该方法计算复杂度高，仅适合于后处理，对工程实践中的实时应用理论指导意义不大。此外，该方法不适用载波多普勒的完好性监测问题，不能解决卫星切换导致的载波多普勒跳变完好性监测问题。

发明内容

本发明要解决的技术问题是：克服现有技术的不足，提供了一种实时多星非高斯载波多普勒异常检测及剔除(Fault Detection and Elimination，FDE)方法，通过最近邻低轨用户钟速约束实现非高斯分布下载波多普勒高完好性监测和故障剔除，可实现多颗卫星切换引起的载波多普勒异常跳变的检测和剔除，从而克服了多颗卫星切换引入的用户钟速序列的异常跳变问题，提高了采用GNSS对低轨卫星时钟驯服时用户钟速的估计精度和完好性，具有广阔的应用前景。

本发明目的通过以下技术方案予以实现：

一种多星非高斯载波多普勒异常检测及剔除方法，包括如下步骤：

S1、根据三个历元的载波相位计算出中间历元的载波多普勒；并将中间历元作为当前历元；

S2、利用当前历元的载波多普勒，解算当前历元用户三维速度和钟速；

S3、利用当前历元用户三维速度和钟速，根据三个历元中的第一个历元对应的L点钟速均值，进行当前历元的最近邻钟速约束的载波多普勒异常检测及剔除；然后获得当前历元的最优载波多普勒RAIM后的钟速结果，和，当前历元的L点钟速均值；然后进行下一个历元的检测及剔除。

优选的，低轨卫星GNSS接收机第n个历元、第f个频点、第s颗GNSS卫星的载波相位

为：

式中，

是卫星天线和接收机天线相位中心之间的几何距离，c是光速,dT^s和dt_r分别是卫星和接收机用户的钟差，

是电离层误差，λ_f和

分别是波长和模糊度，

代表载波相位噪声。

优选的，第n个历元载波多普勒

的计算方法为：

式中，

为低轨卫星GNSS接收机第n个历元、第f个频点、第s颗GNSS卫星的载波相位；T_s为两个相邻历元间的采样间隔

优选的，所述载波多普勒异常检测方法为：

第n个历元，m颗卫星是否存在一个或多个卫星的载波多普勒异常通过将检测统计量T与阈值Thd进行比较获得，若T＞Thd，则第n个历元的m颗卫星中一颗到多颗卫星的载波多普勒测量值存在完好性异常；若T≤Thd，则第n个历元的m颗卫星的载波多普勒测量值均无完好性异常。

优选的，所述多星非高斯载波多普勒异常检测及剔除方法的平滑长度不超过1周。

一种计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时，实现上述多星非高斯载波多普勒异常检测及剔除方法的步骤。

本发明相比于现有技术具有如下有益效果：

(1)本发明采用最近邻钟速约束模型，可对非高斯分布下的多颗卫星载波多普勒异常变化进行检测和剔除，可有效抑制多颗卫星同时切换导致的用户钟速抖动，提升用户钟速估计的中长期稳定度；

(2)本发明方法将非高斯载波多普勒噪声的RAIM问题分解为高斯噪声拟合部分和残余非高斯分布两部分，并分别给出了相应的处理技术方法；

(3)本发明方法采用组合优化的方式处理残余非高斯分布噪声下的多颗GNSS卫星优选问题，算法复杂度低、通用性强，利于工程实际应用；

(4)本发明方法采用最近邻低轨用户钟速约束技术，降低GNSS卫星切换对低轨用户钟速引入的钟速跳变，提高了低轨用户钟速估计的准确度和稳定度；

(5)本发明方法提出一种矢量检验统计量，可一次性检测出多颗精密测速残差较大的卫星，检测效率高。

附图说明

图1为多星非高斯载波多普勒异常检测及剔除方法流程图；

图2为载波多普勒完好性FDE的用户钟速及可见卫星数；

图3为检测剔除一颗故障卫星的用户钟速；

图4为检测剔除多颗故障卫星的用户钟速；

图5(a)为未采取载波多普勒完好性FDE的用户钟速稳定度；

图5(b)为多星载波多普勒完好性FDE后的用户钟速稳定度；

图5(c)为多颗卫星载波多普勒FDE对用户时钟稳定度的改善幅度。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明的实施方式作进一步详细描述。

常规的接收机RAIM方法，假设同一个时刻仅仅有一颗卫星的测量发生完好性故障，在载波多普勒服从高斯分布假设下，用户钟速采用相应的卡方分布的标量检验统计量，若该检验统计量超过阈值门限，每次从观测卫星列表中剔除一颗卫星，重新计算用户钟速，最后挑选出最小的用户钟速作为最优用户钟速。然而，在载波多普勒完好性监测过程中发现以上处理策略效果并不理想，因为同一个历元内可能存在两颗卫星同时切换。本发明提出的载波多普勒完好性异常检测和剔除方法可以同时检测并剔除多个卫星的载波多普勒的异常测量值。在第n个历元，载波多普勒完好性H0和H1假设分别可描述为：

H0:第n个历元的m颗卫星的载波多普勒测量值均无完好性异常；

H1:第n个历元的m颗卫星中一颗到多颗卫星的载波多普勒测量值存在完好性异常。

一种多星非高斯载波多普勒异常检测及剔除方法，其工作流程如图1所示。多星非高斯载波多普勒异常检测及剔除主要包括以下五个部分：(1)根据三个历元的载波相位计算出中间历元的载波多普勒；(2)解算中间历元用户三维速度和钟速；(3)根据上个历元对应的L点钟速均值，进行中间历元的最近邻钟速约束的载波多普勒异常检测及剔除；(4)输出最优载波多普勒RAIM后的钟速结果；(5)更新中间历元L点钟速均值。

该方法将多颗GNSS卫星非高斯分布的载波多普勒异常跳变故障检测及剔除分为以下两步处理：第一步，将所有可观测GNSS卫星的载波多普勒的残差按高斯分布假设，构建矢量检验统计量，检测该矢量检验统计量是否超过高斯分布下的理论门限，给出备选载波多普勒异常卫星列表；第二步，先在备选载波多普勒异常卫星列表中挑选一颗、二颗…多颗同时失效卫星的组合，再将其从当前历元的观测卫星列表中剔除，并分别计算相应的用户钟速，最后将这些备选用户钟速与前L点平均用户钟速之差的绝对值进行比较，采用均值距离最小准则，选择出最优的卫星剔除方案，完成对载波多普勒非高斯分布残差部分的自主完好性检测与剔除。本发明提出的多星非高斯载波多普勒异常检测及剔除工作流程如图1所示，其中L点钟速均值计算中平滑长度L一般不超过1周，本实施例选L＝604800s。

具体的，一种多星非高斯载波多普勒异常检测及剔除方法，包括如下步骤：

1)获取单历元载波多普勒测量值

载波多普勒接收机自主完好性检测的输入是由三个历元中首尾两个历元载波相位的差分计算得到的中间历元的载波多普勒。低轨卫星GNSS接收机第n个历元、第f个频点、第s颗GNSS卫星的载波相位

可以表示为：

其中

是卫星天线和接收机天线相位中心之间的几何距离，c是光速，dT^s和dt_r分别是卫星和接收机用户的钟差，

是电离层误差，λ_f和

分别是波长和模糊度，

代表载波相位噪声。

本发明提出的第n个历元载波多普勒

计算方法如下式所示

其中T_s表示两个相邻历元间的采样间隔。由于

和

中的误差是非高斯分布的，故

中的误差也是非高斯分布的，其接收机自主完好性监测不能套用伪距RAIM算法。

2)解算使用全部卫星载波多普勒时的用户速度和钟速

当同时可观测4颗以上卫星时，用户三维速度v_r(n)和用户钟速

组成的未知矢量

可以通过迭代求解下式(3)获取

b＝Hx+ε_v (3)

其中

s＝1,2,...,m是第n个历元同时观测到m颗GNSS卫星的三维速度。ε_v为噪声残差矢量

为第m颗GNSS卫星的卫星钟速，

为用户与第m颗卫星的星地单位矢量。

卫星s对应的星地单位矢量可表示为：

其中r^s(n)＝[x^s(n),y^s(n),z^s(n)^T]和r_r(n)＝[x_r(n),y_r(n),z_r(n)]^T分别是卫星和接收机用户的三维位置。Sagnac效应的变化率

可以表示为：

其中ω_e是地球转动角速度，v_x,r(n)、v_y,r(n)是用户三维速度v_r(n)中的x、y方向上的分量。

为卫星s在y方向的速度分量，x_r(n)为用户在x方向的位置分量，y^s(n)为卫星s在y方向的位置分量，

为卫星s在x方向的速度分量，x^s(n)为卫星s在x方向的位置分量。

在户三维速度和用户钟速组成的未知矢量x求解迭代收敛后，令使用全部GNSS卫星载波多普勒时，计算得到的用户速度和钟速矢量为x₀。

3)计算高斯分布假设下备选故障卫星列表

对x₀对应的H进行正交三角QR分解得到

H＝QR (8)

其中R是上三角矩阵，Q是正交矩阵，根据Q^-1＝Q^T，将Q^T按下式进行分块

令N＝N_(5:m,:)，用N左乘x求解迭代收敛后的残差b，得到p＝Nb。实时计算出p的标准差σ_d＝p^Tp/(m-4)。按下式组成矢量检测统计量T，并给出高斯分布下载波多普勒异常卫星检测的阈值。

Thd＝μ·std(T),μ∈{1,2,3} (11)

其中.×表示元素相乘，sum(A,1)表示对A的按列求和形成的行向量，μ是阈值乘数因子，std(T)为……计算矢量检测统计量T的标准差p^T为p的转置。第n个历元，m颗卫星是否存在一个或多个卫星的载波多普勒异常通过将检测统计量T与阈值Thd进行比较获得，判决规则为：若T＞Thd，则判决H1假设成立；若T≤Thd，则判H0假设成立。若H1假设成立，将所有T＞Thd的卫星号放入一个集合Ω，其元素个数记为df。

4)构成最优的用户三维速度和用户钟速的估计值的备选集合

记第n时刻的可见卫星列表为Υ，对Ω中元素按照

进行排列分别组成去掉1～df个故障卫星的集合Ω₁、Ω₂、…Ω_df，然后将这些集合中的单个元素从Υ中依次剔除，得到如下的一系列用户三维速度和用户钟速的估计子集合。

……

其中card(·)表示集合的元素。最优的用户三维速度和用户钟速的估计值的备选集合可表示为

X_ensemble＝X₁∪X₂∪...∪X_df∪{x₀} (15)

5)剔除多颗非高斯分布误差引起的故障卫星

记最优的用户三维速度和用户钟速组成的矢量为

并令其初始值为x^*(0)＝x₀。按照下式计算用户钟速的历元平均值

根据第n个历元备选集合X_ensemble中用户钟速与第n-1个历元用户钟速平均值的绝对值差异最小，得到其在X_ensemble元素下标，选择出最优x^*(n)

x^*(n)＝X_ensemble{i} (18)

其中X_ensemble{i}₄表示X_ensemble中第i个元素中的第4个子元素，即用户钟速分量。通过以上途径，可实现从第n个历元观测卫星集合Υ中剔除多颗由于非高斯分布引起的故障卫星。

为了进一步验证本发明所提多星非高斯载波多普勒异常检测及剔除方法的科学性、正确性、有效性，同时分析评估其对接收机钟速估计所能达到的性能指标。本实施例通过Trimble Net R9接收机接收Spirent模拟器产生的GPS射频信号数据，采用氢钟为Spirent模拟器和Trimble接收机提供同源的10MHz参考频率信号。由于低轨接收机采用双频无电离层组合消除99.9％的电离层延迟误差，故Spirent模拟GPS信号时关闭了电离层和对流层延迟。

未采用载波多普勒异常跳变和剔除手段，利用观测时段内载波多普勒计算出的用户钟速以及可见卫星的变化规律如图2所示。

图2中用户钟速的标准差为2.3mm。从图2中可以看出，未经载波多普勒FDE的用户钟速变化起伏较明显，用这样的钟速估计值对用户时钟驯服将直接降低驯服后时钟的准确度和稳定度。理论上由于Trimble接收机和Spirent模拟器的频率是同源在氢钟上的，故接收机的钟速理论上应该稳定在零附近。此外，整个观测过程中卫星数目切换的较频繁，且卫星切换与用户钟速估计值的跳变具有较明显的对应关系。故可以确定导致卫星钟速的变化的一种因素是由于卫星切换导致的。

按照传统的RAIM算法假设，在一个时刻仅仅有一颗卫星的载波多普勒观测值发生RAIM故障。首先在式(11)中令μ＝1，检测并统计出超过门限的卫星列表。将超出阈值绝对值最大的那颗卫星挑选过去，其结果如图3所示。

图3中用户钟速的标准差为3.1mm。从图3可看出，检测剔除一颗故障卫星后，用户钟速的估计曲线总体变得更加平缓，但是有些区域的用户钟速较未采用FDE的结果有恶化，其超过6E-11，故其钟速估计值跳变更大。其原因是卫星切换一个时刻可能超过一颗，仅剔除一颗卫星的测量值没有真实消除多颗切换引入的卫星钟速跳变，这从图2中可见卫星的变化曲线中可得到证实。采用本发明方法检测剔除多颗故障卫星后，用户钟速估计的曲线如图4所示。

图4中用户钟速的标准差为1.6mm。由图4可看出，经过本发明方法检测剔除多颗故障卫星后，用户钟速的精度较未载波多普勒完好性FDE的用户钟速及仅检测剔除一颗故障卫星的用户钟速精度有明显提升，用户钟速的变化在整个观测区间内也更为平稳。图5中给出了多颗卫星载波多普勒完好性FDE对用户时钟稳定度提升效果对比。

从图5(a)～图5(c)可看出，采用本发明多颗卫星载波多普勒完好性FDE方法，64s之前的稳定度有微小的恶化，恶化幅度小于5E-13；64s以后的中长期稳定度均具有较大幅度的明显改善，改善幅度在5E-13～2E-12。64s之前稳定度的微小恶化是由于本发明提出的载波多普勒RAIM算法中利用了前L点平均用户钟差的信息。

本发明说明书中未作详细描述的内容属本领域技术人员的公知技术。

本发明虽然已以较佳实施例公开如上，但其并不是用来限定本发明，任何本领域技术人员在不脱离本发明的精神和范围内，都可以利用上述揭示的方法和技术内容对本发明技术方案做出可能的变动和修改，因此，凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化及修饰，均属于本发明技术方案的保护范围。

Claims (5)

1.一种多星非高斯载波多普勒异常检测及剔除方法，其特征在于：包括如下步骤：

S1、根据三个历元的载波相位计算出中间历元的载波多普勒；并将中间历元作为当前历元；

S2、利用当前历元的载波多普勒，解算当前历元用户三维速度和钟速；

S3、利用当前历元用户三维速度和钟速，根据三个历元中的第一个历元对应的前L点钟速均值，进行当前历元的最近邻钟速约束的载波多普勒异常检测及剔除；然后获得当前历元的最优载波多普勒RAIM后的钟速结果，和，当前历元的前L点钟速均值；然后进行下一个历元的检测及剔除。

2.根据权利要求1所述的一种多星非高斯载波多普勒异常检测及剔除方法，其特征在于：低轨卫星GNSS接收机第n个历元、第f个频点、第s颗GNSS卫星的载波相位

为：

式中，

是卫星天线和接收机天线相位中心之间的几何距离，c是光速,dT^s和dt_r分别是卫星和接收机用户的钟差，

是电离层误差，λ_f和

分别是波长和模糊度，

代表载波相位噪声。

3.根据权利要求1所述的一种多星非高斯载波多普勒异常检测及剔除方法，其特征在于：第n个历元载波多普勒

的计算方法为：

式中，

为低轨卫星GNSS接收机第n个历元、第f个频点、第s颗GNSS卫星的载波相位；T_s为两个相邻历元间的采样间隔。

4.根据权利要求1所述的一种多星非高斯载波多普勒异常检测及剔除方法，其特征在于：所述多星非高斯载波多普勒异常检测及剔除方法的平滑长度不超过1周。

5.一种计算机可读存储介质，其上存储有计算机程序，其特征在于，该程序被处理器执行时，实现权利要求1～4之一所述方法的步骤。

Images