CN109541659A - 一种基于北斗的地基增强***载波相位平滑伪距方法 - Google Patents

Abstract

本发明公开了一种基于北斗的地基增强***载波相位平滑伪距方法，为了探测载波相位是否存在周跳，根据伪距与载波相位观测值计算载波差量，并利用Gaussian膨胀法确定载波差量周跳判决门限值；为了减小异常电离层干扰对载波相位平滑伪距精度影响，通过计算伪码‑载波偏离度以自适应调节Hatch滤波器平滑窗口时间长度，当电离层出现异常时减小载波相位在平滑滤波器中的权重。本发明能够有效提高伪距平滑精度，同时能够提高机载端位置计算精度。

Description

一种基于北斗的地基增强***载波相位平滑伪距方法

技术领域

本发明属于星基导航增强***领域，具体涉及一种基于北斗的地基增强***载波相位平滑伪距方法。

背景技术

飞机的精密进近与着陆是飞行过程中最为关键的阶段，其对导航***精度与完好性的要求较高，将北斗作为一种主要的导航设备应用于飞机的进近与着陆是未来民用航空导航的发展趋势。民用航空导航使用的北斗标准定位服务，难以满足其在定位精度、连续性、完好性和可用性等方面的需求，因此出现了基于北斗的各种增强***，其中地基增强***(Ground Based Augmentation System,GBAS)通过差分定位提高了卫星导航精度，并增加了一系列完好性监测算法，能够为机场覆盖终端区空域范围内的配置相应机载设备的飞机提供I类(CAT-I)甚至更高等级精密进近、着陆引导服务.美国、欧洲、日本等先后开展了基于全球定位***(global positioning system，GPS)的地基增强***研究工作，其中美国Honeywell公司研制的满足CAT-I类精密进近GBAS于2009年获得首次***设计认证，澳大利亚悉尼机场于2009年安装了Honeywell研制的SLS-4000型GBAS。我国的北斗卫星导航***(Beidou Navigation System,BDS)截止2012年已经具备向亚太地区提供无源定位、导航与授时服务，2014年我国正式启动基于BDS的地基增强***工程建设，目前已经投入使用，在全国建设了150个框架基站，能够提供米级精度的定位服务。目前国内科研院所正在积极开展基于北斗的差分定位与完好性监测等相关研究工作。

地基增强***基准站接收机能够实时接收卫星信号并根据卫星的轨道参数计算得到卫星的空间位置坐标，利用卫星的空间位置坐标以及基准站的精确位置信息可以得到基准站与卫星的几何距离值，与卫星播发的伪距信息进一步可以得到伪距差分修正值。由于一个C/A码的波长约300m，其测量噪声误差一般为波长的1％，即为3m，同时多路径误差会影响伪距测量精度，而载波相位测量值虽然含有整周模糊度，但是其非常平滑，其精度高达毫米量级。载波相位平滑伪距将伪距无整周模糊度与载波相位观测量高精度的特点相结合，能够有效提高伪距测量精度。目前用于伪距平滑的方法主要有Hatch滤波器和Kalman滤波器，这两种方法都能实现伪距的平滑，但是对平滑过程中电离层异常以及相位周跳的研究较少。

发明内容

发明目的：针对现有技术的不足，本发明提出一种基于北斗的地基增强***载波相位平滑伪距方法，利用伪距与载波相位观测值计算载波差量，从而判断载波相位是否存在周跳，进一步当载波相位无周跳发生时，根据伪距与载波相位观测值计算伪码-载波偏离度，以此判断电离层是否存在异常，如果电离层发生异常，则自适应调节Hatch滤波器平滑窗口时间大小，克服异常电离层干扰对伪距平滑精度影响。

技术方案：本发明所述的一种基于北斗的地基增强***载波相位平滑伪距方法，包括以下步骤：

步骤1，根据北斗卫星导航***(Beidou Navigation System,BDS)接收机伪距与载波相位观测值计算载波差量，并与计算得到的载波差量门限值进行比较，如果载波差量大于载波差量门限值，说明载波相位发生周跳，则初始化滤波器，否则根据伪距与载波相位观测值计算伪码-载波偏离度值，并与设置的伪码-载波偏离度门限值进行比较；

步骤2，如果伪码-载波偏离度大于门限值，则说明电离层发生异常，此时平滑窗口时间为T1，如果伪码-载波偏离度不大于门限值，说明电离层较为平稳，此时设置平滑窗口时间为T2，利用Hatch滤波器对伪距进行平滑；

步骤1中，根据如下公式计算载波差量temp：

temp＝|λ(φ_k-φ_k-1)-(ρ_k-ρ_k-1)|，

式中，λ为载波波长，φ_k和φ_k-1分别为历元k时刻的载波相位观测值和k-1时刻的载波相位观测值，ρ_k和ρ_k-1分别为k时刻的伪距观测值和k-1时刻的伪距观测值。

步骤1中，根据如下方法计算得到载波差量门限值：

步骤1-1，计算得到每颗卫星每个历元时刻的载波差量；

步骤1-2，按卫星仰角划分载波差量区间，以10°为一个区间进行划分，并分别计算每个区间内载波差量的均值μ_elev与标准差σ_elev；

步骤1-3，利用5阶多项式拟合其它角度均值与标准差，其表达式为：

y＝A₅x⁵+A₄x⁴+A₃x³+A₂x²+A₁x+A₀，

式中，A_i为需要拟合的第i个参数，i取值为0～5，x为输入卫星仰角，其取值分别为5°、15°、25°、35°、45°、55°、65°、75°、85°，y为需要拟合的载波差量均值或者标准差；

步骤1-4，对步骤1-1得到的载波差量进行归一化；

步骤1-5，统计归一化后的载波差量，绘制归一化后载波差量分布直方图，并根据每个间隔点内的载波差量数据个数与载波差量总个数计算其概率密度值，即利用间隔内载波差量数据个数除以载波差量总个数，由此得到其概率密度值；

步骤1-6，计算膨胀系数f；

步骤1-7，根据得到的膨胀系数f以及均值μ_elev与标准差σ_elev获得该组样本数据的载波差量门限值Threshold_elev：

Threshold_elev＝μ_elev±6fσ_elev。

步骤1-4中，采用如下公式对样本参数进行归一化：

式中，x表示载波差量采样数据，表示归一化后的载波差量值；

步骤1-6包括：根据步骤1-5计算得到的概率密度分布情况绘制样本数据的概率密度分布图，并在其中绘制标准正态分布概率密度曲线，然后对标准正态分布的概率密度曲线进行膨胀，直至膨胀所得到的曲线完全包络样本数据的概率密度曲线两侧的样本数据，从而得到膨胀系数f。

步骤1中，根据如下公式计算伪码-载波偏离度：

ρ_k＝ρ_IF,k+I_k+ε_ρ,k，

φ_k＝ρ_IF,k-I_k+λN+ε_φ,k，

ρ_IF,k＝r_k+c[δt_u,k-δt_s,k]+T_k，

式中，ρ_IF,k为k时刻无电离层影响的伪距，r_k为k时刻星站真实的几何距离，c为光速，δt_u,k和δt_s,k分别为k时刻接收机钟差和卫星钟差，I_k和T_k分别为k时刻电离层延迟和对流层延迟，λ为波长，N为整周模糊度，ε_ρ,k和ε_φ,k分别为k时刻接收机端伪距和相位观测误差，包括多路径误差；

定义z_k表示k时刻的伪距观测值ρ_k与载波相位观测值φ_k之差：

z_k＝ρ_k-φ_k

＝2I_k+ε_ρ,k-ε_φ,k-λ·N

设定观测时段内，接收机一直锁定，无载波失锁和失周，整周模糊度N值保持不变，同时令δε_k＝ε_ρ,k-ε_φ,k，则有：

式中，dz_k表示待估计伪码-载波偏离度，δε_k表示伪距噪声与载波相位噪声差分量，表示估计得到的电离层延时速率，其大小与卫星仰角、空间电子浓度值等因素相关，T_s＝1s，表示相邻两个历元间隔时间；

为了抑制伪距噪声、多路径测量误差等对伪码-载波偏离度估计精度影响，采用二阶线性时不变低通滤波器滤除dz_k中高频噪声，计算公式如下：

式中，d_1,k与d_1,k-1为过程中间变量，Dvgc_k与Dvgc_k-1表示k时刻经过低通滤波后的伪码-载波偏离度与k-1时刻经过低通滤波后的伪码-载波偏离度；τ₁＝τ₂＝30s，表示单位观测时刻的长度。

步骤1中，根据大量实验，设置伪码-载波偏离度门限值为0.02m/s。

步骤2中，根据大量实验，设置T1为20s，T2为100s。

有益效果：

1、本发明提出的载波相位周跳判决门限值利用Gaussian膨胀法对样本数据进行包络，减小完好性虚警率；

2、本发明提出了利用伪码-载波偏离度监视电离层异常，从而实现了Hatch滤波器平滑窗口时间长度的自适应调节；

3、本发明为了减小异常电离层干扰对伪距平滑精度影响，提出了一种平滑窗口时间长度可变的设计思路，通过减小时间长度值以达到减小载波相位在滤波器中的权重。

本发明提出了一种基于伪码-载波偏离度的自适应Hatch滤波器，该方法根据伪码-载波偏离度值判断电离层是否存在异常，如果电离层发生异常，则通过减小Hatch滤波器平滑窗口时间以达到减小载波相位在滤波器中的权重，从而能够克服异常电离层干扰对伪距平滑精度影响，同时本发明对载波相位周跳探测及其门限值确定方法进行了研究，根据伪距与载波相位观测值计算载波差量并利用Gaussian膨胀法确定得到载波差量门限值。

附图说明

下面结合附图和具体实施方式对本发明做更进一步的具体说明，本发明的上述或其他方面的优点将会变得更加清楚。

图1是载波差量门限值生成流程。

图2是载波差量随卫星仰角变化曲线。

图3是载波差量标准差随卫星仰角变化曲线。

图4是样本数据经过归一化后随卫星仰角分布曲线。

图5是归一化后样本数据分布直方图。

图6是载波差量分布概率密度曲线。

图7是载波差量及其阈值曲线。

图8是自适应平滑窗口计算流程。

图9是载波相位平滑伪距算法流程。

图10a为4号卫星实验结果。

图10b为7号卫星实验结果。

具体实施方式

下面结合附图及实施例对本发明做进一步说明。

如图1所示，本发明提供了一种基于北斗的地基增强***载波相位平滑伪距方法，图1为载波差量门限值生成流程，包括以下步骤：

步骤1，接收机采集卫星信号，由此计算得到每颗卫星的载波差量，图2为采集一段时间计算得到的载波差量随卫星仰角变化曲线，图中黑色点线表示卫星载波差量，黑色虚线表示利用Gaussian膨胀法计算得到的门限值；

步骤2，按卫星仰角划分参数区间，以10°为一个区间划分采样数据，并分别计算不同区间内样本参数的均值μ_elev与标准差σ_elev，具体区间划分为(0°,10°]∪(10°,20°]∪(20°,30°]…(80°,90°]；

步骤3，利用5阶多项式拟合其它角度均值与标准差，其表达式为：

y＝A₅x⁵+A₄x⁴+A₃x³+A₂x²+A₁x+A₀

式中，A_i(i＝0,1,2,3,4,5)为需要拟合的参数，x为输入卫星仰角，其取值分别为5°、15°、25°、35°、45°、55°、65°、75°、85°，y为需要拟合的载波差量均值或者标准差。

图3为载波差量标准差随卫星仰角变化曲线，图中黑色圈为每个区间的标准差，黑色实线为经过高阶拟合后得到的载波差量标准差随仰角变化曲线，黑色虚线为经过Gaussian膨胀后得到的载波差量标准差曲线，膨胀系数f＝1.4256.

步骤4，对样本参数进行归一化，具体计算过程为将采样数据减去步骤3中求得的均值μ_elev，然后除以相应角度下的标准差σ_elev

式中，x表示载波差量采样数据，μ_elev和σ_elev由步骤3拟合得到，表示归一化后的载波差量值。

图4为样本数据经过归一化后随卫星仰角分布曲线，由图可以看出，经过归一化后的载波差量幅值变化平稳，其均值μ＝-0.028，σ＝1.056，这也进一步说明步骤3中高阶多项式拟合结果的正确性。

步骤5，归一化参数分布统计，绘制归一化数据的分布直方图，并根据每个间隔点内的样本数据与总样本数量计算其概率密度值，图5为归一化后样本数据分布直方图；

步骤6，计算膨胀因子，根据步骤5计算得到的概率密度分布情况绘制样本数据的概率密度分布图，并在其中绘制标准正态分布概率密度曲线，然后对标准正态分布的概率密度曲线进行膨胀，直至膨胀所得到的曲线完全包络样本数据的概率密度曲线两侧的样本数据，从而得到膨胀系数为f，如图6所示，图中黑色点线为样本实际概率分布曲线，黑色实线为服从N(0,1)的高斯分布概率密度曲线，黑色虚线为经过高斯膨胀后服从N(0,1.4256)的高斯分布概率密度曲线，从而得到膨胀系数为f＝1.4256；

步骤7，步骤7生成示警阈值，根据得到的膨胀系数f以及步骤3中多项式插值计算得到的均值μ_elev与标准差σ_elev获得该组样本数据的门限值Threshold_elev

Threshold_elev＝μ_elev±6fσ_elev

图7为GBAS工作一段时间14号卫星的载波差量随观测时间变化曲线，其示门限值为利用本文所提出的方法计算得到。

如图8所示，本发明提供了一种基于北斗的地基增强***载波相位平滑伪距方法，图8为自适应平滑窗口计算流程，包括以下步骤：

(1)根据BDS接收机输出伪距与载波相位计算伪码-载波偏离度，监测电离层是否存在异常，如果电离层发生异常，设置平滑窗口时间为20s，否则平滑窗口时间设置为100s；

(2)利用Hatch滤波器对伪距进行平滑，得到较为光滑的伪距值；

(3)对一个历元内所有可见卫星进行完好性监测，形成可用卫星集合并计算得到伪距差分校正值。

如图9所示，本发明提供了一种基于北斗的地基增强***载波相位平滑伪距方法，图9为载波相位平滑伪距算法流程，包括以下步骤：

(1)根据BDS接收机伪距与载波相位观测值计算载波差量，并与设置的门限值进行比较，如果载波相位发生周跳，则初始化滤波器；

(2)计算伪码-载波偏离度值，并与其设置的示警阈值进行比较；

(3)如果电离层发生异常则平滑窗口时间选择为20s，否则设置平滑窗口时间为100s，利用Hatch滤波器对伪距进行平滑；

(4)对一个观测历元内所有卫星的完好性进行监测，形成可用卫星集合及其伪距差分校正值。

实施例

下面以中国电科某研究所建立的试验型GBAS基准站为实验对象，验证本发明提出的差分定位方法的准确性。该试验型GBAS地面***包括4台基准接收机、1台GNSS信号源、1台数据处理***、VDB(Very High Frequency,VHF Data Broadcasting)电台以及天线，其中，GNSS信号源能够模拟生成GPS/BDS信号，并能够模拟飞机进近航迹；架设在楼顶天线可用来接收GPS/BDS信号；VDB电台用来广播与接收地面站播发的差分数据与完好性参数，同时，用户***包括1台机载接收机及相应的天线，实验步骤如下：

步骤1，将4台接收机设置为基准接收机，另一台接收机设置为机载端移动接收机，天线放置于楼顶位置处；

步骤2，分别采集地面站与机载端接收机数据，包括：星历数据、观测数据；

步骤3，利用地面站软件对4台基准接收机数据进行处理，保存可观测卫星的伪距差分修正值，同时进行机载端定位解算，读取伪距差分修正值进行误差校正，输出最终的位置信息。

采集了多组实验数据，选择其中一组数据进行了详细分析，图10a和图10b为地面站伪距差分值对比曲线，其中黑色虚线为无载波平滑时的伪距差分值，黑色实线为有载波平滑时的伪距差分值，图10a为4号卫星实验结果，图10b为7号卫星实验结果。

本发明提供了一种基于北斗的地基增强***载波相位平滑伪距方法，具体实现该技术方案的方法和途径很多，以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。本实施例中未明确的各组成部分均可用现有技术加以实现。

Claims (8)

1.一种基于北斗的地基增强***载波相位平滑伪距方法，其特征在于，包括以下步骤：

步骤1，根据北斗卫星导航***接收机伪距与载波相位观测值计算载波差量，并与计算得到的载波差量门限值进行比较，如果载波差量大于载波差量门限值，说明载波相位发生周跳，则初始化滤波器，否则根据伪距与载波相位观测值计算伪码-载波偏离度值，并与设置的伪码-载波偏离度门限值进行比较；

步骤2，如果伪码-载波偏离度大于门限值，则说明电离层发生异常，此时平滑窗口时间为T1，如果伪码-载波偏离度不大于门限值，说明电离层较为平稳，此时设置平滑窗口时间为T2，利用Hatch滤波器对伪距进行平滑。

2.根据权利要求1所述的方法，其特征在于，步骤1中，根据如下公式计算载波差量temp：

temp＝|λ(φ_k-φ_k-1)-(ρ_k-ρ_k-1)|，

式中，λ为载波波长，φ_k和φ_k-1分别为历元k时刻的载波相位观测值和k-1时刻的载波相位观测值，ρ_k和ρ_k-1分别为k时刻的伪距观测值和k-1时刻的伪距观测值。

3.根据权利要求2所述的方法，其特征在于，步骤1中，根据如下方法计算得到载波差量门限值：

步骤1-1，计算得到每颗卫星每个历元时刻的载波差量；

步骤1-2，按卫星仰角划分载波差量区间，以10°为一个区间进行划分，并分别计算每个区间内载波差量的均值μ_elev与标准差σ_elev；

步骤1-3，利用5阶多项式拟合其它角度均值与标准差，其表达式为：

y＝A₅x⁵+A₄x⁴+A₃x³+A₂x²+A₁x+A₀，

式中，A_i为需要拟合的第i个参数，i取值为0～5，x为输入卫星仰角，其取值分别为5°、15°、25°、35°、45°、55°、65°、75°、85°，y为需要拟合的载波差量均值或者标准差；

步骤1-4，对步骤1-1得到的载波差量进行归一化；

步骤1-5，统计归一化后的载波差量，绘制归一化后载波差量分布直方图，并根据每个间隔点内的载波差量数据个数与载波差量总个数计算其概率密度值，即利用间隔内载波差量数据个数除以载波差量总个数，由此得到其概率密度值；

步骤1-6，计算膨胀系数f；

步骤1-7，根据得到的膨胀系数f以及均值μ_elev与标准差σ_elev获得该组样本数据的载波差量门限值Threshold_elev：

Threshold_elev＝μ_elev±6fσ_elev。

4.根据权利要求3所述的方法，其特征在于，步骤1-4中，采用如下公式对样本参数进行归一化：

式中，x表示载波差量采样数据，表示归一化后的载波差量值。

5.根据权利要求4所述的方法，其特征在于，步骤1-6包括：根据步骤1-5计算得到的概率密度分布情况绘制样本数据的概率密度分布图，并在其中绘制标准正态分布概率密度曲线，然后对标准正态分布的概率密度曲线进行膨胀，直至膨胀所得到的曲线完全包络样本数据的概率密度曲线两侧的样本数据，从而得到膨胀系数f。

6.根据权利要求5所述的方法，其特征在于，步骤1中，根据如下公式计算伪码-载波偏离度：

ρ_k＝ρ_IF,k+I_k+ε_ρ,k，

φ_k＝ρ_IF,k-I_k+λN+ε_φ,k，

ρ_IF,k＝r_k+c[δt_u,k-δt_s,k]+T_k，

式中，ρ_IF,k为k时刻无电离层影响的伪距，r_k为k时刻星站真实的几何距离，c为光速，δt_u,k和δt_s,k分别为k时刻接收机钟差和卫星钟差，I_k和T_k分别为k时刻电离层延迟和对流层延迟，λ为波长，N为整周模糊度，ε_ρ,k和ε_φ,k分别为k时刻接收机端伪距和相位观测误差；

定义z_k表示k时刻的伪距观测值ρ_k与载波相位观测值φ_k之差：

z_k＝ρ_k-φ_k

＝2I_k+ε_ρ,k-ε_φ,k-λ·N

设定观测时段内，接收机一直锁定，无载波失锁和失周，整周模糊度N值保持不变，同时令δε_k＝ε_ρ,k-ε_φ,k，则有：

式中，dz_k表示待估计伪码-载波偏离度，δε_k表示伪距噪声与载波相位噪声差分量，表示估计得到的电离层延时速率，T_s＝1s，表示相邻两个历元间隔时间；

采用二阶线性时不变低通滤波器滤除dz_k中高频噪声，计算公式如下：

式中，d_1,k与d_1,k-1为过程中间变量，Dvgc_k与Dvgc_k-1表示k时刻经过低通滤波后的伪码-载波偏离度与k-1时刻经过低通滤波后的伪码-载波偏离度；τ₁＝τ₂＝30s，表示单位观测时刻的长度。

7.根据权利要求6所述的方法，其特征在于，步骤1中，设置伪码-载波偏离度门限值为0.02m/s。

8.根据权利要求7所述的方法，其特征在于，步骤2中，设置T1为20s，T2为100s。