王韦舒，上官伟,2,3，刘 江,2,3，姜 维,2,3
(1.北京交通大学 电子信息工程学院，北京 100044；
2.北京交通大学 轨道交通控制与安全国家重点实验室，北京 100044；
3.北京市电磁兼容与卫星导航工程技术研究中心，北京 100044)
随着卫星导航系统的发展完善，以卫星定位为主要特征的下一代列车运行控制系统已成为新的重点研究方向。卫星定位技术的引入，旨在实现列车自主感知，减少对轨道电路或应答器等地面设备的依赖，提高车载设备的智能化水平，降低运营维护成本[1-2]。通过卫星导航技术进行列车位置估计和列车完整性检查，有助于推动既有列车运行控制系统中闭塞制式的转变，由固定闭塞向移动闭塞的跨越将有利于实现列车动态高效追踪，进一步缩短列车间隔，从而提高运输效率。在列车动态追踪运行中，精准可靠的列车位置将是列控系统生成列车许可和速度控制曲线的重要前提，直接影响着列车的安全运行。
目前，扩展卡尔曼滤波(Extended Kalman Filter, EKF)是实施卫星定位计算的常用方法，其估计准则是最小均方误差(Minimum Mean Square Error, MMSE)准则，在量测噪声服从高斯分布的条件下能够实现列车位置状态的最优估计。但卫星观测量容易受到钟跳、轨道参数建模误差、多径效应等的影响产生阶跃故障或者瞬时粗大偏差[3]，量测噪声不再服从高斯分布。如果不对异常量测进行处理，将会造成列车定位性能下降，影响列车控制决策。针对MMSE准则下残差二次型对于故障的敏感性问题，故障检测与排除(Fault Detection and Exclude, FDE)和鲁棒滤波估计是两种典型的解决方案。FDE的核心是排除所有可能会对定位结果产生不利影响的量测，即尽最大可能检测出异常量测并剔除。当可见卫星数量足够多时，FDE算法能够保障定位估计过程的稳定性。但在受限观测环境或可见卫星较少条件下，排除故障量测所导致的空间几何结构的改变，同样会直接影响到定位解算性能[4]。不同于FDE，鲁棒估计方法能够对未知统计特性的有界噪声进行处理，重在通过权函数降低异常数据的不利影响。基于Huber的鲁棒滤波已经应用于列车定位中用于抑制环境中多源噪声对定位性能的影响[5-6]。但该方法受到权函数的约束，对于明显异常的量测仍旧会赋予其一定的权值，使得滤波精度由于错误的量测而出现下降；
并且对于近似正常的量测其权值较小，未能充分利用量测。
目前，基于信息理论学习的优化准则逐渐受到研究人员广泛的关注[7]。其中，相关熵描述了两个随机变量之间的广义相似性，相关熵越大表示估计值与真实值越接近[8]。与MMSE不同的是，相关熵能够捕获更高阶的统计量，提取由于高斯性偏离引起的各种信息，包含了二阶统计量没有的大量有效信息。由相关熵特性产生的最大相关熵准则作为一种局部相似性度量，对量测噪声非零均值、非高斯或者有较大异常值等情况都具有较强的抑制能力，开始逐渐引入到滤波理论中，增强既有滤波算法的鲁棒性[9]。在此思路下，本文对传统列车定位解算进行改进，推导最大相关熵准则下的扩展卡尔曼滤波，通过重建EKF的量测噪声矩阵，实现对滤波增益中先验信息的权重调整；
在此基础上，结合M估计理论提出一种自适应核宽度的扩展卡尔曼滤波定位方法，并利用现场数据对该方法进行验证和分析。
要实现列车自主定位，仅仅依靠卫星导航定位系统难以满足安全需求。通常采用融合算法将不同类型的定位传感器信息进行联合处理，如惯性导航[10]、轮速传感器[11]、激光雷达[12]、光电传感器[13]等。其中惯性导航系统(Inertial Navigation System, INS)作为一种自主式导航系统，不依赖于外部信息，具有采样频率高、隐蔽性强、输出信息丰富等优点，已成为解决自主定位问题的重要选择[14]。根据信息融合的深度，可选择松组合或者紧组合的方式将卫星导航系统与惯性器件测量的姿态、加速度相融合。与松组合相比，紧组合在伪距域将卫星原始观测数据和INS计算的预测伪距/伪距率相结合，能够克服可观测卫星数量的约束，在观测条件欠佳下仍能实现有效定位。GNSS/INS紧组合处理流程见图1。
图1 GNSS/INS紧组合处理流程
列车组合定位系统的性能主要取决于数据融合算法。对于一个非线性系统而言，其状态方程和量测方程可表示为
( 1 )
式中：xk为k时刻列车状态向量；
zk为k时刻量测向量；
f(·)、h(·)分别为二阶可微的状态方程和观测方程；
wk-1、vk分别为相互独立且服从零均值高斯分布的过程噪声、量测噪声，相应的协方差矩阵分别为Qk-1、Rk，系统的先验误差协方差矩阵记为Pk-1。
通过最小均方误差准则，损失函数JMMSE可表示为
( 2 )
求解argminJ(xk)能够获取系统状态量统计意义上的最优估计，其前提在于系统噪声和量测噪声为服从零均值正态分布的高斯白噪声。然而在实际应用中，由于运行环境等不确定因素的影响有可能导致量测故障或产生粗大偏差，此时如果继续使用MMSE准则的滤波算法，将使得滤波估计性能显著下降甚至滤波发散，造成较大定位偏差。因此，如何设计一个鲁棒性强的滤波算法来保证列车定位性能是需要解决的关键问题。
2.1 最大相关熵基本原理
相关熵是两个随机变量间的广义相似度量。对于给定的随机变量X和Y，其联合分布函数为FX,Y(X,Y)，则熵V(X,Y)可以表示为
V(X,Y)=E[φ(X,Y)]=∬φ(X,Y)dFX,Y(X,Y)
( 3 )
式中：E为计算期望；
φ(X,Y)为核函数，通常选择高斯核函数，表示为
( 4 )
其中，e=X-Y；
σ为核宽度。
在实际场景中，只有有限的数据是已知的，并且联合分布函数FXY(X,Y)大多数情况是未知的。在这种情况下，相关熵通过样本均值估计得到，即
( 5 )
式中：e(i)=X(i)-Y(i)；
N为抽取的总样本数。
将高斯核函数按照泰勒级数展开，可得熵为
( 6 )
由式(6)可以看出，熵为变量X-Y各阶项的加权和；
核宽度为二阶及更高阶项的权重参数；
核宽度越大，高阶项的影响越小，熵值大小将主要由二阶展开项所决定。
最大相关熵损失函数JMCC可以写为
( 7 )
假设W为自适应系统中需要求解的参数向量，X(i)和Y(i)分别为系统模型输出和期望响应，则通过最大熵得到的最优解可以表示为
( 8 )
2.2 基于相关熵的自适应扩展卡尔曼滤波
基于最大相关熵准则的特性，重新推导新准则下的扩展卡尔曼滤波。将式(1)线性化可得
xk=Fk-1xk-1+wk-1
( 9 )
zk=Hkxk+vk
(10)
将式(9)、式(10)重构成非线性递推方程的形式，可表示为
(11)
可进一步得到
(12)
(13)
Lk=Γkxk+ek
(14)
为解决最小均方误差准则对于故障的敏感性问题，不再采用式(2)作为代价函数，而是引入最大相关熵准则降低伪距量测故障所带来的不利影响。重新构建的代价函数为
(15)
状态量的最优解可表示为
xk=argmaxJMCC(xk)
(16)
通过对式(16)进行一阶求导得
(17)
定义权重矩阵为
(18)
CP,k=diag(Gσ(e1,k),…,Gσ(em,k))
(19)
CR,k=diag(Gσ(em+1,k),…,Gσ(en,k))
(20)
(21)
将式(13)、式(18)代入式(21)可得
(22)
修改后的先验误差协方差矩阵、量测噪声协方差矩阵为
(23)
(24)
(25)
(26)
(27)
上述过程推导了基于固定核宽度的相关熵扩展卡尔曼滤波(Fixed Maximum Correntropy Extended Kalman Filter, FMCEKF)，选用高斯核函数代替MMES估计，能够获取更高阶的量测统计信息，降低含故障量测信息的权重，从而保证滤波估计性能的稳健性。当量测中夹杂有较大的伪距偏差时，滤波器的ek将大于正常值，通过核函数所获得的权重值将非常小。不同核宽度d下权函数变化情况见图2。
图2 不同核宽度下权函数变化情况
由图2可以进一步获知，不同核宽度下对于残差的权重分配差异较大。当核宽度较小时，相同量测误差下所设定的权重更小，更有利于增强系统的稳健性。因此，如何选择合适的核宽度将影响到滤波算法对异常点的削弱能力。
受到鲁棒估计思想的启发，将Pseudo-Huber权函数引入到FMCEKF来调整核宽度大小，提出一种自适应核宽度的相关熵扩展卡尔曼滤波(Adaptive Maximum Correntropy Extended Kalman Filter, AMCEKF)。根据量测质量的不同而设定可变的核宽度，最大程度上利用观测到的信息，并削弱异常数据对估计结果的影响。
(28)
第i颗卫星观测量对应的核宽度di,k可表示为
(29)
(30)
式中：ri为第i颗卫星的新息，由k时刻量测信息和状态预测值之差得到；
ξr,k为k时刻新息标准差；
ωi,k为第i颗卫星的标准化新息；
τ为调节因子，通常取1.210 7[13]。当一个或者多个量测受到外界扰动发生异常时，其标准化新息将显著增大，所对应的核宽度将选择一个相对较小的值，进而在增益矩阵中分配给其更小的权重比例，以提高状态估计的鲁棒性。
2.3 基于伪距/伪距率的组合定位系统
2.3.1 状态方程
基于伪距、伪距率的GNSS/INS组合状态方程由INS误差方程和GNSS误差方程组成，可以表示为
xk=Fk-1xk-1+wk-1
(31)
(32)
式中：Fk-1为状态转移矩阵，由k-1时刻GNSS、INS的误差传播模型组成。xk-1为系统的17维状态估计量。
由惯性导航系统的姿态、速度、位置、加速度计零偏、陀螺漂移误差，以及接收机的钟差、钟漂误差组成，具体表示为
X=[φxφyφzδvxδvyδvzδxδyδzεxεyεz
(33)
式中：(φx,φy,φz)为姿态误差；
(δvx,δvy,δvz)为速度误差；
(δx,δy,δz)为位置误差；
(εx,εy,εz)为陀螺仪的三轴漂移误差；
(x,y,z)为加速度计的三轴随机误差；
分别为等效钟差距离误差和等效钟漂距离误差。
2.3.2 量测方程
系统的量测方程可以表示为
zk=Hkxk+vk
(34)
式中：Hk为组合系统的量测矩阵；
zk为系统量测值，由惯性导航推算得到的预测伪距、伪距率，以及定位接收机测量得到的伪距、伪距率所确定，表示为
(35)
(36)
同理，接收机与第i颗卫星的测量伪距为
δb+η+δe
(37)
式中：δb为接收机钟差造成的距离误差；
η为卫星时钟偏差、电离层误差、对流层误差的总和，能够通过近似模型加以估计；
δe为残留的随机误差。
INS估计的伪距率可以表示为
(38)
(39)
同理，可将GNSS伪距率建模为
(40)
综合式(36)～式(40)可得量测矩阵为
(41)
2.4 基于AMCEKF的组合定位算法
本文所提出的算法以EKF为基础，根据最大相关熵准则和Pseudo-Huber权函数重新构造量测噪声矩阵，实现自适应调节量测所占配比，充分利用所观测到的量测信息，并降低故障量测对状态估计的影响。具体步骤如下：
Step1初始化滤波器参数，包括状态估计量X、系统噪声、量测噪声等。
Step2构建新息，根据式(35)～式(40)计算伪距和伪距率量测误差，确定量测矩阵，并对式(13)采用Cholesky分解误差协方差矩阵，进而计算系统新息。
Step3权重更新，根据式(28)对每一个量测计算特定核宽度，并结合式(18)重构权重矩阵。
Step4确定增益矩阵，根据式(23)～式(25)，由新的权重矩阵更新先验误差协方差矩阵和量测噪声协方差矩阵，并计算各状态的增益。
Step5滤波输出，基于增益矩阵更新状态估计量和误差协方差矩阵。
本文采用2019年于拉萨—林芝铁路采集的列车行驶数据对所提出的定位算法进行测试和验证。主要分析瞬时故障、单颗卫星阶跃故障、多颗卫星阶跃故障场景下算法的实际性能。
所选实验数据位于日喀则—拉萨，列车上安装诺瓦泰SPAN-FSAS分体式组合导航定位系统，用于采集卫星原始观测量和惯性导航数据，数据更新频率分别为1、200 Hz。在列车从日喀则向拉萨行驶过程中，卫星观测条件良好，水平精度因子在0.6左右，可见卫星数不低于10颗。
3.1 瞬时故障下定位性能验证
为测试组合定位算法在瞬时故障下的定位性能，需仿真故障卫星所对应的量测数据。选择PRN 26号卫星，在原始观测量中每隔50 s加入50 m的伪距偏差。分别用EKF、FMCEKF、AMCEKF算法对列车进行定位解算，瞬时故障下位置误差见图3。
图3 瞬时故障下位置误差
从图3可以看出：EKF对于瞬时故障异常敏感，滤波估计始终处于跳变和收敛过程中，最大定位误差达到了4.8 m。对于FMCEKF而言，不同的核宽度下定位结果差异较大；
当核宽度选择5.0时，在故障发生时会发生突变，定位结果与EKF相似，位置误差在5.0 m左右波动；
但不同于EKF，在无故障后FMCEKF能够快速收敛；
进一步当核宽度选取较小的值0.5时，大部分情况下定位误差趋于稳定，能够较好地抑制脉冲异常的影响。相比EKF和FMCEKF，AMCEKF估计结果更加平稳，均值为1.675 m，能够最大程度上减少瞬时故障对于定位性能的影响，对所注入的伪距偏差具有较低的敏感性。
与FMCEKF相比，AMCEKF具有较强鲁棒性的原因在于核宽度不再是恒定不变的值，而是随观测情况进行调整。瞬时故障中不同卫星量测的核宽度变化过程见图4。当PRN 26存在故障时，在自适应调整策略下，核宽度均被设置成更小的值为0.32；
无故障时，核宽度在0.85～0.96之间波动。这取决于不同时刻各个卫星的伪距新息变化。新息越大，所确定的核宽度越小，从而获得较小的相对熵，以达到降低异常影响的目的。
图4 瞬时故障中不同卫星量测的核宽度变化过程
为了对估计误差对比情况进行量化分析，选择均方根误差RMSE和距离均方根误差DRMS评价定位结果[16]，瞬时故障下不同算法定位误差统计结果见表1。从均方根误差来看，与EKF相比，核宽度=0.5的FMCEKF算法在北向、东向、地向上分别增加了26.8%、41.7%、13.3%；
AMCEKF则增加了52.6%、78%、23%；
对水平误差，FMCEKF和AMCEKF分别增加了29.4%和56.5%。
表1 瞬时故障下不同算法定位误差 m
3.2 阶跃故障下定位性能验证
为进一步检验所提算法的实际性能，对单颗卫星阶跃故障和多颗卫星阶跃故障两种场景下的定位结果进行对比分析。单故障时，选择PRN 25，在100～150 s注入50 m的伪距偏差；
多故障时，在350～400 s对PRN 25、26号卫星注入30 m的伪距偏差。阶跃故障下不同算法对列车位置估计的结果见图5。
图5 阶跃故障下位置误差
从图5可以得出：单星故障时，EKF估计误差会随着时间不断增大，北向、东向、地向最大误差分别达到11、15、45 m。对于FMCEKF而言，如果核宽度选择不当，其对异常的敏感性将会增强，如图5中核宽度为5.0时所示；
当核宽度为0.8时，在正常观测范围时，其估计性能与EKF基本一致；
但在异常区域时，仍能保持稳定的精度水平，受到故障值的影响较小。AMCEKF的估计误差略微有所增大，但北向、东向误差最大不超过2 m，估计准确性和鲁棒性均优于FMCEKF和EKF。多故障下，FMCEKF在核宽度为5.0时，对于东向、地向上的估计性能退化较为严重；
而核宽度为0.8时，仍能够保持稳定的估计性能。AMCEKF虽然略微有所浮动，但整体而言，无论是在正常还是故障条件下均能够实现更高的定位精度。
阶跃故障下不同卫星量测的核宽度变化见图6，示出了实验过程中PRN 22、25、26、31号卫星所设定的核宽度变化情况。由图6可知，在100～150、350～400 s时刻，AMCEKF对故障卫星(PRN 25、26)的核宽度自适应调整到了0.3～0.4，通过较小的核宽度获得增益矩阵中较小的权重，以减弱故障量测的影响。对于正常卫星(PRN 22、31)的量测，AMCEKF也能够根据其观测质量动态选择相应的核宽度，从而实现充分利用有效量测的目的。
图6 阶跃故障下不同卫星量测的核宽度变化
阶跃故障下不同算法定位误差见表2。从表2可以看出，EKF算法由于无法抑制故障值的影响，其水平位置的均方根误差较大；
d=5.0时的FMCEKF在北向、东向的RMSE分别增加12.4%、7.9%，d=0.8时则增加了20.6%、62.9%；
从DRMS来看，AMCEKF相比EKF、FMCEKF(d=0.8)，分别增加了61.2%、35%。总体可得，FMCEKF的估计性能受核宽度影响较大，选择一个合适的核宽度将直接决定着定位算法的鲁棒性能。AMCEKF能够合理选择相应的核大小，无论是在瞬时故障还是阶跃故障情况下，都能表现出较好的鲁棒性。
表2 阶跃故障下不同算法定位误差 m
以卫星导航技术为核心的列车位置感知已经成为下一代列控系统实际工程应用中必须攻克的关键技术。面对复杂多变的环境，如何保障外界干扰条件下列车位置估计的准确性和稳定性，是列车组合定位融合估计中面临的重大难题。
(1)针对外界异变条件下滤波精度下降，以及固定核宽度对估计性能的约束问题，将鲁棒估计思想和信息理论学习相结合，融入到组合定位系统的状态估计中，设计了基于最大相关熵准则的自适应扩展卡尔曼滤波定位方法。
(2)以拉萨—林芝铁路某区段实测数据对所提出的算法进行了实际验证。结果表明，所提方法在瞬时故障、阶跃单故障和阶跃多故障3种不同的场景中，均能够有效抑制故障检测对列车定位估计准确性的不利影响，能够随时变观测质量动态调整各个卫星量测在定位估计中所占的比重，有效解决了故障量测带来的估计失准问题，提高了整体定位精度，具有较强的鲁棒性和适应性，体现出该方法在未来实际应用中的巨大潜力。
猜你喜欢 伪距阶跃宽度 BDS-3载波相位平滑伪距单点定位性能分析地理空间信息(2022年10期)2022-10-31综采液压支架大流量水压比例阀阶跃响应特性研究＊中国科技纵横(2021年22期)2022-01-08BDS 三频数据周跳探测与修复方法研究导航定位学报(2021年5期)2021-10-13北斗二号伪距偏差特性分析及其对定位的影响导航定位与授时(2021年3期)2021-05-18两种伪距定位精度分析及计算程序的实现导航定位学报(2021年2期)2021-04-22浅析Lexus车系λ传感器工作原理及空燃比控制策略(三)汽车维修与保养(2021年11期)2021-03-11阶跃响应在系统超调抑制中的运用电子技术与软件工程(2020年17期)2020-02-02基于SRD的梳状谱发生器设计中文信息(2016年7期)2016-05-14孩子成长中，对宽度的追求更重要人生十六七(2015年5期)2015-02-28你有“马屁股的宽度”吗？销售与市场·管理版(2009年21期)2009-09-03