基于抗差方差分量估计的码伪距/多普勒联合测速方法及其应用与流程

本发明属于gnss(全球导航卫星系统)定位与导航技术领域，特别涉及基于抗差方差分量估计的码伪距/多普勒联合测定接收机的运动速度。

背景技术：

速度是航空飞行、智能导航、无人驾驶、海上交通等领域的重要运动参数之一，gnss测速的高精度、实时、廉价等特点使得其在这些领域得到了广泛的应用。国外学者通过差分数据模拟测试，得到了测速误差与载体的运动速度和加速度变化率成正比的关系，误差范围从亚毫米每秒到几米每秒不等。国内学者也曾对伪距、载波、多普勒及其差分观测值的测速模型和测速精度进行了系统的剖析和评价。详细地，伪距、载波和多普勒观测值均可测定接收机速度，直接利用gnss定位结果进行位置差分就能求取平均速度，但这非常依赖于定位结果的精度，实用性不强。另外，通过历元间伪距差分可求得差分时间内的平均速度，但受码伪距精度影响，测速误差较大；利用历元间相位差分也能求得差分时间内的平均速度，并且精度较高，但载波相位易受干扰而周跳频繁，甚至无输出，可靠性较低；多普勒观测值可以获得高精度的瞬时速度，对轨道误差、接收机钟差、大气误差不敏感，不存在周跳，伪距单点定位误差对测速精度的影响也在毫米级别，是一类比较实用可靠的观测值。

纵观有诸多种观测信息或测速方法，但是如何融合不同类别的观测信息，最大化地实现高精度稳定测速是工程应用中面临的迫切问题。

技术实现要素：

技术问题：针对上述现有技术，借鉴多源信息融合的思想，提出采用单频码伪距和多普勒两类不同精度的观测值联合解算接收机速度的方法，同时基于抗差估计与方差分量估计相结合的方式，探索建立更为合理有效的随机模型，以最大化地发挥联合模型的效果。

技术方案：基于抗差方差分量估计的码伪距/多普勒联合测速方法，首先采用码伪距、多普勒、载噪比、卫星高度角的观测信息分别建立测速的函数模型和随机模型；然后采用最小二乘抗差估计的方法获得能够抑制粗差影响的抗差随机模型；最后将码伪距和多普勒的测速函数模型与抗差随机模型组成联合测速模型，采用方差分量估计的方法迭代求解卫星接收机的运动速度。

进一步地，基于抗差方差分量估计的码伪距/多普勒联合测速方法包括如下具体步骤：

步骤1)，采用码伪距、多普勒、载噪比、卫星高度角观测信息分别建立测速的函数模型和随机模型，包括如下具体步骤：

a)，通过历元间差分的方式求得码伪距和载波的变化率，如式(1)所示：

式中，ρ和分别为码伪距及其变化率，和分别为载波相位及其变化率，d为多普勒观测值，δt为差分时间，其中各项的下标k、k+1分别为第k时刻和第k+1时刻；

b)，利用码伪距和多普勒观测值建立测速的函数模型，一般地gnss观测方程如式(2)所示：

式中，ρ为接收机码伪距观测值，为接收机载波观测值，n0为整周模糊度，λ为载波波长，r为接收机与卫星间的真实几何距离，c为光速，δt为钟差，δρion为电离层延迟，δρtrop为对流层延迟，ερ和为包含轨道误差、多路径效应、观测噪声在内的其他误差项，其中各项的上标s和下标r分别代表卫星和接收机。

根据gnss观测方程对时间t求导并进行线性化，如式(3)所示：

式中，各项表示随时间t的变化率，且

式中，r^s、为ecef(earth-centered,earth-fixed)坐标系下的卫星位置、速度列向量，rr、为ecef坐标系下接收机位置、速度列向量，r0为由接收机概略坐标rr0求得的接收机与卫星间的几何距离，为接收机概略速度，x、y和z为位置矢量在ecef坐标系下的投影，vx、vy和vz为速度矢量在ecef坐标系下的投影。联合公式(1)、(3)和(4)可得到相应的测速函数模型；

c)，利用载噪比、卫星高度角等观测信息建立测速的随机模型，如式(5)所示：

式中，σ为协方差，e表示卫星高度角，下标i表示卫星编号，s为缩放因子，常数项a0、a1和e0由表1定义：

其中，缩放因子s由卫星载噪比定义，如式(6)所示：

式中，c/n0表示卫星载噪比，int(*)为取整函数。由公式(5)、(6)和表1可确定测速的先验随机模型。

步骤2)，采用最小二乘抗差估计的方法获得能够抑制粗差影响的抗差随机模型，包括如下具体步骤：

a)，通过最小二乘估计求得残差向量和相应的协因数，如式(7)所示：

式中，b为设计矩阵，p为先验权，l为观测值向量，q为观测值协因数，待估参数，v为残差向量，qvv为残差协因数。

b)，通过周江文提出的iggiii等价权方案获得能够抑制粗差影响的抗差随机模型，如式(8)所示：

式中，为标准化残差，k0和k1为常量，一般k0∈[1.0～1.5]，k1∈[2.5～8.0]，为抗差等价权，下标i代表第i个观测值。

步骤3)，将码伪距和多普勒的测速函数模型与抗差随机模型组成联合测速模型，采用方差分量估计的方法迭代求解卫星接收机的运动速度，包括如下具体步骤：

a)，将码伪距和多普勒的测速函数模型与抗差随机模型组成联合测速模型，如式(9)所示：

式中，各项下标1和2分别代表码伪距和多普勒，n和w分别代表的相应组合项；

b)，采用方差分量估计的方法迭代调整各类观测值的权重，如式(10)所示：

式中，tr(*)表示矩阵求迹，e(*)表示取期望，n为观测值个数，为单位权方差的估计值，v为残差向量，为抗差等价权，各项的下标1和2分别代表码伪距和多普勒。通过求解方程组(10)式，得到单位权方差的估值代入式(11)获得调整后的各类观测值权重：

式中，表示单次迭代调整得到的权重，c为常数，可固定选取中任一个，各项的下标1和2分别代表码伪距和多普勒。

c)，利用经方差分量估计调整后的联合测速模型求解卫星接收机的运动速度，如式(12)所示：

式中，为估计的速度矢量，上标-1代表矩阵求逆运算；

重复步骤a)、b)、c)直到各类单位权方差的估值相等或经假设检验其相等为止，即可求得接收机最终的运动速度。

进一步地，当进行步骤3)时，若经过4～5次循环迭代时仍无法满足单位权方差估值相等或通过假设检验，则跳出循环，直接采用抗差解作为最终的解算结果。

此外，上述引用了周江文，该人为大地测量学家，iggiii方案是周江文于1989年根据测量误差有界性提出来的第3套抗差估计方法，即抗差的等价权由公式(8)获得。这是抗差估计领域普遍采用且抗差效果较好的抗差估计方法，已被业界广泛应用和推崇，成为了公知的实用抗差方法。公式(8)是其提出的各套抗差方法当中运用最广的一套方法中的核心公式。本发明也是采用其抗差方法进行粗差的抑制。

有益效果：本发明所提出的一种基于抗差方差分量估计的码伪距/多普勒联合测速方法，通过利用码伪距和多普勒观测信息联合解算接收机速度，突破目前单频卫星接收机仅用单一观测值进行测速的局限；针对联合测速模型易受粗差等不良观测值影响的问题，采用等价权抗差估计原理抑制同类观测值内部粗差对模型的影响；对于融合多源信息进行联合测速时存在各类信息间不相容，精度不一致的问题，采用方差分量估计方法平衡不同类观测值之间的权重比例，自适应地调节联合测速的随机模型。使用本发明所提出的方法，可实现单频观测信息的充分有效利用，显著提升联合测速方法的鲁棒性、无偏性及稳定性。

附图说明

附图是用来提供对本发明的进一步理解，并且构成说明书的一部分，与下面的具体实施方式一起用于解释本发明，但并不构成对本发明的限制。在附图中：

图1为本发明一实施例中基于抗差方差分量估计的码伪距/多普勒联合测速方法的流程图；

图2为本发明一实施例中各解算方案测速真误差；

图3为本发明一实施例中单点定位n、e方向真误差；

图4为本发明一实施例中各解算方案x轴测速真误差；

图5为本发明一实施例中各解算方案合成速度误差；

具体实施方式

以下结合附图对本发明的具体实施方式进行详细说明。应当理解的是，此处所描述的具体实施方式仅用于说明和解释本发明，并不用于限制本发明。

一种基于抗差方差分量估计的码伪距/多普勒联合测速方法，使用码伪距、多普勒、载噪比、卫星高度角等观测信息分别建立测速的函数模型和随机模型，采用最小二乘抗差估计的方法获得能够抑制粗差影响的抗差随机模型，将码伪距和多普勒的测速函数模型与抗差随机模型组成联合测速模型，基于方差分量估计的方法迭代求解卫星接收机的运动速度。

在gnss测速解算当中，是通过建立观测方程并对方程中的未知参数(如位置、速度)进行求解(估计)，从而获得位置或速度信息。观测方程可分为两部分：观测值的函数模型与随机模型。函数模型刻画了观测值与待估计间的数学关系(即自然物理规律用数学方程表达，如位移＝速度*时间等，其中位移为各时间节点的观测值，速度为待估计的参数)；随机模型则描述观测值的噪声(即误差)大小及观测值间的相关性。要使待估计参数获得准确有效的求解(估计)，需要建立准确的函数模型和与之相匹配的随机模型。在实际解算当中，函数模型由确定的物理含义能够较为准确地获得，而由于观测设备和环境的随机的或非随机的综合影响，往往难以保证获得的观测值准确误差，也即随机模型难以准确获取建立。所以为建立准确的随机模型，本发明通过由观测值的主要特征信息(即信噪比、卫星高度角等)，联合建立先验的随机模型，然后通过抗差方差分量估计方法，利用函数模型和随机模型进行细致的调节，尽可能地使得函数模型与随机模型相匹配，且能准确地反映实际的测量情况。本发明主要特点在于观测值不仅仅采用单一的常规码伪距观测值，还融合了多普勒观测值，这是两种不同类型的观测值，观测值的物理含义和观测误差水平都不一致。所以不仅需要建立联合测速的函数模型，还需建立准确的随机模型，从而提高速度信息的估计精度，采用的方法就是具有抗差能力的方差分量估计方法。

实施例1

一种基于抗差方差分量估计的码伪距/多普勒联合测速方法包括如下具体步骤：

步骤1)，采用码伪距、多普勒、载噪比、卫星高度角等观测信息分别建立测速的函数模型和随机模型，包括如下具体步骤：

a)，通过历元间差分的方式求得码伪距和载波的变化率，如式(1)所示：

式中，ρ和分别为码伪距及其变化率，和分别为载波相位及其变化率，d为多普勒观测值，δt为差分时间，其中各项的下标k、k+1分别为第k时刻和第k+1时刻。

b)，利用码伪距和多普勒观测值建立测速的函数模型，一般地gnss观测方程如式(2)所示：

式中，ρ为接收机码伪距观测值，为接收机载波相位观测值，n0为整周模糊度，λ为载波波长，r为接收机与卫星间的真实几何距离，c为光速，δt为钟差，δρion为电离层延迟，δρtrop为对流层延迟，ερ和为包含轨道误差、多路径效应、观测噪声在内的其他误差项，其中各项的上标s和下标r分别代表卫星和接收机。

根据gnss观测方程对时间t求导并进行线性化，如式(3)所示：

式中，各项表示随时间t的变化率，且

式中，r^s、为ecef(earth-centered,earth-fixed)坐标系下的卫星位置、速度列向量，rr、为ecef坐标系下接收机位置、速度列向量，r0为由接收机概略坐标rr0求得的接收机与卫星间的几何距离，为接收机概略速度，x、y和z为位置矢量在ecef坐标系下的投影，vx、vy和vz为速度矢量在ecef坐标系下的投影。联合公式(1)、(3)和(4)可得到相应的测速函数模型；

c)，利用载噪比、卫星高度角等观测信息建立测速的随机模型，如式(5)所示：

式中，σ为协方差，e表示卫星高度角，下标i表示卫星编号，s为缩放因子，常数项a0、a1和e0由表1定义：

其中，缩放因子s由卫星载噪比定义，如式(6)所示：

式中，c/n0表示卫星载噪比，int(*)为取整函数。由公式(5)、(6)和表1可确定测速的先验随机模型。

步骤2)，采用最小二乘抗差估计的方法获得能够抑制粗差影响的抗差随机模型，包括如下具体步骤：

a)，通过最小二乘估计求得残差向量和相应的协因数，如式(7)所示：

式中，b为设计矩阵，p为先验权，l为观测值向量，q为观测值协因数，待估参数，v为残差向量，qvv为残差协因数。

b)，通过周江文提出的iggiii等价权方案获得能够抑制粗差影响的抗差随机模型，如式(8)所示：

式中，为标准化残差，k0和k1为常量，一般k0∈[1.0～1.5]，k1∈[2.5～8.0]，为抗差等价权，下标i代表第i个观测值。

步骤3)，将码伪距和多普勒的测速函数模型与抗差随机模型组成联合测速模型，采用方差分量估计的方法迭代求解卫星接收机的运动速度，包括如下具体步骤：

a)，将码伪距和多普勒的测速函数模型与抗差随机模型组成联合测速模型，如式(9)所示：

式中，各项下标1和2分别代表码伪距和多普勒。

b)，采用方差分量估计的方法迭代调整各类观测值的权重，如式(10)所示：

式中，tr(*)表示矩阵求迹，e(*)表示取期望，n为观测值个数，为单位权方差的估计值，v为残差向量，为抗差等价权，各项的下标1和2分别代表码伪距和多普勒。通过求解方程组(10)式，得到单位权方差的估值代入式(11)获得调整后的各类观测值权重：

式中，表示单次迭代调整得到的权重，c为常数，可固定选取中任一个，各项的下标1和2分别代表码伪距和多普勒。

c)，利用经方差分量估计调整后的联合测速模型求解卫星接收机的运动速度，如式(12)所示：

式中，为估计的速度矢量，上标-1代表矩阵求逆运算；

重复步骤a)、b)、c)直到各类单位权方差的估值相等或经假设检验其相等为止，即可求得接收机最终的运动速度。其中，当进行步骤3)时，若经过4～5次循环迭代时仍无法满足单位权方差估值相等或通过假设检验，则跳出循环，直接采用抗差解作为最终的解算结果。

在本实施例中设定粗差判定阈值k0为1.2、k1为5.5，最大循环阈值为4。

实施例2

下面对基于抗差方差分量估计的码伪距/多普勒联合测速方法的应用场景进行举例说明：

为综合比较本发明所提方法的实际解算效果，设计了静态、静态模拟动态以及车载动态等三种测试方式，分4种解算方案进行比较：

方案a：最小二乘估计，即直接采用先验随机模型和联合测速的函数模型进行最小二乘平差解算。

方案b：方差分量估计，即在方案a的基础上直接采用常规方法分量估计的方法调节伪距和多普勒的先验随机模型。

方案c：最小二乘抗差估计，即分别对伪距和多普勒的先验权阵，采用抗差最小二乘估计的方式调节不良观测值的权重，然后组合成联合测速模型进行最小二乘平差解算。

方案d：抗差方差分量估计，即在方案c的基础上加入方差分量估计，从而形成抗差与方差分量相结合的估计方法对伪距和多普勒联合测速模型进行平差解算。

1)、静态测试

图2为各解算方案测速真误差图。选取的是澳洲curtinuniversitycut0站年积日为148d全天的静态观测数据，具体为trimblenetr9接收机采集的gpsl1和bdsb1上伪距和多普勒观测信息进行解算。常规最小二乘估计的速度真误差有较明显的波动，并且由于升降卫星等不良观测值的影响，造成速度真误差存在一定程度的跳动。通过抗差估计后，这种波动和跳跃得到了一定程度的抑制，另外常规的方差分量估计也能达到类似的效果，但是基于抗差估计与方差分量相结合的估计方法效果最佳。

下表1为各解算方案真误差统计表。各解算方案的均值中，抗差方差分量估计的速度均值最接近于速度真值零，方差分量估计和抗差估计依次次之，最小二乘估计最差。依据参数估计无偏性可知，抗差方差分量的速度估值无偏性最好，由此可认为抗差方差分量估计方法能有效平衡不同精度观测值之间的权重，纠正系统性偏差的影响，这一点也通过其内外符合的统计量基本一致得到验证。另外，抗差方差分量估计的内外符合精度均小于其他方案，依据参数估值的有效性可以得到，抗差方差分量估计是这四种方案中最优的。

表1

2)、静态模拟动态测试

图3为单点定位n、e方向真误差。采用的是江苏某cors站年积日为152d，时长为12h，频率为1hz的gps和bds单频观测数据进行静态模拟动态测试。除了少量粗差等不良观测值造成尖刺和跳动之外，单点定位各方向的定位误差都在5m之内。

图4为各解算方案x方向测速真误差。从图中可以看出方差分量估计受到不良观测值的影响，出现了明显偏离真值的速度估值，通过抗差估计能够较大程度地予以抑制，使得抗差方差分量的估值得到了改善。由此说明，方差分量估计对不良观测值较为敏感，通过与抗差估计结合鲁棒性能有所提升。

下表2为各方案真误差统计表。其中各解算方案统计量的表现与表1基本一致，抗差方差分量估计仍是四种参数估计方法中的最优估计。不同的是与表1相比，速度误差量有所增大，这表明测速精度与原始观测值的质量紧密相关，因此在方差分量估计中加入抗差估计是有必要的。

表2

3)、运动测试

图5为各解算方案合成的测速误差。通过在机动车上安装加拿大novatel公司的span系列高精度光纤闭环惯导组合导航系统和ubloxneo-m8t单频接收机模块，两者共用同一卫星天线，在南京市区进行动态数据采集测试。采用高精度后处理inertialexplorer软件解算组合导航的速度作为参考值与各解算方案计算ublox采集的单频数据的速度作差。从图5的测速误差图中可以看出，除了553～581历元的隧道和起始停车期间，卫星信号受到遮挡干扰外，其他路段所得的速度估值误差均能保证在1之内。

下表3是各解算方案误差统计表。从表3中可以看出，各统计量的表现与表1、表2基本一致，从而也进一步验证了其中的一般性规律。

实施例3

一种全球导航卫星系统，包括接收机，其特征为：所述系统包括基于抗差方差分量估计的码伪距/多普勒联合测速方法。通过接收机采集卫星数据，当速度解算时，通过接收机输出的观测数据，在计算平台(如电脑、手机等)进行解算并输出或显示速度信息。

可以理解的是，以上实施方式仅仅是为了说明本发明的原理而采用的示例性实施方式，然而本发明并不局限于此。对于本领域内的普通技术人员而言，在不脱离本发明的精神和实质的情况下，可以做出各种变型和改进，这些变型和改进也视为本发明的保护范围。