BDS广播星历A1型故障完好性监测方法

2021-03-01刘媛媛

导航定位学报 2021年1期

刘媛媛，李亮，李慧，程春

（哈尔滨工程大学智能科学与工程学院，哈尔滨 150001）

0 引言

北斗卫星导航系统（BeiDou navigation satellite system, BDS）是我国自行研制的全球卫星导航系统（global navigation satellite system, GNSS），可全天候、全天时为各类用户提供高精度、高可靠的定位、导航及授时服务，具备支持民航的导航能力[1]。将卫星导航应用于民航系统，可有效地提高飞行安全和效率。服务于民航的导航系统，其性能需求主要体现在精度、完好性、连续性和可用性4个方面[2]。飞机精密进近的完好性和连续性，是航空用户最关注的性能需求，其中卫星广播星历故障，是影响飞机精密进近的主要风险源之一。为保障与生命安全相关用户的定位完好性，对广播星历故障的监测是必不可少的。

卫星广播星历故障可根据是否与卫星机动有关而分为两大类型[3]：

1）A型故障。发生卫星机动后，广播星历参数错误；包括地面站知道卫星机动发生的A1型故障和地面站不知道卫星机动发生的A2型故障。

2）B型故障。没有发生卫星机动，但广播星历参数错误。

广播星历B型故障可根据预测星历减今日星历监测的方法进行监测。本文主要研究利用广播星历监测A型故障的方法，由于A2型故障发生的概率低于系统要求的误警率，所以文中仅对广播星历中A1型故障监测方法展开研究。

对广播星历A1型故障的监测，目前应用最多的方法是接收机自主完好性监测（receiver autonomous integrity monitoring, RAIM）算法[4]。大多数RAIM算法采用伪距观测量进行完好性监测，由于其观测噪声较大、定位精度不高，只能满足航路飞行和非精密进近应用，不能满足精密进近服务的要求[5]。基于载波相位观测量的实时导航技术，能够满足飞机精密进近的精度要求[6]，在短基线条件下，观测量的双差能够消除大部分测量误差。虽然基于载波相位观测量完好性监测的精度较高，但存在整周模糊度可靠解算等问题，增加了解算难度[7]。部分学者也采用伪距和载波相位观测量联合的方式来实现完好性监测[8]，但仍存在大量的模糊度解算等问题，增加了系统算法的复杂度。

考虑到伪距和载波相位观测量在进行广播星历A1型故障监测时存在的问题，本文提出顾及参数灵敏度的广播星历A1型故障完好性监测方法。本文以北斗卫星导航（区域）系统即北斗二号（BeiDou navigation satellite （regional） system,BDS-2）、第三代北斗卫星导航系统即北斗三号（BeiDou navigation satellite system with global coverage,BDS-3）的C01～C30卫星为研究对象，采用伪距和载波相位观测量联合的方式测量轨道误差，并与卫星位置误差对广播星历参数变化的灵

敏度相结合，将关于卫星位置误差的观测方程转化为关于广播星历参数偏差的观测方程，在提高观测精度的同时，降低载波相位观测中模糊度解算的复杂度。本文首先分析卫星位置误差对广播星历参数变化的灵敏度，获取每个广播星历参数的不确定性对卫星位置精度的影响大小及位置误差函数。其次，采用伪距差分和双差载波相位测量轨道误差，与灵敏度分析中的位置误差函数结合构建观测方程。最后，根据系统的完好性风险要求计算最小可检测误差（minimum detectable error,MDE），并与满足地基增强系统（ground-based augmentation systems, GBAS）完好性风险要求的MDE检测门限比较，实现对BDS广播星历A1型故障的完好性监测。

1 广播星历参数灵敏度分析

通过卫星位置误差对广播星历参数变化的灵敏度进行分析，可获得每个广播星历参数的不确定性对卫星位置精度的影响大小以及关于广播星历参数偏差的位置误差函数，将其和测量轨道误差的观测方程结合，可转化为关于广播星历参数偏差的观测方程，以降低模糊度解算的复杂度。本文采用2019年全年BDS-2、BDS-3 C01～C30卫星的广播星历进行灵敏度分析。

广播星历正常的情况下，卫星位置误差对广播星历参数变化的灵敏度是可确定的。广播星历参数描述了在一定拟合间隔下的卫星轨道，它包括15个广播星历参数和1个星历参考时间toe。为书写方便，本文用矩阵p表示15个广播星历参数，即

式中：M0为平近点角；Δn为平均角速度校正值；e为偏心率；a为长半径；0Ω为升交点赤经；0i为轨道倾角；ω为近地角距；Ω˙为升交点赤经变化率；i˙为轨道倾角变化率；Cuc为升交点角距余弦调和校正振幅；Cus为升交点角距正弦调和校正振幅；为轨道半径余弦调和校正振幅；Crs为轨道半径正弦调和校正振幅；Cic为轨道倾角余弦调和校正振幅；Cis为轨道倾角正弦调和校正振幅。

根据15个广播星历参数和星历参考时间toe，可以在任意时刻t，使用BDS接口控制文件中定义的卫星位置非线性函数，计算卫星位置，计算公式为

式中：x、y、z分别为3维卫星位置的计算函数。卫星位置误差对广播星历参数变化的灵敏度可以通过卫星位置函数对每个广播星历参数的微分体现，具体公式为

式中：等式右边的第1个矩阵为3×15的灵敏度矩阵，该矩阵中的每个项是15个卫星广播星历参数的函数。可将式（3）的位置误差函数简写为

将式（4）与广播星历参数偏差的标准差σ( δpi)结合，可以得到每个广播星历参数变化导致的卫星位置变化，即

根据式（5）可以得到广播星历正常情况下卫星位置对地球静止轨道（geostationary Earth orbit,GEO）、倾斜地球同步轨道（inclined geosynchronous orbits, IGSO）及中圆地球轨道（medium Earth orbit, MEO）卫星每个广播星历参数变化的灵敏度曲线，如图1所示。

图1 GEO、IGSO、MEO卫星位置误差对星历参数变化的灵敏度曲线

图1中从左至右依次为GEO、IGSO、MEO卫星位置误差对每个广播星历参数变化的灵敏度曲线。其中M0和ω参数的灵敏度曲线重合，由于部分广播星历参数卫星位置误差较大，导致部分位置误差较小的Cic和Cis灵敏度曲线在图1中并不明显，GEO、IGSO、MEO卫星参数Cic和Cis的灵敏度曲线如图2所示。

图2 GEO、IGSO、MEO卫星Cic和Cis的灵敏度曲线

由图1和图2可以看出，GEO、IGSO、MEO卫星大部分广播星历参数的灵敏度曲线存在一定的周期性。对于GEO卫星，参数M0、0Ω和ω的变化对卫星位置影响较大；对于IGSO卫星，参数和ω的变化对卫星位置影响较大；对于MEO卫星，参数M0、a、0Ω、0i和ω的变化对卫星位置影响较大。而参数Cuc、Cus、Crc、Crs、Cic和Cis的变化，对GEO、IGSO、MEO三类卫星位置的影响都比较小，其中Cic和Cis对卫星位置的影响最小，其影响量均不超过3 m。因此，在构建广播星历A1型故障监测模型时，可忽略Cic和Cis对卫星位置的影响，从而引入2个载波相位观测的整周模糊度。

2 A1型故障完好性监测

广播星历以及根据卫星信号获得的伪距观测量或载波相位测量值，均可用于计算用户位置。因此，当广播星历发生故障后，将对定位结果产生很大影响。差分是消除许多误差源的有力手段。GBAS是基于差分技术的局域卫星增强系统，将地面基准站全部可视卫星的差分修正量播发给机载用户，机载用户利用接收的差分修正量对自身观测量进行校正，采用修正后的观测量进行定位解算，从而能够提高用户定位精度，这有助于提升飞机精密进近的安全性能[9]。GBAS使用多个具有高性能、多径限制天线的接收机，以满足机载用户对导航精度、完好性、连续性和可用性的严格要求。

本文基于短基线参考接收机的地面设施，根据系统的完好性风险要求推导了MDE：采用伪距和载波相位观测量联合的方式测量轨道误差，构建关于位置误差的观测方程；结合广播星历参数灵敏度分析中的位置误差函数，将关于位置误差的观测方程转化为关于广播星历参数偏差的观测方程；根据最小二乘法得到卫星位置误差的协方差矩阵；最后计算MDE并与检测门限比较，实现广播星历A1型故障的完好性监测。

2.1 完好性风险约束

垂直保护水平是垂直方向的定位误差边界，与满足完好性风险要求的漏检率有关，飞机用户不断生成垂直保护水平来反映卫星位置误差[10]。当垂直保护水平超过GBAS最低运行性能标准中规定的固定警报限值时，将会在飞机上发出警报[11]。

垂直方向的卫星位置误差[12]为

式中：Svert,i为加权最小二乘投影矩阵中第i颗卫星的垂直系数；ie为卫星视距方向单位向量；x为飞机天线与地面站天线的矢量位移；iρ为卫星i到地面站接收机的标量范围；νρ为误差源投影到用户位置域的差分测距误差。

假设存在检验统计量te，则检测门限eT可以定义为

式中：tσ为检验统计量误差tν的标准差；kFFA为确保无故障报警概率的参数。假设tν服从零均值的正态分布，则有

结合式（6）至式（9），可以得到

式中jσ为νp的标准差，即第j颗卫星无故障差分测距误差的标准差，且有

定义1个基于漏检率的漏检检测参数kmd，则垂直保护水平的公式为

MDE的计算公式为

地面站向每个卫星的空中用户广播的去相关参数iP与EMD的关系[13]为

所以飞机在星历失效下的垂直保护水平VPL可以表示为

从式（15）可以看出，MDE的值越小，则垂直保护水平越小。MDE为广播星历故障的检测能力，全球定位系统（global positioning system,GPS）广播星历A1型故障监测方法中，MDE的检测门限为3 500 m[14]，由于BDS目前的发展还没有GPS成熟，在定位精度等方面劣于GPS，采用GPS的检测门限能够满足BDS的完好性风险要求[15]，所以文本采用3 500 m作为BDS的MDE检测门限。

星历异常情况下，检验统计量具有非中心的卡方分布，该分布的最小非中心性参数λ与完好性风险要求的漏检率PMD保持一致[16]，相关公式为

2 × 10-7的漏检率可以满足GBAS的完好性风险要求，根据式（16）可得到sλ的值。星历异常情况下，最小非中心性参数的值小于检验统计量的值，关系为

式中：δer为星历异常下的卫星位置误差；Cδr为卫星位置误差的协方差矩阵。

根据公式（17）可以得到

式中 eigCδr为协方差矩阵的特征值。

MDE的计算公式为

若根据式（19）计算的卫星MDE值均小于BDS的MDE检测门限，则满足GBAS的完好性风险要求。

2.2 误差观测模型构建

单差观测方程是在同一时刻、2台接收机对相同卫星的观测值相减得到的观测方程，即站间差分。对于在同一时刻的同一卫星，卫星具有的卫星钟差大小完全相同；因此站间差分可以用来消除和卫星有关的误差项，且在短基线情况下也可以通过站间差分的方式消除大气延迟误差[17]。对于差分用户，定位精度主要受正交于视距方向的轨道误差影响，须采用精度较高的双差载波相位来测量正交于视距方向的轨道误差。投影到视距方向的轨道误差直接导致伪距测量误差，对于视距方向的轨道误差，可采用伪距差分测量，在保证定位精度的同时降低模糊度解算难度。

采用伪距差分zkρ测量卫星视距方向ek的轨道误差，采用双差载波相位zkφ测量正交于视距方向的轨道误差，并结合广播星历参数灵敏度分析中的位置误差函数，即得到观测方程为

式中：l是基线长度；xb是接收机基线方向单位向量；I是大小为3×3的单位矩阵；N是整固模糊度；sλ是信号波长；Rk是卫星到接收机的距离；kφυ是载波相位的测量误差；kρυ是伪距测量误差。测量误差kρυ和kφυ分别服从以下分布：

式中：k为不同时刻，k=1,2, …,24；σρ是伪距差分测量误差的标准偏差；φσ是差分载波相位测量误差的标准偏差。

将式（20）表示为

式中：Z为所选取地面站的观测量矩阵；H为观测矩阵；υ为观测误差矩阵。

根据式（20）可知观测矩阵H为

式中O为矩阵元素为0的矩阵。

根据卫星位置误差对广播星历参数的灵敏度分析可知，Cic和Cis对GEO、IGSO、MEO卫星位置的影响最小，所以在构建A1型故障监测模型时，可忽略这2个参数。下文中的δaP不包含Cic和Cis，加入2个接收机的未知整周模糊度。

测量误差υ服从分布

式中

式中ρV为伪距差分测量误差矩阵。

由于kφυ主要由接收机白噪声和多路径造成，且多路径噪声的时间长度通常为1 min，采样时间大于2 min，上述多路径噪声通常不相关[18]。而伪距差分测量误差随时间的变化可以描述为一阶高斯马尔可夫过程为

式中：Δt是采样时间间隔；τ是伪距测量误差的时间常数。

通过加权最小二乘法得到广播星历参数偏差的最佳估计值及其协方差矩阵分别为：

式中W为加权矩阵。

结合式（4）可得到卫星位置误差的协方差矩阵为

将式（30）代入式（19）中可以计算卫星的MDE。MDE的值与伪距测量误差的时间常数τ、伪距差分测量误差的标准差ρσ、双差载波相位测量误差的标准差φσ有关，仅改变上述参数中其中的某个参数，令其他参数的值保持不变，分析MDE对每个参数的灵敏度。若计算的MDE值均小于MDE检测门限，则满足系统的完好性风险要求，完成A1型故障监测。

3 实验与结果分析

本文选取香港大地测量中心的HKKS、HKSS、HKWS 三个观测站，2020-07-01 BDS-2、BDS-3 C01～ C30卫星的观测数据和广播星历进行A1型故障监测。本文选择GPS局域增强系统原型中，多径抑制天线性能中规定的σρ=10m、σφ=0.05m作为初始值。

首先改变τ的值，分别选取τ的值为3、12和24 h计算每个卫星的MDE值，部分卫星MDE实验结果如图3所示。

图3 τ分别为3、12、24 h时计算的MDE

从图3可以看出，τ取不同值时，计算的MDE结果近似，即τ的大小对MDE的结果值影响很小，且每个卫星的MDE值均小于3 500 m的检测门限，满足GBAS的完好性风险要求，所以在下面的实验中，选取τ=12 h。

其次改变ρσ的值，分别选取ρσ的值为6、20和1×1010m，计算各个卫星的MDE，部分卫星结果如图4所示。根据图4可知，当ρσ值小于20 m时，ρσ的改变对MDE的影响很小，每个卫星的MDE值均小于3 500 m，满足GBAS的完好性风险要求。

图4 ρσ分别为6、20、1× 1010 m时计算的MDE

最后改变φσ的值，分别选取φσ的值为0.04和0.06 m时计算MDE，部分卫星结果如图5所示。由图5可以看出，当φσ的值为0.06 m时，C20卫星某时刻的MDE结果值略大于3 500 m，即φσ的大小对MDE影响较大。所以，当φσ小于0.06 m时，满足GBAS的完好性风险要求。

图5 φσ分别为0.04、0.06 m时计算的MDE

为验证上述监测方法的有效性，选取上述初始值进行实验，在0时刻引入C01号卫星广播星历参数误差后，计算的MDE结果趋近于无穷大，远大于3 500 m的MDE检测门限。综上所述，选定天线性能要求的初始值后，仅改变影响MDE的某个参数后计算得到的MDE数值均小于检测门限，采用伪距和载波相位观测量联合的方式测量轨道误差，与卫星位置误差对广播星历参数变化的灵敏度结合，构建A1型故障监测模型是有效的。