一种基于非差无电离层的北斗钟差估计方法

1.本发明属于卫星定位中的卫星钟差估计技术领域，尤其涉及一种基于非差无电离层的北斗钟差估计方法。


背景技术：

2.目前国内外对于钟差预报及卫星钟性能分析的研究主要集中于gps，验证基于可靠模型得到的自估钟差与igs数据中心提供的精密钟差产品的精度比较，分析其可靠性与稳定性。通过对灰色预测模型和二次多项式模型分析比较，发现两个模型短期预报精度相当，对于长期预报，灰色模型优于二次多项式模型。二次多项式和线性模型对卫星钟差进行分析，总结出不同类型的卫星钟差与灰色模型指数的规律关系。小波神经网络钟差预报算法，结果表明该预报算法精度得到提高。改进的小波神经网络钟差预报，提出“shannon熵值-能量比”优选小波函数，提高了小波神经网络钟差预报的精度和可靠性。
3.对于实时钟差估计研究，2012年葛茂荣等人通过混合观测模型实时钟差估计，与非差观测模型实时钟差估计对比分析，结果表明精度相当，混合观测模型解决了初始偏差问题，最终计算速度能够提升25％。2015年武汉大学赵齐乐等基于非差观测值实现了北斗卫星实时钟差估计，通过与武汉大学提供的事后钟差作对比，实验解算得到的gps实时钟差与事后产品相当，bds实时钟差精度略低于gps；将得到的钟差带入ppp进行结算得到的定位结果与事后产品相当，表明实时钟差也具有一定可靠性。2020年黄观文等基于非差模型研究北斗二号和北斗三号卫星的实时钟差估计，并通过ppp进行验证定位精度，实验表明北斗全系统较北斗二号钟差精度std能够提升18％，获取更好的定位性能。
4.上述方法，其钟差是基于全球范围内的大量测站或部分区域网进行联合解算进行估计得到，其大量数据和时效性及精度等方面对于钟差估计具有一定的局限性。
5.针对上述技术问题，本发明基于事后精密钟差估计数学模型及参数解算理论，研究精密钟差估计的算法和处理流程，通过区域网解算进行钟差估计，与武汉大学导航中心(wum)提供的精密钟差产品进行对比验证其可靠性，为后续bds全系统实现高精度定位提供可靠的钟差数据支持。


技术实现要素：

6.为实现上述目的，本发明提供了一种基于非差无电离层的北斗钟差估计方法，包括如下步骤：
7.s1：获取测站观测值、轨道产品、igs sinex文件：利用卫星接收机观测值文件，选择多系统伪距观测值数据，武汉大学提供的轨道产品和sinex文件；
8.s2：对数据进行预处理：对数据进行粗差、周跳探测；
9.s3：针对具体误差进行改正：根据误差改正模型对相位缠绕、地球自转、相对论效应、潮汐误差，进行误差改正；
10.s4：进行参数估计：卫星轨道采用wum事后精密轨道产品，对测站坐标和接收机钟
差进行参数估计，电离层延迟使用无电离层组合进行消除，对流层延迟使用saastamoinen模型消除，模糊度采用lambda算法固定，bds系统偏差采用伪距估计和白噪声进行修正；
11.s5：残差检验、输出精密钟差产品：如果残差超出规定则重新回到上述s4步骤，若残差达标则输出精密钟差产品；
12.s6：钟差精度评估：选取一颗卫星作为参考星，通过二次做差法将其余所有卫星与该参考卫星的相对估计钟差值做差进而消除由基准钟带来的影响，最后统计std精度。
13.进一步的，步骤s2对数据进行预处理主要包括以下步骤：
14.步骤2-1，剔除粗差、探测周跳，利用mw组合进行周跳探测，其计算公式为：
[0015][0016]
式中：f1、f2代表各个波段的频率，λ1、λ2为各个频率的波长，为载波相位观测量，p1、p2为伪距载波观测量，n
wl
、λ
wl
分别为宽巷模糊度和宽巷波长，其组合观测模糊度为：
[0017][0018]
利用递推公式，求出每个历元i的平均宽巷模糊度及均方根：
[0019][0020][0021]
式中：表示前i个历元宽巷模糊度平均值；σ2(i)表示前i个历元方差，若所求结果满足式(5)周跳检测方程条件，则认为当前历元存在周跳：
[0022][0023]
步骤2-2，进行周跳修复，在gnss方程中设一个新的模糊度未知参数，原先的模糊度参数n(1)和新的模糊度参数n(2)之间存在一个条件，其公式为：
[0024]
n(1,i,j,k)＝n(2,i,j,k)+i(i,j,k)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0025]
式中：i为一个整周常数，i、j和k分别为接收机、卫星和观测频率，对n(1)和n(2)进行运算可以得到整周常数，如果i等于零，则没发生周跳。
[0026]
进一步的，步骤s3针对具体误差进行改正主要包括以下步骤：
[0027]
步骤3-1，根据改正模型对相位缠绕进行改正，计算公式为：
[0028]
δψ＝2nπ+δφ
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)
[0029][0030]
d＝x-(k
·
x)
·
k-k
×
y (9)
[0031]
d'＝x-(k
·
x')
·
k+k
×
y' (10)
[0032]
n＝int[(δψ
pre-δφ)/2π] (11)
[0033]
式中：k为卫星至测站接收机方向的单位矢量，δψ为卫星天线相位缠绕的改正
数，x和y为卫星天线的姿态矢量，x'和y'为接收机天线的姿态矢量，int[
·
]为取整函数符号；
[0034]
步骤3-2，根据改正模型对地球自转进行改正，地球内部瞬时自转轴的极移转动角度可达到0.8角秒，使相应的测站坐标在高程方向偏差能够达到25mm，在水平方向偏差能够达到7mm，可通过相应的改正模型进行消除，采用二阶love和shida数对经纬度和高程方向进行改正,计算公式为：
[0035][0036][0037][0038]
式中：δφ、δλ、δh分别为测站纬度、经度和高程方向改正数，φ、λ分别为相应测站坐标纬度和经度，为平均地极坐标，为极点变化量(rad/s)；
[0039]
步骤3-3，根据改正模型对相对论效应进行改正，卫星在轨道上运动，假设卫星钟在轨道所处惯性空间中处于匀速运动状态，在狭义相对论中，数值表达式为：
[0040]
δf1＝-0.835
×
10-10
f (15)
[0041]
式中：f为卫星钟的标准频率，卫星钟与接收机钟之间相对钟频偏差为δf1，在广义相对论中，数值表达式如下：
[0042]
δf2＝-5.284
×
10-10
f (16)
[0043]
式中：卫星钟和接收机钟之间的相对钟频偏差为δf2，结合狭义和广义相对论的总量影响，数值表达式如下：
[0044]
δf＝δf
1-δf2＝4.449
×
10-10
f (17)
[0045]
式中：卫星钟和接收机钟之间产生的相对钟频偏差为δf，卫星绕地球转动的轨道为椭圆轨道，并非为卫星轨道半径近似为a的圆轨道，且卫星在转动时在多方面因素的影响下，其卫星运动状态也随着时间改变，非相对论效应理论的匀速运动状态，其卫星钟与接收机钟之间相对钟频偏差并非为δf，依然存在较大的误差，其中较为主要的误差因素由轨道偏心率e引起，偏心率e数值虽然较小，但是其产生的周期项变化和受卫星钟频率误差影响的卫星信号传播时间的误差，需对该误差进行改正，公式如下：
[0046][0047][0048]
式中：为卫星的位置矢量，为卫星的速度矢量，c为光速传播速度，gm为地球引力场数，m为地球质量，r为卫星至地心距离，r为接收机至地心距离，ρ为卫星和接收机之间的空间距离；
[0049]
步骤3-4，根据iers2010协议改正对潮汐误差进行改正。
[0050]
进一步的，步骤s4进行参数估计主要包括以下步骤：
[0051]
步骤4-1，卫星轨道使用wum事后精密轨道产品，对测站坐标和接收机钟差进行参数估计，接收机钟差可通过不同卫星之间卫星观测值做差消除，每一个观测历元之间接收
机钟差相互独立，作为未知参数解算通过白噪声处理与测站xyz位置坐标、对流层延迟误差未知参数一起求解，bds系统偏差采用伪距估计和白噪声进行修正；
[0052]
步骤4-2，电离层延迟使用无电离层组合进行消除，基于appleton-hartree公式得到电离层延迟的计算公式如下：
[0053][0054]
载波相位观测值中和伪距观测值表示的卫星信号传播路经长度与真实站星距离之间的关系如下式：
[0055][0056][0057]
式中：n
φ
为电离层相折射率，n
p
为群折射率，ne为电子密度，tec总电子含量(total electeon content，tec)为卫星信号在传播过程中受到影响电粒子的密度积分，ρ为传播路径几何距离，ion为电离层延迟。从上式中可发现，在实际计算过程中，电离层延迟效应对于伪距和载波影响的总量级相同，符号相反；
[0058]
步骤4-3，对流层延迟使用saastamoinen模型消除，其改正模型如下：
[0059][0060]
f(φ，h)＝1-0.0026
×
cos(2
×
φ)-0.00028h (22)
[0061]
式中，p为相应测站地理位置的气压，t为相应测站地理位置的温度，e为相应测站地理位置的水气压，f(φ，h)为相应测站经纬度φ和相应测站高程h的数学函数；
[0062]
步骤4-4，对模糊度采用lambda算法固定，该方法主要由两部分组成，一是为降低模糊度参数之间相关性而进行的多维整数变换，二是在转换后的空间内进行模糊度搜索，然后再将结果转换回模糊度空间中，进而求得模糊度整数解，lambda算法通过对模糊度参数及其方差-协方差阵进行整数高斯变换(也称z变换)，将他们从原空间变换到新的空间中去，实现模糊度的降相关，计算过程为：
[0063][0064]
式中：n、n0、qn为原空间内的实数模糊度、整数模糊度以及方差-协方差阵，为新空间内的实数模糊度、整数模糊度以及方差-协方差阵，z矩阵中所有元素均为整数，行列式绝对值为1，其逆矩阵所有元素也皆为整数，z中包含多个整数变换矩阵，
然后在新的空间中进行模糊度搜索，把搜索到的模糊度组合通过逆变换变回原空间中，整数模糊度搜索条件(目标函数f)：
[0065][0066]
式中：为已知浮点的模糊度向量，为浮点模糊度的方差-协方差阵，其逆为浮点模糊度的权矩阵，为整数模糊度向量，代入目标函数可得：
[0067][0068]
求得各模糊度备选组合后，将备选组中的整数进行组合，并将各组组合带入目标函数f中，使目标函数f取得最小值的组合为最佳模糊度组合；
[0069]
步骤4-5，bds系统偏差采用伪距估计和白噪声进行修正。
[0070]
进一步的，步骤s5残差检验、输出精密钟差产品主要包括以下步骤：
[0071]
根据残差检验，如果残差超出规定则重新回到上述s4步骤，若残差达标则输出精密钟差产品。
[0072]
进一步的，步骤s6钟差精度评估主要包括以下步骤：
[0073]
步骤6-1，实际得到的钟差值为相对估计钟差，在进行钟差精度评估时，我们首先需要选取一颗卫星作为参考星，通过二次做差法将其余所有卫星与该参考卫星的相对估计钟差值做差进而消除由基准钟带来的影响，最后统计std精度，具体公式如下：
[0074][0075][0076][0077][0078]
上式中：i、s、n分别表示为对应历元时刻、卫星号和对应历元时刻总数，t、t0分别表示为估计钟差值和事后精密产品钟差值，δ、δ0分别表示为做差消除基准钟后的估计钟差值和事后精密产品钟差值，δ为做差消除基准钟的差值，为求平均值符号，δ
std
为最终钟差评估精度。
[0079]
有益效果：一种基于非差无电离层的北斗钟差估计方法与现有技术相比，具有如下优势：
[0080]
(1)本发明只需区域内较少参考站快速获取的卫星钟差的bds钟差估计算法和数据处理流程，提高了时效性；
[0081]
(2)本发明弥补了较大区域内生成的钟差精度问题，可以在一定区域内提供更加精确的钟差，为bds钟差获取提供了一种良好的估计方法。
附图说明
[0082]
图1为本发明一种基于非差无电离层的北斗钟差估计法的流程图。
具体实施方式
[0083]
按图1所示步骤，对本发明一种基于非差无电离层的北斗钟差估计方法进行详细说明。
[0084]
步骤1，获取测站观测值、轨道产品、igs sinex文件，包括以下具体步骤：
[0085]
利用卫星接收机观测值文件，选择多系统伪距观测值数据，武汉大学提供的轨道产品和sinex文件。
[0086]
步骤2，对数据进行粗差、周跳探测等预处理，包括以下具体步骤：
[0087]
步骤2-1，剔除粗差、探测周跳，利用mw组合进行周跳探测，其计算公式为：
[0088][0089]
式中：f1、f2代表各个波段的频率，λ1、λ2为各个频率的波长，为载波相位观测量，p1、p2为伪距载波观测量，n
wl
、λ
wl
分别为宽巷模糊度和宽巷波长，其组合观测模糊度为：
[0090][0091]
利用递推公式，求出每个历元i的平均宽巷模糊度及均方根：
[0092][0093][0094]
式中：表示前i个历元宽巷模糊度平均值；σ2(i)表示前i个历元方差，若所求结果满足式(5)周跳检测方程条件，则认为当前历元存在周跳：
[0095][0096]
步骤2-2，进行周跳修复，在gnss方程中设一个新的模糊度未知参数，原先的模糊度参数n(1)和新的模糊度参数n(2)之间存在一个条件，其公式为：
[0097]
n(1,i,j,k)＝n(2,i,j,k)+i(i,j,k)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0098]
式中：i为一个整周常数，i、j和k分别为接收机、卫星和观测频率，对n(1)和n(2)进行运算可以得到整周常数，如果i等于零，则没发生周跳。
[0099]
步骤3，针对具体误差进行改正，包括以下具体步骤：
[0100]
步骤3-1，根据改正模型对相位缠绕进行改正，计算公式为：
[0101]
δψ＝2nπ+δφ
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)
[0102][0103]
d＝x-(k
·
x)
·
k-k
×
y (9)
[0104]
d'＝x-(k
·
x')
·
k+k
×
y' (10)
[0105]
n＝int[(δψ
pre-δφ)/2π] (11)
[0106]
式中：k为卫星至测站接收机方向的单位矢量，δψ为卫星天线相位缠绕的改正
数，x和y为卫星天线的姿态矢量，x'和y'为接收机天线的姿态矢量，int[
·
]为取整函数符号；
[0107]
步骤3-2，根据改正模型对地球自转进行改正，地球内部瞬时自转轴的极移转动角度可达到0.8角秒，使相应的测站坐标在高程方向偏差能够达到25mm，在水平方向偏差能够达到7mm，可通过相应的改正模型进行消除，采用二阶love和shida数对经纬度和高程方向进行改正,计算公式为：
[0108][0109][0110][0111]
式中：δφ、δλ、δh分别为测站纬度、经度和高程方向改正数，φ、λ分别为相应测站坐标纬度和经度，为平均地极坐标，为极点变化量(rad/s)；
[0112]
步骤3-3，根据改正模型对相对论效应进行改正，卫星在轨道上运动，假设卫星钟在轨道所处惯性空间中处于匀速运动状态，在狭义相对论中，数值表达式为：
[0113]
δf1＝-0.835
×
10-10
f (15)
[0114]
式中：f为卫星钟的标准频率，卫星钟与接收机钟之间相对钟频偏差为δf1，在广义相对论中，数值表达式如下：
[0115]
δf2＝-5.284
×
10-10
f (16)
[0116]
式中：卫星钟和接收机钟之间的相对钟频偏差为δf2，结合狭义和广义相对论的总量影响，数值表达式如下：
[0117]
δf＝δf
1-δf2＝4.449
×
10-10
f (17)
[0118]
式中：卫星钟和接收机钟之间产生的相对钟频偏差为δf，卫星绕地球转动的轨道为椭圆轨道，并非为卫星轨道半径近似为a的圆轨道，且卫星在转动时在多方面因素的影响下，其卫星运动状态也随着时间改变，非相对论效应理论的匀速运动状态，其卫星钟与接收机钟之间相对钟频偏差并非为δf，依然存在较大的误差，其中较为主要的误差因素由轨道偏心率e引起，偏心率e数值虽然较小，但是其产生的周期项变化和受卫星钟频率误差影响的卫星信号传播时间的误差，需对该误差进行改正，公式如下：
[0119][0120][0121]
式中：为卫星的位置矢量，为卫星的速度矢量，c为光速传播速度，gm为地球引力场数，m为地球质量，r为卫星至地心距离，r为接收机至地心距离，ρ为卫星和接收机之间的空间距离；
[0122]
步骤3-4，根据iers2010协议改正对潮汐误差进行改正。
[0123]
步骤4，进行参数估计，包括以下具体步骤：
[0124]
步骤4-1，卫星轨道使用wum事后精密轨道产品，对测站坐标和接收机钟差进行参数估计，接收机钟差可通过不同卫星之间卫星观测值做差消除，每一个观测历元之间接收
机钟差相互独立，作为未知参数解算通过白噪声处理与测站xyz位置坐标、对流层延迟误差未知参数一起求解，bds系统偏差采用伪距估计和白噪声进行修正；
[0125]
步骤4-2，电离层延迟使用无电离层组合进行消除，基于appleton-hartree公式得到电离层延迟的计算公式如下：
[0126][0127]
载波相位观测值中和伪距观测值表示的卫星信号传播路经长度与真实站星距离之间的关系如下式：
[0128][0129][0130]
式中：n
φ
为电离层相折射率，n
p
为群折射率，ne为电子密度，tec总电子含量(total electeon content，tec)为卫星信号在传播过程中受到影响电粒子的密度积分，ρ为传播路径几何距离，ion为电离层延迟。从上式中可发现，在实际计算过程中，电离层延迟效应对于伪距和载波影响的总量级相同，符号相反；
[0131]
步骤4-3，对流层延迟使用saastamoinen模型消除，其改正模型如下：
[0132][0133]
f(φ，h)＝1-0.0026
×
cos(2
×
φ)-0.00028h (22)
[0134]
式中，p为相应测站地理位置的气压，t为相应测站地理位置的温度，e为相应测站地理位置的水气压，f(φ，h)为相应测站经纬度φ和相应测站高程h的数学函数；
[0135]
步骤4-4，对模糊度采用lambda算法固定，该方法主要由两部分组成，一是为降低模糊度参数之间相关性而进行的多维整数变换，二是在转换后的空间内进行模糊度搜索，然后再将结果转换回模糊度空间中，进而求得模糊度整数解，lambda算法通过对模糊度参数及其方差-协方差阵进行整数高斯变换(也称z变换)，将他们从原空间变换到新的空间中去，实现模糊度的降相关，计算过程为：
[0136][0137]
式中：n、n0、qn为原空间内的实数模糊度、整数模糊度以及方差-协方差阵，为新空间内的实数模糊度、整数模糊度以及方差-协方差阵，z矩阵中所有元素均为整数，行列式绝对值为1，其逆矩阵所有元素也皆为整数，z中包含多个整数变换矩阵，
然后在新的空间中进行模糊度搜索，把搜索到的模糊度组合通过逆变换变回原空间中，整数模糊度搜索条件(目标函数f)：
[0138][0139]
式中：为已知浮点的模糊度向量，为浮点模糊度的方差-协方差阵，其逆为浮点模糊度的权矩阵，为整数模糊度向量，代入目标函数可得：
[0140][0141]
求得各模糊度备选组合后，将备选组中的整数进行组合，并将各组组合带入目标函数f中，使目标函数f取得最小值的组合为最佳模糊度组合；
[0142]
步骤4-5，bds系统偏差采用伪距估计和白噪声进行修正。
[0143]
步骤5，：残差检验、输出精密钟差产品，包括以下具体步骤：
[0144]
根据残差检验，如果残差超出规定则重新回到上述s4步骤，若残差达标则输出精密钟差产品。
[0145]
步骤6，钟差精度评估，包括以下具体步骤：
[0146]
实际得到的钟差值为相对估计钟差，在进行钟差精度评估时，我们首先需要选取一颗卫星作为参考星，通过二次做差法将其余所有卫星与该参考卫星的相对估计钟差值做差进而消除由基准钟带来的影响，最后统计std精度，具体公式如下：
[0147][0148][0149][0150][0151]
上式中：i、s、n分别表示为对应历元时刻、卫星号和对应历元时刻总数，t、t0分别表示为估计钟差值和事后精密产品钟差值，δ、δ0分别表示为做差消除基准钟后的估计钟差值和事后精密产品钟差值，δ为做差消除基准钟的差值，为求平均值符号，δ
std
为最终钟差评估精度。
[0152]
本发明一种基于非差无电离层的北斗钟差估计方法，针对钟差产品具有延时，较大区域内精度不高，基于北斗观测数据，提出了非差无电离层的北斗钟差估计方法，该方法可以在一定的区域内根据较少的参考站快速获取卫星钟差，提高了时效性和精准度，为bds高精度定位提供了基础。
[0153]
以上所述仅为本发明的最佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。