GNSS电离层、对流层延迟计算程序设计

电离层延迟计算原理：Klobuchar模型、双频改正

对流层延迟计算原理：Saastamoinen模型

本文从上述模型出发，主要阐述程序设计思路、预期功能、算例及结果分析的阐述，以及笔者对编程过程中一些常见问题和注意事项的总结。

相关代码详见：

GNSS电离层延迟、对流层延迟计算matlab程序设计实验资源-CSDN文库https://download.csdn.net/download/m0_58307078/88615019?spm=1001.2014.3001.5503

1 实验任务及目的

1. 了解常见电离层、对流层延迟消除方法；
2. 掌握模型消除、双频改正消除电离层延迟的原理和计算方法，进行程序设计、分析消除效果；
3. 掌握模型消除对流层延迟的原理和计算方法，进行程序设计、分析消除效果；
4. 以双频改正值作为真值，对比分析GPS和北斗Klobuchar模型消除效果。

2 数据来源及编程测试环境

1．算例数据来源：PPP精密单点定位原理实验平台

2．编程环境：MATLAB 2016a

3．测试环境：MATLAB 2016a及PPP精密单点定位原理实验平台

3 实验内容及成果

3.1 电离层延迟计算程序设计

3.1.1 Klobuchar模型

1．程序设计

预期功能：基于几个关键参数（经纬度、方位角、仰角、天内秒），利用Klobuchar模型，分别实现GPS卫星和北斗卫星的电离层延迟计算。

设计思路：GPS卫星和北斗卫星的电离层延迟计算几乎相同，仅有两点区别：

1）GPS电离层模型基于地磁坐标系，需要将地理纬度转换为地磁纬度，北斗电离层模型基于地理坐标系，直接采用地理纬度；

2）北斗电离层模型最后还需根据频率进行修正：

基于上述情况，可通过选取恰当的计算过程设置if分支，分别处理两者的区别，其余步骤可共用，因此该函数的输入除了需要四个必要参数（经纬度、方位角、仰角、天内秒）外，还需输入卫星类型，用于if语句条件判断。

程序的具体计算流程如图3-1所示：

图3-1 Klobuchar模型计算流程示意

2．遇到的问题及解决过程

（1）角度参数的单位转换

问题描述：忽略了将以角度为单位的参数转化为以半圆为单位的单位转换问题，导致一系列计算错误，一定程度上造成debug时思绪混乱。

解决过程：第一次调试的直观现象是结果错误和AMP计算错误，查看工作区AMP变量值发现：AMP计算结果为负值。由于AMP只与电离层延迟参数和地磁纬度（此次调试的算例为GPS卫星）有关，首先查证电离层延迟参数输入正确，因此注意力回到纬度计算上，地磁纬度由地理纬度转换而来，要么地理纬度错误，要么转换出错，要么两者都有；进一步检查后发现：程序中一些角度参数的单位转换方式不统一：直接全部转化为弧度计算导致角度参数与有些公式单位不适配。

经查阅及询问发现的解决方法为：将角度参数的单位转化为半圆，设计三角函数计算时通过*pi化为弧度计算。

图3-2 debug思路1

值得肯定的一点是，编写程序时，在AMP计算部分增添了计算错误提示（如图3-3所示），有助于快速锁定debug方向。

图3-3 代码截图示意

（2）编写程序时只计算到延迟时间，没有计算延迟距离。

解决过程：对比计算结果和PPP实验平台上输入形式，发现大概相差8-9个数量级，仔细回想，发现自己算的延迟是时间，PPP平台上检查的是距离，刚好是光速能解决的数量级差异。

图3-4 debug思路2

3．算例结果

利用Klobuchar模型计算135周518448秒的电离层延迟算例结果如下表所示：

表3-1 Klobuchar模型电离层延迟算例结果

卫星号	仰角E	方位角A	修正前伪距	天内秒	K电离层延迟	修正后伪距
'gps2'	57.76	79.28	21457996.543	30	1.926	21457998.469
'gps6'	16.83	111.23	23989819.717	30	4.891	23989824.608
'gps13'	74.72	144.74	20175583.918	30	1.677	20175585.595
'gps29'	42.56	306.63	21916662.553	30	2.107	21916664.660
'bds2'	48.17	237.81	37172108.285	30	1.969	37172110.254
'bds5'	24.48	256.37	39035001.681	30	3.023	39035004.704
'bds8'	48.66	2.08	37397257.752	30	1.956	37397259.708
'bds28'	23.95	321.58	24717961.188	30	3.056	24717964.244

对比表中GPS与北斗的电离层延迟改正，并未发现明显差别。

3.1.2 双频改正

1．程序设计

预期功能：基于4个参数（f1、f2、ρ1、ρ2），将两个不同频率的观测量进行组合，消除电离层延迟，输出改正后的伪距。

设计思路：按组合原理编写函数，计算流程如下图所示

图3-5 双频改正计算流程示意

2．算例结果

利用双频组合计算135周518448秒的电离层延迟算例结果如下表所示：

表3-2 双频改正电离层延迟算例结果

卫星号	f1	f2	伪距1	伪距2	修正后伪距
'gps2'	1575.420	1227.600	21457994.053	21457989.389	21458001.262
'gps6'	1575.420	1227.600	23989825.057	23989827.844	23989820.749
'gps13'	1575.420	1227.600	20175582.378	20175580.412	20175585.417
'gps29'	1575.420	1227.600	21916662.031	21916659.243	21916666.341
'bds2'	1561.098	1207.140	37172116.292	37172110.002	37172125.646
'bds5'	1561.098	1207.140	39035009.182	39035007.145	39035012.211
'bds8'	1561.098	1207.140	37397269.073	37397263.605	37397277.205

3.1.3 GPS和北斗的Klobuchar模型改正对比

将双拼改正的结果作为真值，对比GPS与北斗卫星在135周518448秒的Klobuchar模型电离层延迟，算例如表4-1所示：

表3-3 GPS与北斗算例结果对比

卫星号	Klobuchar模型	双频改正（真值）	误差
'gps2'	21457998.469	21458001.262	-2.793
'gps6'	23989824.608	23989820.749	3.859
'gps13'	20175585.595	20175585.417	0.178
'gps29'	21916664.660	21916666.341	-1.681
'bds2'	37172110.254	37172125.646	-15.392
'bds5'	39035004.704	39035012.211	-7.507
'bds8'	37397259.708	37397277.205	-17.497

显然，GPS卫星的Klobuchar模型电离层延迟改正效果要优于北斗卫星。若仅从Klobuchar模型电离层延迟改正数对比GPS和北斗，两者无明显差别（见表3-1），因此，可以初步判断，GPS和北斗卫星的双频改正值差异较大。

3.2对流层延迟计算程序设计

1．程序设计

预期功能：基于4个参数（f1、f2、ρ1、ρ2），将两个不同频率的观测量进行组合，消除电离层延迟，输出改正后的伪距。

设计思路：按组合原理编写函数，程序计算流程及步骤如图3-6所示：

图3-6 对流层延迟计算流程示意

2．算例结果

2021.11.13日00：00：30时刻对流层算例结果如下表所示：

表3-4 2021.11.13 00：00：30对流层延迟算例结果

卫星号	位置	仰角E	修正前伪距	对流层延迟	修正后伪距
'gps2'	经度：113.91168121 纬度：22.52138978 高度：23.36	57.76	21457995.6	2.866118	21457998.467
'gps6'		16.83	23989816.23	8.372870	23989824.604
'gps13'		74.72	20175583.08	2.513060	20175585.594
'gps29'		42.56	21916661.08	3.584209	21916664.659
'bds1'		46.97	37381379.79	3.316325	37381383.108
'bds3'		65.08	36244560.71	2.673096	36244563.387
'bds5'		24.48	39034998.85	5.850302	39035004.701
'bds9'		34.46	338372765.8	4.284359	338372770.095
'bds33'		51.99	22892356.88	3.076802	22892359.960

代码详见：

GNSS电离层延迟、对流层延迟计算matlab程序设计实验资源-CSDN文库https://download.csdn.net/download/m0_58307078/88615019?spm=1001.2014.3001.5503

附录1 .m函数文件说明

文件名（函数名）	输入参数	函数说明
I_delay	Type,E,A,rou0,t_gps,pos	Klobuchar模型计算电离层延迟
ffv	f1,f2,rou1,rou2	双频改正计算电离层延迟
T_delay	E,rou0	对流层延迟计算