【保姆级教程/代码】使用PRIDE PPP-AR软件解算ZTD、反演GNSS PWV,并评估精度
目录
- 1 准备数据
-
- 1.1 IGS ZTD数据
- 1.2 ERA5气象数据、ERA5 PWV产品
- 1.3 GNSS观测数据
- 2 PRIDE PPP-AR软件精密单点定位
-
- 2.1 PRIDE PPP-AR软件下载和安装
- 2.2 PRIDE PPP-AR软件解算ZTD
- 3 PRIDE ZTD精度评价
- 4 基于PRIDE ZTD的GNSS PWV反演
1 准备数据
1.1 IGS ZTD数据
下载网址和批量提取方法可参考我的另一篇文章:
https://blog.csdn.net/weixin_43910227/article/details/146906263?spm=1001.2014.3001.5501
在本算例中,我下载了BJFS站2024年的IGS ZTD产品,用于评估PRIDE PPP-AR软件解算出ZTD结果的精度。
1.2 ERA5气象数据、ERA5 PWV产品
下载网址:https://cds.climate.copernicus.eu/
使用ERA5的气温、气压参与GNSS PWV反演,使用ERA5 PWV评估GNSS PWV的精度。再使用双线性插值算法得到站点处的以上值。
1.3 GNSS观测数据
下载网址:http://www.igs.gnsswhu.cn/index.php
注意,精密星历等数据不需要自己下载,PRIDE PPP-AR软件可以自动下载。若下载的观测数据是.d文件,需要将其转换为.o文件,不然PRIDE PPP-AR软件无法跑,方法可以见我的另一篇文章:https://blog.csdn.net/weixin_43910227/article/details/149309339?spm=1001.2014.3001.5501
2 PRIDE PPP-AR软件精密单点定位
2.1 PRIDE PPP-AR软件下载和安装
下载网址:https://github.com/PrideLab/PRIDE-PPPAR
软件介绍和教程:http://pride.whu.edu.cn
注意,该软件分Linux命令行版本和GUI两个版本,其中Linux版本是最新版v3.1,而GUI版本是2022年发布的v2.2版本,v3.1需要安装Ubuntu配置Linux虚拟环境。在本算例中,我使用的是window GUI版本的软件。
运行pride_pppar_winGUI.exe完成安装
2.2 PRIDE PPP-AR软件解算ZTD
使用该软件进行精密单点定位的操作步骤可参考该软件附带的使用说明PRIDE PPP-AR v3.1 manual-ch.pdf。
输入数据目录,配置各项参数,开始解算
绘制解算结果
完成PPP后,可以得到以下文件,提取出每个文件中包含的ztd文件
3 PRIDE ZTD精度评价
使用IGS ZTD评估PRIDE ZTD的精度
部分代码:
clearIgsZTDpath=\'F:\\pride精度评估\\ZTD\\BJFS\\\'; %修改为你自己IGS文件的路径PrideZTDpath=\'F:\\pride精度评估\\bjfs\\2024\\\';%读取基于ZTDIGSZTD=ReadIgsZTD(IgsZTDpath);PrideZTD=ReadPrideZTD(PrideZTDpath);IGSZTD(:,2)=IGSZTD(:,2);PrideZTD(:,2)=PrideZTD(:,2);[R,RMSE,y,a,b,Bias,STD,mae]=accuracy_comparison(IGSZTD(ia,2),PrideZTD(ib,2));figurescatter(IGSZTD(ia,2),PrideZTD(ib,2),70,\'.\');box on;
RMSE为5.6mm,说明该结果精度还不错
4 基于PRIDE ZTD的GNSS PWV反演
基于PRIDE ZTD、ERA5气象数据,完成GNSS PWV反演,原理可以参考我的另一篇文章:https://blog.csdn.net/weixin_43910227/article/details/146504760?spm=1001.2014.3001.5501
部分代码:
clear%% 要修改的部分IgsZTDpath=\'F:\\pride精度评估\\ZTD\\BJFS\\\'; %修改为你自己IGS文件的路径ERA5path2=\'F:\\pride精度评估\\era5\\2024-1y_t2m_p.nc\';ERA5path3=\'F:\\pride精度评估\\era5\\2024-1y-h.nc\';ResultPath=\'BJFS反演结果.xlsx\'; %结果保存的路径ResultPathmat=\'BJFS_PWV.mat\';%BJFSX=-2148744.689; %修改为你自己测站的坐标Y=4426641.149; %修改为你自己测站的坐标Z=4044655.758; %修改为你自己测站的坐标%% 以下不用修改%读取基于IGS ZTDdisp(\'正在读取GNSS ZTD......\')% IGSZTD=ReadIgsZTD(IgsZTDpath);PrideZTDpath=\'F:\\pride精度评估\\bjfs\\2024\\\';IGSZTD=ReadPrideZTD(PrideZTDpath);IGSZTD(:,2)=IGSZTD(:,2)/1000;%测站的空间直角坐标转大地坐标[GNSSLon,GNSSLat,GNSSH]=xyz2blh(X,Y,Z,\'WGS84\');%双线性插值出GNSS测站处的ERA5 PWV,气温、气压、并做垂直转换disp(\'正在提取测站处的ERA5......\')ERA5data=readERA52(ERA5path1,ERA5path2,ERA5path3,GNSSLat,GNSSLon,GNSSH);disp(\'开始反演GNSS PWV......\')%% 画时间序列figurescatter(IGSZTD(ia,1),GNSSPWV,5,\'filled\');hold on;dateaxis(\'x\',6);scatter(ERA5data(ib,1),ERA5data(ib,2),5,\'filled\');grid on;xlabel(\'日期\');ylabel(\'PWV (mm)\');legend(\'GNSS-PWV\',\'ERA5-PWV\',\'Location\',\'northwest\');hold offset(gca,\'LineWidth\',1,\'FontSize\',12);xlim([IGSZTD(ia(1),1),IGSZTD(ia(end),1)]);set(gcf,\'color\',\'white\',\'Position\',[787.4000000000001,645,749.5999999999999,256.8]);box on;set(gca,\'LineWidth\',1,\'FontSize\',14);%% 精度评价[R,RMSE,y,a,b,Bias,STD,mae]=accuracy_comparison(GNSSPWV,ERA5data(ib,2));figure%画散点图scatter(GNSSPWV,ERA5data(ib,2),70,\'.\');box on;xlabel(\'GNSS PWV (mm)\'); ylabel(\'ERA5 PWV (mm)\'); hold on;
完整代码、技术指导,可咸鱼搜:GNSS_MET