CN104020482A - 一种高动态卫星导航接收机精确测速方法 - Google Patents

Abstract

本发明公开了卫星导航接收机领域的一种高动态卫星导航接收机精确测速方法，用以解决目前卫星导航接收机领域研究中存在的问题。该方法包括：首先，提取卫星载波与本地载波的频率差值；然后，对提取的载波频率差值进行滤波处理；接着，利用最小二乘法求解接收机载体速度；最后，运行速度RAIM算法。本发明针对传统高动态卫星导航接收机测速技术中存在的卫星信号跟踪失锁、频率抖动、测量速度算法复杂、测量结果随机误差大等实际应用的局限性，创造性地提出了一种高动态卫星导航接收机精确测速方法。

Description

一种高动态卫星导航接收机精确测速方法

技术领域

本发明涉及卫星导航接收机领域，特别涉及一种用于高动态卫星导航接收机精确测速方法。

背景技术

北斗卫星导航系统(BeiDou Navigation Satellite System,BDS)是中国正在实施的自主发展、独立运行的全球卫星导航系统。北斗卫星导航系统和美国全球定位系统(Global Positioning System，GPS)、俄罗斯格洛纳斯(Global Navigation Satellite System，GLONASS)及欧盟伽利略(GALILEO)系统一起，是联合国卫星导航委员会已认定的供应商。

卫星导航接收机可以确定载体的位置、速度和时间，广泛应用于导航、测量与授时等领域。其中，测速技术需要提取卫星信号载波多普勒频率值。提取的卫星多普勒频率值会由于载波跟踪环路的热噪声而存在随机抖动；另外由于卫星运动、卫星健康状况或者信号干扰、遮挡等原因会引起卫星信号跟踪失锁、频率抖动，提取的多普勒频率会有较大误差和错误，从而引起计算的速度误差变大。现在有些导航接收机采用一些滤波算法平滑处理，但会影响接收机的动态性能，使接收机在高速、高加速度运动场景时，测速误差变大。

测速是卫星导航接收机的重要功能，是移动目标导航基本需求，尤其对于高动态飞行器导航，常用惯性导航与卫星导航接收机组合的工作模式。在这种应用场合卫星导航接收机的测速结果实时修正惯性导航的结果，消除惯性导航的积累误差，对接收机的测速精度具有较高的要求，一般要求测速误差小于0.2m/s。

发明内容

本发明的目的在于，提出一种高动态卫星导航接收机精确测速方法，用以解决目前导航系统测速研究存在的问题。

为实现上述目的，本发明提出了一种高动态卫星导航接收机精确测速方法。卫星导航接收机包括射频部分和基带部分：射频部分完成导航信号的放大、滤波及变频功能；基带部分又可分为基带信号处理部分和基带信息处理部分，信号处理部分包括对各个导航卫星信号的捕获、跟踪、解调功能，信息处理部分根据信号处理电路得到的信息完成位置、速度与时间的解算，其中基带部分包括一个嵌入式CPU来运行信息处理程序，其特征是所述方法包括下列步骤：

步骤1：通过直接读取载波跟踪环频率误差量的方法提取卫星载波与本地载波的频率差值；

步骤2：设计一种多档α-β-γ滤波器对提取的卫星载波频率差值进行滤波处理；

步骤3：将滤波后的各频率差值，组成求解接收机速度的方程组，利用最小二乘法求解出接收机速度；

步骤4：设计一种针对速度求解的接收机自主完好性监测(Receiver Autonomous Integrity Monitoring，RAIM)算法，排除速度解算方程组中的错误方程，确保解算速度不出现错误。

所述多档α-β-γ滤波器由α-β-γ滤波器模块、状态判决模块、档位切换模块、输出切换模块组成；

所述α-β-γ滤波器模块输入端口与所述档位切换模块输出端口相连，所述α-β-γ滤波器模块输出端口与所述输出切换模块输入端口相连，用于对提取的卫星载波频率差值进行滤波处理；

所述状态判决模块输入端口与所述输出切换模块相连，用于反馈滤波结果；

所述状态判决模块输出端口分别与所述档位切换模块、输出切换模块相连，用于判断接收机状态是否满足切换到其他档位进行滤波的条件。

所述使用多档α-β-γ滤波器的具体步骤：

步骤21：对卫星载波频率差值，通过第1档α-β-γ滤波器开始滤波；

步骤22：分析滤波结果及导航结算结果，判断接收机状态是否满足切换到其他档位进行滤波的条件，并将其反馈至档位切换模块；

步骤23：根据状态判决模块的反馈，在档位切换模块估算当前信息的方差，判决适合使用哪个档位的滤波器；

步骤24：将卫星载波频率差值，通过档位切换模块确定的合适档位滤波器进行滤波处理；

步骤25：循环步骤22至步骤24，直至收敛到稳定工作状态时，由输出切换模块输出滤波结果。

所述速度求解的RAIM算法的具体步骤：

步骤41：先用位置解算使用的那几颗卫星的多普勒观测量做一次速度最小二乘解算，然后开始执行速度RAIM算法；

步骤42：如果解算的卫星数小于或等于4颗，则速度RAIM算法无法检测故障星，直接输出计算的速度结果；

步骤43：若卫星数超过4颗，则根据速度最小二乘法的计算结果，算出多普勒残余误差向量ε及速度RAIM算法统计检测量T_X；

步骤44：将T_X与对应的检测门限T_D比较，若T_X≤T_D，没有故障星，RAIM算法完成，并输出定位测速结果；

步骤45：若T_X>T_D，则判断卫星数是否超过5颗。若星数未超过5颗，无法排除故障星，RAIM算法结束，同时输出测速结果不可用的标志，提醒使用者；

步骤46：若星数超过5颗，则可以找到发生故障的卫星，多普勒残余最大的方程判断为错误方程；

步骤47：将错误方程剔除后，用剩余的卫星做速度最小二乘解算，然后重复前面的步骤，直至速度RAIM算法统计检测量T_X小于检测门限T_D，或方程组的个数小于5，输出解算结果。

本发明的有益效果为：

1.找到一种适合于高动态运动物体的卫星载波多普勒提取的滤波算法：

通过此方法可以滤除热噪声引起的随机抖动，实现对多普勒频率的平滑，同时可以适应高动态移动引起的频率变化，并且该算法要尽量简洁，节约资源，便于嵌入式中央处理器(Central Processing Unit,CPU)实现。

2.对速度解算方程进行检测，排除错误方程，使解算速度不出现错误跳动：

接收机输出的解算速度出现大的错误，一般是由于定位位置错误或则提取的卫星多普勒频率错误引起，本发明设计一种针对速度的RAIM算法，排除速度解算方程组中的错误方程，使解算速度不出现错误。

3.利用本发明测速性能可以达到以下指标：

测速范围：0～3000m/s；加速度：<40g；加加速度：<10g/s；测速误差：<0.2m/s。

附图说明

图1是多档位α-β-γ滤波器结构框图。

图2是没有速度RAIM算法的测试结果。

图3是加上速度RAIM算法的测试结果。

具体实施方式

下面结合附图，对优选实施例作详细说明。应该强调的是，下述说明仅仅是示例性的，而不是为了限制本发明的范围及其应用。

本发明解决问题的思路是：首先，提取卫星载波与本地载波的频率差值；然后，对提取的载波频率差值进行滤波处理；接着，利用最小二乘法求解接收机载体速度；最后，运行速度RAIM算法。

下面结合附图说明本发明的具体实现方式。该发明包括如下的步骤：

步骤1：通过直接读取载波跟踪环频率误差量的方法提取卫星载波与本地载波的频率差值。

目前的卫星导航系统信号均为直序扩频二进制相移键控(BinaryPhase Shift Keying,BPSK)调制方式的信号，在接收机基带信号处理部分需要实现对各颗卫星信号的载波跟踪和伪码跟踪。其中在载波跟踪环可以得到本地合成载波与接收载波的频率误差。此误差包含接收机移动、卫星本身移动、地球自转和接收机钟漂等引起的误差。由载波跟踪环路读取的载波频率误差值，具有实时准确的特点，适合于高动态场景应用，本算法就是直接采用从载波跟踪环提取的方法得到卫星载波与本地载波的频率差值。通过导航电文读取卫星自身的钟差修正与地球自转引起的频差修正。

步骤2：设计一种多档α-β-γ滤波器对提取的卫星载波频率差值进行滤波处理。

由于载波环路热噪声的存在会引起抖动，为了提升测量精度，对观测量进行平滑滤波是必不可少的环节。而通过滤波器的平滑提高精度，总会以参数变化率的失真作为代价。对于高动态接收机来说，精度与变化速度同等重要，故需要寻求一种，对减小抖动和快速变化兼顾的滤波方法，同时还要保证良好的实时性和尽量少的计算量。本算法基于α-β-γ滤波器，尝试通过扩展改进，得到满足快速收敛速度与平滑精度双重要求的滤波器。

α-β-γ滤波器由α-β滤波器扩展而来，其实质上是运动方程为匀加速直线运动的卡尔曼滤波器的稳态形式。其计算形式十分简单，需要存储的数据量也很少，有着很高的实时性，滤波效果也比较理想，因此在各领域都有着广泛的应用，尤其适合应用于高动态导航接收机的观测量平滑。

(1)α-β-γ滤波器基本原理

假设变量x作匀加速变化，则状态量X是包含x及其一二阶导数的三维向量：

$X={(x,\stackrel{\&CenterDot;}{x},\stackrel{\&CenterDot;\&CenterDot;}{x})}^T---\left(1\right)$

状态方程为：

X_k+1＝ΦX_k+Γq_k    (2)

其中：

$\mathrm{\&Phi;}=\left[\begin{array}{ccc}1& T& T^2/2\\ 0& 1& T\\ 0& 0& 1\end{array}\right]$ $\mathrm{\&Gamma;}=\left[\begin{array}{c}{T}^3/6\\ T^2/2\\ T\end{array}\right]$

q_k为均值为零的高斯白噪声，T为采样间隔。

观测方程为：

y_k＝HX_k+v_k    (3)

其中，H＝(1,0,0)，v_k为测量噪声，是均值为零、标准差为的高斯白噪声。

α-β-γ滤波器实质上是运动方程为匀加速的卡尔曼滤波器的稳态形式。当卡尔曼滤波递推足够多的步数后，增益矩阵将趋于常数矩阵K，由此得到α-β-γ滤波器公式：

$\left\{\begin{array}{c}{\hat{X}}_{k+1/k}=\mathrm{\&Phi;}{\hat{X}}_k\\ {\hat{X}}_{k+1/k}={\hat{X}}_{k+1/k}+K(y_k-H{\hat{X}}_{k+1/k})\end{array}\right.---\left(4\right)$

其中， $K=\left[\begin{array}{ccc}\mathrm{\&alpha;}& \frac{\mathrm{\&beta;}}{T}& \frac{2\mathrm{\&gamma;}}{T^2}\end{array}\right],$ α、β、γ为常系数。

α-β-γ滤波器的平滑精度和收敛速度取决于α、β、γ这三个常系数的选取。选取的方法主要有：临界阻尼选择法、最佳选择法以及卡尔曼滤波稳态增益选择法等。本算法选用最佳选择法，选取规则如公式(5)。

$\left\{\begin{array}{c}\mathrm{\&beta;}=\frac{(2{\mathrm{\&alpha;}}^3-4{\mathrm{\&alpha;}}^2)+\sqrt{4{\mathrm{\&alpha;}}^6-64{\mathrm{\&alpha;}}^5+64{\mathrm{\&alpha;}}^4}}{8(1-\mathrm{\&alpha;})}\\ \mathrm{\&gamma;}=\frac{\mathrm{\&beta;}(2-\mathrm{\&alpha;})-{\mathrm{\&alpha;}}^2}{\mathrm{\&alpha;}}\end{array}\right.---\left(5\right)$

其中α的值自由选取，α越大，滤波器平滑精度越低，但收敛速度越快，反之则平滑精度高，收敛速度慢。

(2)多档位α-β-γ滤波器设计应用

由于α-β-γ滤波器是基于观测量匀加速这个假设而得到的，在实际应用中，就会遇到增大系数取值会造成平滑精度低，减小系数取值又会造成平滑收敛速度慢的两难境地。尤其对于高动态接收机而言，面向高动态应用，就意味着接收机的速度、加速度、乃至加加速度等，都有很大的变动范围，如果仅使用单一α-β-γ滤波器，必定不能满足精度和收敛速度的要求，而使用卡尔曼滤波等自适应能力强的滤波器，又会使得计算量大为增加，造成接收机计算负担过重，或者使得实时性不能得到保障。

由此很自然的想到采用多档位α-β-γ滤波器切换的方式进行滤波处理，当接收机运动状态稳定时，系统将滤波器的α、β、γ系数切换到值较小的档位上，使得测量值的平滑精度得到保障；反之，当接收机运动状态变化很剧烈时，则系统将滤波器α、β、γ系数切换到值较大的档位，使得平滑的收敛速度得到保障。

图1是多档位α-β-γ滤波器结构框图。

其中，α-β-γ滤波器组模块中存有N组α、β、γ参数值，每组参数都满足公式(5)的条件，参数值由大到小，分为第1档到第N档，分别记为F₁,F₂,…,F_N，每一档位的参数分别适应不同的工作范围，F₁的α、β、γ参数值最大，收敛速度最快，平滑精度最低，而F_N则相反，α、β、γ参数值最小，收敛速度最慢，平滑精度最高。

在开始阶段，先由F₁开始滤波，状态判决模块分析滤波结果及导航结算结果，判断接收机状态是否满足切换到其他档位进行滤波的条件，满足则控制切换α、β、γ参数值，转到滤波器的另一档位。

状态判决模块的判决依据主要有滤波结果的反馈和其他导航信息，其他导航信息包括位置解算结果、速度解算结果、卫星伪距、卫星信号载波相位等。

档位切换方法，通过估算滤波信息的方差，来判决当前状态适合使用哪个档位的滤波器，具体方法如下。引入统计量：

$L\left(k\right)=\underset{i=k-L+1}{\overset{k}{\mathrm{\&Sigma;}}}\frac{e^2\left(i\right)}{{\mathrm{\&sigma;}}_e^2}---\left(6\right)$

其中e(i)为i时刻的信息，即有：

$e\left(k\right)=y_k-{H\hat{X}}_{k+1/k}---\left(7\right)$

其方差：

${\mathrm{\&sigma;}}_e^2={\mathrm{\&sigma;}}_v^2/(1-\mathrm{\&alpha;})---\left(8\right)$

不同的滤波器档位，有着不同的信息方差，即第n档位的滤波器F_n的信息方差为：

${\mathrm{\&sigma;}}_n^2={\mathrm{\&sigma;}}_v^2/(1-{\mathrm{\&alpha;}}_n)---\left(9\right)$

其中α_n为F_n对应的α系数。

则如果当前处于F_m滤波器，判定接收机状态适合F_n的判决公式为：

${\mathrm{\&chi;}}_0^2\frac{{\mathrm{\&sigma;}}_{n+1}^2}{{\mathrm{\&sigma;}}_m^2}<L\left(k\right)<{\mathrm{\&chi;}}_0^2\frac{{\mathrm{\&sigma;}}_n^2}{{\mathrm{\&sigma;}}_m^2}---\left(10\right)$

其中是自由度为L的χ²变量，在显著水平为p₀时的临界值，即满足：

$P({\mathrm{\&chi;}}^2<{\mathrm{\&chi;}}_0^2)=1-p_0---\left(11\right)$

一般的取L＝30，p₀＝0.05，由χ²分布表查得

除此之外，还可引入其他导航信息来提高判决的准确性，如解算的位置信息、伪距信息、载波相位信息等。由于这些信息依赖于接收机的其他功能的准确性，故这里不做研究探讨。

实验测试发现，如果各档位滤波器的α、β、γ参数跳度过小会造成切换档位次数过多，收敛不够迅速，而参数跳度过大时，又会造成切换过程中滤波输出出现较大跳变，这是由于滤波器中含有保存的之前档位滤波器的预测状态。当提高滤波器档位时，其预测状态含有的噪声分量过大，不利于高档位滤波器的收敛；而当降低滤波器档位时，由于高档位的收敛相对低档位的滤波器慢，会造成预测值偏移真值过大，同样对滤波器的收敛造成影响。

故适当选取较大的参数跳度，每次仅提高或降低一个档位，并采取同时运算相邻档位滤波器的方式，事先对可能要切换的滤波器进行收敛计算，这样当满足切换条件切换时，滤波器能更快的收敛到稳定的工作状态，这样做可大为减少切换时带来的输出跳变。

有上述介绍可总结，使用多档α-β-γ滤波器的具体步骤包括：

步骤21：对卫星载波频率差值，通过第1档α-β-γ滤波器开始滤波；

步骤22：分析滤波结果及导航结算结果，判断接收机状态是否满足切换到其他档位进行滤波的条件，并将其反馈至档位切换模块；

步骤23：根据状态判决模块的反馈，在档位切换模块估算当前信息的方差，判决适合使用哪个档位的滤波器；

步骤24：将卫星载波频率差值，通过档位切换模块确定的合适档位滤波器进行滤波处理；

步骤25：循环步骤22至步骤24，直至收敛到稳定工作状态时，由输出切换模块输出滤波结果。

步骤3：将滤波后的各频率差值，组成求解接收机速度的方程组，利用最小二乘法求解出接收机速度。

接收机实现对信号稳定跟踪，本地载波频率可由载波NCO的频率控制字得到，将其减去设定的中频频率就可得到载波频率差值得观测量，数学模型可以写成如式(12)所示：

$f_d=\mathrm{DC}*(v_s-v_u)+c\stackrel{\&CenterDot;}{b}+n_d---\left(12\right)$

其中，是卫星和用户之间的径向方向余弦矢量；ν_s＝[ν_sx,ν_sy,ν_sz]^T是卫星的速度矢量；ν_u＝[ν_ux,ν_uy,ν_uz]^T是用户的速度矢量；b是钟漂；n_d是总的噪声项；[x_s,y_s,z_s]是卫星位置坐标；[x_u,y_u,z_u]是用户位置坐标。

假设接收机能够获取m颗卫星的多普勒观察量，根据式(12)可得如下方程组：

$\begin{array}{c}{f}_{d1}={\mathrm{DC}}_1*(v_{s1}-v_u)+c\stackrel{\&CenterDot;}{b}+n_{d1}\\ f_{d2}={\mathrm{DC}}_2*(v_{s2}-v_u)+c\stackrel{\&CenterDot;}{b}+n_{d2}\\ .\\ .\\ .\\ f_{\mathrm{dm}}={\mathrm{DC}}_m*(v_{\mathrm{sm}}-v_u)+c\stackrel{\&CenterDot;}{b}+n_{\mathrm{dm}}\end{array}---\left(13\right)$

式中，f_di是第m颗卫星的多普勒观测量，单位是m/s；DC_m是第m颗卫星的方向余弦矢量；v_si是第m颗卫星的速度，可由导航电文计算得到。

上组式中共有4个待解的未知量，分别为用户的速度矢量和钟漂为以下描述方便，统一用一个矢量表示。将上式中的卫星速度稍做整理，移到等号左边得到：

$\begin{array}{c}{f}_{d1}-{\mathrm{DC}}_1*v_{s1}={\mathrm{DC}}_1*v_u+c\stackrel{\&CenterDot;}{b}+n_{d1}\\ f_{d2}-{\mathrm{DC}}_2*v_{s2}={\mathrm{DC}}_2*v_u+c\stackrel{\&CenterDot;}{b}+n_{d2}\\ .\\ .\\ .\\ f_{\mathrm{dm}}-{\mathrm{DC}}_m*v_{\mathrm{sm}}={\mathrm{DC}}_m*v_u+c\stackrel{\&CenterDot;}{b}+n_{\mathrm{dm}}\end{array}---\left(14\right)$

如果定义：

$f_d=\left[\begin{array}{c}{f}_{d1}-{\mathrm{DC}}_1*v_{s1}\\ f_{d2}-{\mathrm{DC}}_2*v_{s2}\\ .\\ .\\ .\\ f_{\mathrm{dm}}-{\mathrm{DC}}_m*v_{\mathrm{sm}}\end{array}\right],$ $H=\left[\begin{array}{c}-\frac{x_{s1}-x_u}{r_1},-\frac{y_{s1}-y_u}{r_1},-\frac{z_{s1}-z_u}{r_1},1\\ -\frac{x_{s2}-x_u}{r_2},-\frac{y_{s2}-y_u}{r_2},-\frac{z_{s2}-z_u}{r_2},1\\ .\\ .\\ .\\ -\frac{x_{\mathrm{sm}}-x_u}{r_m},-\frac{y_{\mathrm{sm}}-y_u}{r_m},-\frac{z_{\mathrm{sm}}-z_u}{r_m},1\end{array}\right],$ $n_d=\left[\begin{array}{c}{n}_{d1}\\ n_{d2}\\ .\\ .\\ .\\ n_{\mathrm{dm}}\end{array}\right]$

则式可以写成矩阵的形式，即：

f_d＝HX_u+n_d    (15)可以直接用最小二乘方法求解：

${\hat{X}}_u={\left(H^{T}H\right)}^{-1}{H}^T{f}_d---\left(16\right)$

在接收机内部利用LS方法求解速度，一般的顺序先求出位置，卫星运动方向余弦矢量形成H矩阵，即可根据上式得到用户速度。

步骤4：发明设计一种针对速度求解的RAIM算法，排除速度解算方程组中的错误方程，使解算速度不出现错误。

(1)速度RAIM算法概念

RAIM技术是设在用户接收机中的一种算法，它利用导航卫星的冗余信息，对多个导航解进行一致性检验，从而达到完好性监测的目的。

RAIM的理论基础是相差的探测和分离理论，它需要解决两个问题：卫星是否存在故障和故障存在于哪颗卫星。当可视卫星数大于4颗时，可以进行故障监测，当可视卫星数大于5颗时，可以进行故障识别。目前RAIM算法主要用于位置解算中，当用于求解位置的伪距存在较大误差时，通过检测伪距位置方程可以排除故障卫星，保证定位结果正确。

由上式可以看出，卫星多普勒频率值的误差对测速精度有较大的影响。而卫星载波多普勒频率值是从载波跟踪环读取，环路的失锁或则抖动将得到错误的多普勒频率，从而导致速度解算错误。与位置解算类似，速度解算也采用最小二乘法解算，同样可以构造一种RAIM算法对速度解算方程进行检测排除错误方程。下面本文对RAIM技术用于卫星多普勒频率的故障检测和识别的进行分析。

(2)基于残差平方和的RAIM故障检测：

${\hat{f}}_d=H*{\hat{X}}_u---\left(17\right)$

用ω∈R来表示残差向量，则有：

$\begin{array}{c}\mathrm{\&omega;}=f_d-{\hat{f}}_d=f-H*{\hat{X}}_u\\ =[I-H{\left(H^{T}H\right)}^{-1}{H}^T]*f_d\\ =[I-H{\left(H^{T}H\right)}^{-1}{H}^T]*f_d\\ =[I-H{\left(H^{T}H\right)}^{-1}{H}^T]*[{\mathrm{HX}}_u+n_d]\\ =[H-H{\left(H^{T}H\right)}^{-1}{H}^{T}H]*X_u+[I-H{\left(H^{T}H\right)}^{-1}{H}^T]*n_d\\ =[H-H]*X_u+[I-H{\left(H^{T}H\right)}^{-1}{H}^T]*n_d\\ =[I-H{\left(H^{T}H\right)}^{-1}{H}^T]*n_d\\ =S*n_d\end{array}---\left(18\right)$

由上式知，向量ω中包含了卫星多普勒频率误差信息，可知作为判断有无故障卫星的依据。定义它表示各个卫星伪距残差的平方和，它是一个自由度为(n－4)的χ²分布，在实际应用中，通常使用作为统计检测量的，在每次定位时将它与门限值T_D比较，若T_X≥T_D，有故障星；若T_X<T_D，无故障星。

当统计检测量超过门限时，就认为存在故障星。此时，可以再根据最大似然估计方法，找到出现故障可能性最大的那颗卫星。因为ω＝S*n_d，如果某一颗卫星的值达到最大，则认为该星出现故障的可能性最大。

(3)速度RAIM检测门限的确定

根据统计理论，假设卫星测量值服从零均值的高斯白噪声分布，则在卫星无故障时，若给定虚警概率值P_FA，就可以将检测门限T_D计算出来。

$T_D={\mathrm{\&sigma;}}_n^2*Q^{-1}\left(P_{\mathrm{FA}}\right|n-4)---\left(19\right)$

式中，表示测量的等效测距误差；Q(χ²|r)＝1-P(χ²|r)；这里P(χ²|r)为自由度是r的χ²分布，它的概率密度函数为：

$f\left(x\right)=\left\{\begin{array}{cc}\frac{1}{2^{\frac{r}{2}}\mathrm{\&Gamma;}\left(\frac{r}{2}\right)}{x}^{\frac{r}{2}-1}{e}^{-\frac{r}{2}}& x>0\\ 0& x<0\end{array}\right.---\left(20\right)$

由式(10)可知，门限值T_D与卫星的伪距测量误差、卫星数和允许的虚警概率有关。现以观测6颗卫星为例来说明如何计算相应的门限。当n＝6时,系统自由度为r＝n-4＝2,此时相应的分布密度函数为：

$f\left(x\right)=\frac{1}{2}{e}^{-x/2},x\&GreaterEqual;0---\left(21\right)$

假设系统等效测距误差为σ_n＝20m，虚警概率取PFA＝1/15000，则有：

${\&Integral;}_0^{\mathrm{\&alpha;}}f\left(x\right)\mathrm{dx}=1-\frac{1}{15000}---\left(22\right)$

于是，可以计算出积分上限和检测门限为：

α＝2ln15000＝19.2316

$T_D=\sqrt{\mathrm{\&alpha;}*{\mathrm{\&alpha;}}_n^2/n-4}=\sqrt{(19.2316\&times;{20}^2)/2}=62.02m---\left(23\right)$

同理，不同卫星数的检测门限也可按此方法计算。

因此，针对速度求解的RAIM算法具体步骤包括：

步骤41：先用位置解算使用的那几颗卫星的多普勒观测量做一次速度最小二乘解算，然后开始执行速度RAIM算法；

步骤42：如果解算的卫星数小于或等于4颗，则速度RAIM算法无法检测故障星，直接输出计算的速度结果；

步骤43：若卫星数超过4颗，则根据速度最小二乘法的计算结果，算出多普勒残余误差向量ε及速度RAIM算法统计检测量T_X；

步骤44：将T_X与对应的检测门限T_D比较，若T_X≤T_D，没有故障星，RAIM算法完成，并输出定位测速结果；

步骤45：若T_X>T_D，则判断卫星数是否超过5颗。若星数未超过5颗，无法排除故障星，RAIM算法结束，同时输出测速结果不可用的标志，提醒使用者；

步骤46：若星数超过5颗，则可以找到发生故障的卫星，多普勒残余最大的方程判断为错误方程；

步骤47：将错误方程剔除后，用剩余的卫星做速度最小二乘解算，然后重复前面的步骤，直至速度RAIM算法统计检测量T_X小于检测门限T_D，或方程组的个数小于5，输出解算结果。

基于上述算法，对GPS接收机进行12小时连续测量，记录测试结果进行分析如图2和3所示。其中，图2为没有速度RAIM算法的结果，在速度解算时只执行前3个步骤。测试过程如下，接收机接收实际信号，连续工作记录12小时，接收机处于静止状态，因此记录的速度结算结果为接收机的测速误差。图3为加上RAIM算法的结果，即，在速度解算时完整的执行4个步骤。在静态下，GPS接收机的实际速度为零，因此图2和图3测试的结果均为接收机的测速误差。从图2可以看出，在不用速度RAIM算法的情况下，测速误差最大可达1.6m/s，多处误差大于0.2m/s。从图3可以看出，加上速度RAIM算法的情况下，测速误差最大为0.2m/s，通常小于0.2m/s。因此，从图2和图3的对比中可以看出速度RAIM算法通过排除错误卫星参与计算提高了测速的精度，使测速误差小于0.2m/s。

以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。

Claims (4)

1.一种高动态卫星导航接收机精确测速方法，其中，卫星导航接收机包括射频部分和基带部分：射频部分完成导航信号的放大、滤波及变频功能；基带部分又可分为基带信号处理部分和基带信息处理部分，信号处理部分包括对各个导航卫星信号的捕获、跟踪以及解调功能，信息处理部分根据信号处理电路得到的信息完成位置、速度与时间的解算，其中基带部分包括一个嵌入式CPU来运行信息处理程序，

其特征是所述方法包括下列步骤：

步骤1：通过直接读取载波跟踪环频率误差量的方法提取卫星载波与本地载波的频率差值；

步骤2：设计一种多档α-β-γ滤波器对提取的卫星载波频率差值进行滤波处理；

步骤3：将滤波后的各频率差值，组成求解接收机速度的方程组，利用最小二乘法求解出接收机速度；

步骤4：设计一种RAIM算法，排除速度解算方程组中的错误方程，确保解算速度不出现错误。

2.根据权利要求1所述的一种高动态卫星导航接收机精确测速方法，其特征是所述多档α-β-γ滤波器由α-β-γ滤波器模块、状态判决模块、档位切换模块、输出切换模块组成；

所述α-β-γ滤波器模块输入端口与所述档位切换模块输出端口相连，所述α-β-γ滤波器模块输出端口与所述输出切换模块输入端口相连，用于对提取的卫星载波频率差值进行滤波处理；

所述状态判决模块输入端口与所述输出切换模块相连，用于反馈滤波结果；

所述状态判决模块输出端口分别与所述档位切换模块、输出切换模块相连，用于判断接收机状态是否满足切换到其他档位进行滤波的条件。

3.根据权利要求1所述的一种高动态卫星导航接收机精确测速方法，其特征是所述使用多档α-β-γ滤波器的具体步骤：

步骤21：对卫星载波频率差值，通过第1档α-β-γ滤波器开始滤波；

步骤22：分析滤波结果及导航结算结果，判断接收机状态是否满足切换到其他档位进行滤波的条件，并将其反馈至档位切换模块；

步骤23：根据状态判决模块的反馈，在档位切换模块估算当前信息的方差，判决适合使用哪个档位的滤波器；

步骤24：将卫星载波频率差值，通过档位切换模块确定的合适档位滤波器进行滤波处理；

步骤25：循环步骤22至步骤24，直至收敛到稳定工作状态时，由输出切换模块输出滤波结果。

4.根据权利要求1所述的一种高动态卫星导航接收机精确测速方法，其特征是所述速度求解的RAIM算法的具体步骤：

步骤41：先用位置解算使用的那几颗卫星的多普勒观测量做一次速度最小二乘解算，然后开始执行速度RAIM算法；

步骤42：如果解算的卫星数小于或等于4颗，则速度RAIM算法无法检测故障星，直接输出计算的速度结果；

步骤43：若卫星数超过4颗，则根据速度最小二乘法的计算结果，算出多普勒残余误差向量ε及速度RAIM算法统计检测量T_X；

步骤44：将T_X与对应的检测门限T_D比较，若T_X≤T_D，没有故障星，RAIM算法完成，并输出定位测速结果；

步骤45：若T_X>T_D，则判断卫星数是否超过5颗。若星数未超过5颗，无法排除故障星，RAIM算法结束，同时输出测速结果不可用的标志，提醒使用者；

步骤46：若星数超过5颗，则可以找到发生故障的卫星，多普勒残余最大的方程判断为错误方程；

步骤47：将错误方程剔除后，用剩余的卫星做速度最小二乘解算，然后重复前面的步骤，直至速度RAIM算法统计检测量T_X小于检测门限T_D，或方程组的个数小于5，输出解算结果。