基于多模PPP-AR 技术的大气水汽探测性能分析

2024-01-30郭秋英薛学龙张海平李国伟毕京学

导航定位与授时 2023年6期

郭秋英,薛学龙,余淼,张海平,李国伟,毕京学

(1.山东建筑大学测绘地理信息学院,济南 250101;2.山东省国土测绘院,济南 250102)

0 引言

大气中水汽含量虽然很少,却是大气中最活跃的成分之一[1]。大气水汽的时空变化是灾害性气候形成的原因之一。因此,快速、精确地获取大气水汽对于气象预报具有重要的意义[2]。

全球卫星导航系统(global navigation satellite system,GNSS)具有高时间分辨率、高精度和全天候的探测大气水汽能力,弥补了传统水汽探测技术的不足,成为探测水汽的重要技术手段。地基GNSS反演大气可降水量(precipitable water vapor,PWV)主要有相对定位和精密单点定位(precise point positioning,PPP)两种模式,都能得到精度较高的对流层延迟(zenith total delay,ZTD)[3]。相对定位需要引入远距离测站进行联合解算,才能有效消除大部分GNSS观测误差,从而提高GNSS反演PWV 的精度。而PPP不需引入远距离测站进行解算,具有成本低、效率高及处理快等优势[4]。此外,整周模糊度固定模式下的PPP(PPP ambiguity resolution,PPP-AR)技术的发展使得PPP 收敛时间进一步缩短[5],从而实现PWV 估值精度的提升。

目前,世界上有四个全球卫星导航系统,包括中国的北斗(BeiDou navigation satellite system,BDS)、美国的GPS、俄罗斯的GLONASS和欧盟的伽利略(Galileo),这些系统相互兼容、相互发展。利用多模GNSS组合可以增加观测卫星数量,从而提高ZTD解算精度[6]。李星星等[7]使用单系统GNSS和多系统GNSS估计ZTD,实验表明,与单一GNSS相比,多模GNSS实时ZTD估计可以实现毫米级的精度和高可靠性。LU[8]基于多GNSS(GLONASS+GPS+伽利略+北斗),利用当前GNSS的所有可用观测结果,得到实时的ZTD/PWV 估值,可以实现毫米级的精度。李宏达等[9]利用GPS/BDS/GLONASS/Galileo组合获取PWV 估值,实验表明,多系统组合PPP估计ZTD的初始化时间缩短,且能够获取更高精度的ZTD和PWV。胡鹏[10]利用实时PPP 技术,分析了GNSS四系统反演水汽的精度差异,结果显示,不同系统得到的ZTD/PWV 结果存在明显的差异,多系统组合观测获取的水汽序列结果最为稳健,但没有对北斗三号获取的PWV 进行评估。夏鹏飞等[11]基于GNSS观测数据与ERA-5资料,构建对流层模型;新的ZTD 模型优于后处理GNSS-ZTD 估值,研究结果表明,对浮点PPP求解施加ZTD 约束后,垂直方向的收敛时间也有显著的改善。

以上研究结果表明,与单一或双系统相比,多系统组合显著提高了PWV 解算精度。由于北斗三号的B1C/B2a相较于其他频点具有更好的信号强度[12],因此本文在利用北斗探测水汽时选择B1C/B2a频点。此外PPP-AR 在传统PPP 的基础上加入了模糊度固定的方法,可以提高定位精度和缩短收敛时间。因此本文基于最终精密星历、采用不同方案下获取ZTD 并进一步转换为PWV,从PPPZTD 与IGS-ZTD 的相关性、PPP-ZTD 的收敛时间、PPP-ZTD 的精度以及PPP-PWV 的精度4个方面评价PPP-AR 探测大气水汽的性能,对于GNSS气象学应用具有重要的意义。

1 PPP反演大气水汽的原理与方法

1.1 PPP基本观测模型

GNSS PPP中伪距和载波相位观测方程可表示为

式中,Pj,Lj分别表示伪距和载波相位观测量;上标s为观测的卫星号;G,C,R,E 分别表示GPS,BDS,GLONASS和Galileo;下标j为载波的频率;ρs为卫星s至接收机的几何距离;c为真空中光速;dtr,dts分别为接收机和卫星钟差;d0为卫星轨道误差;tSTD为对流层斜路径延迟;为电离层延迟;分别为接收机端和卫星端的伪距硬件延迟;分别为接收机端和卫星端载波相位硬件延迟;为伪距观测值残差;为载波相位观测噪声和其他误差。

BDS-3中B1I/B3I伪距偏差可以通过无电离层组合消除,而BDS-3中B1C/B2a无法消除,为了进行改正,可以使用差分码偏差(differential code bias,DCB)产品进行改正,公式为[13]

在PPP解算中,由于未检校的相位延迟(uncalibrated phase delay,UPD)被模糊度吸收,导致模糊度不具备整数特性,无法直接固定为整数。因此,精确地将相位偏差进行分离并改正是实现PPP模糊度固定的关键,目前通常采用整数钟差法或UPD 法来消除相位偏差带来的影响。

1.2 模糊度固定方法

1.2.1 宽巷模糊度固定

模糊度固定一般采用宽窄巷(wide lane-narrow lane,WL-NL)的一个固定的过程[14],melbournewubbena(MW)组合观测值可以转换为

为了避免伪距噪声的影响,在计算LsMW时需要对伪距载波观测值进行平滑处理。具体如下

对平滑后的宽巷模糊度进行星间单差计算,首先,根据宽巷模糊度固定残差和高度角的关系,选择高度角最大的卫星作为基准星,通过计算得到星间单差宽巷模度。然后,对星间单差模糊度取整,在进行取整前后,宽巷模糊度的绝对值差值不大于0.25,则进一步进行下列验证

其中

式中,P表示模糊度固定成功率;为星间单差模糊度的实数解;为单差模糊度整数解;σ为求得的模糊度中误差。

1.2.2 窄巷模糊度固定

经过宽巷模糊度固定后,开始进行窄巷模糊度固定。根据估计的模糊度实数解和星间单差宽巷模糊度可以得出星间窄巷模糊度实数解,公式如下

其对应的协方差阵为

在LAMBDA 算法中,精密单点定位的窄巷模糊度搜索时,其输入量分别为浮点解参数个数、固定解参数个数、具有整数特性的实数解、对应的协方差矩阵以及固定解和固定解的残差平方和。整数最小二乘被认为是模糊度取整方法中最严谨的方法,获取的解为最优解。一般将最优解和次优解进行Ratio检验

1.3 PPP水汽反演的原理

斜路径对流层延迟(slant total delay,STD)可以由映射函数、水平梯度、干延迟(zenith hydrostatic delay,ZHD)与湿延迟(zenith wet delay,ZWD)构成

式中,Md为ZHD 映射函数;MW为ZWD 映射函数;ΔLg(e,a)为大气在水平方向梯度变化引起的误差;e为卫星的截止高度角;a为卫星的方位角。

一般采用Saastamoinen模型来计算ZHD

式中,p为测站下气压,单位为hPa;H为测站高度,单位为m;φ为以弧度为单位的纬度。

GNSS卫星信号STD 映射到观测站的天顶方向上得到ZTD,ZTD 减去ZHD 可以得到ZWD,然后利用ZWD 与PWV 的转换关系获得PWV。1994年,Bevis等[15]给出了ZWD 到PWV 的转换公式,计算方法为

式中,Π为水汽转换系数;ρ为液态水密度;k1=77.604 K/hPa,k2=64.790 K/hPa,k3=3.776×105K2/h Pa分别为大气折射率试验常数;k'2=k2-16.52;Rv=0.462 J·g-1·K-1为水汽气体常数;Rd=0.287 J·g-1·K-1为干空气气体常数;Tm为大气加权平均温度。其中Tm采用Bevis经验公式。基于上述方法就可以获取测站上方的PWV。

2 基于PPP-AR获取ZTD,PWV的性能分析

2.1 数据获取与实验方案

为评估多模PPP-AR在反演水汽方面的能力,选取全球范围16个MGEX 观测站,2022年3月11日到3月17日(年积日70—76)、5月7日到5月13日(年积日126—132)、7月8日到7月14日(年积日189—195)、10月23日到10月29日(年积日294—300)4个时间段的观测数据。采用GPS(G)、BDS-3(C)、GPS+BDS-3(GC)、GPS+GLO+GAL+BDS-3(GREC)以及多系统浮点解(float-GREC)五种方案对ZTD与PWV进行估值,其中GPS采用L1/L2、BDS-3 采用B1C/B2a、GLONASS 采用L1/L2、Galileo采用E1/E5a双频频点,前四种方案为PPP-AR,最后方案为多系统浮点解。实验选取全球范围16个MGEX 观测站分布图见图1。

图1 全球范围16个MGEX观测站分布图Fig.1 Distribution map of 16 MGEX observation stations around the world

实验采用武汉大学IGS数据中心提供的最终精密星历产品。该产品具有较高的精度和时效性,可用于水汽的探测。本文使用Net_Diff 1.16软件进行解算,该软件能够支持所有GNSS系统和所有频率信号、GNSS观测数据和产品下载、双频和三频PPP与PPP-AR 的解算等;且在数据处理过程中,该软件可以更加直观地选择不同的频点及其组合。具体解算策略如表1所示。

表1 数据解算策略Tab.1 Data solution strategy

自1998年以来,国际GNSS服务(international GNSS service,IGS)采用地面观测站数据并结合最终精密卫星轨道、卫星截止高度角为7°,基于PPP技术获取每日时间分辨率为5 min的ZTD 产品,精度约为4 mm[17]。这些对流层产品以每日文件的形式提供,覆盖了IGS网络中的350多个全球导航卫星系统(GNSS)站点。因此,一般将IGS-ZTD 作为参考值。

2.2 多模GNSS探测水汽的性能分析

为分析多模PPP-AR 在探测水汽方面的性能,对不同时间段和不同卫星组合获取ZTD 与PWV的估值进行了比较。从PPP-ZTD 与IGS-ZTD 的相关性、PPP-ZTD 收敛时间、PPP-ZTD 估值精度及PPP-PWV 估值精度4个方面评价多模PPP-AR 探测水汽的性能。

2.2.1 GNSS PPP-ZTD与IGS-ZTD的相关性分析

以IGS 发布的对流层产品(IGS-ZTD)作为参考值,评估不同方案下得到ZTD 估值与IGS-ZTD的相关性。图2展示了2022年3月11日到3月17日(年积日70—76)ALIC 测站5 种方案下获取ZTD 估值与IGS-ZTD 值的相关性。

图2 ALIC测站不同星座组合的相关性Fig.2 Correlation of different constellation combinations of ALIC stations

由图2可以看出,单系统获取ZTD与IGS-ZTD的相关系数(0.902)与单系统C双频B1C/B2a的相关系数(0.895),两者相关系数相当。双系统GC相关系数(0.936)高于单系统G 和单系统C双频B1C/B2a的相关系数。这表明使用两个卫星系统进行组合可以提高相关性。此外,多系统的相关系数(0.938)比双系统的相关系数更高,意味着使用更多卫星系统进行组合时,相关性会进一步增强。相较于float-GREC组合,GREC也展示了更好的相关性。

2.2.2 GNSS PPP-ZTD 收敛时间分析

以IGS 发布的对流层产品(IGS-ZTD)作为参考值,评估ZTD 序列的收敛时间。当GNSS PPPZTD 与IGS-ZTD 差值小于20 mm,并且连续20个历元差值保持在20 mm 内,认为收敛完成,这个过程需要的时间称为收敛时间。图3为DGAR 测站4个不同时间段、5种方案下的收敛时间。

图3 DGAR 测站不同时间段、不同组合下的收敛时间Fig.3 Convergence time of DGAR station under different time periods and combinations

图3显示了DGAR 测站4个不同时间段的收敛过程。在PPP收敛的过程中,刚开始各方案的初始ZTD 与IGS-ZTD 之间存在较大的差值,随着时间的推移,GNSS PPP-ZTD 与IGS-ZTD 之间的差值逐渐减小,并且变化趋势逐渐趋于一致。可以看出,双系统的收敛时间小于单系统的收敛时间,而多系统的收敛时间更加迅速,固定解在浮点解的基础之上实现了更加迅速的收敛。多系统PPP 具有更好的收敛性能,可显著提高ZTD 估计的准确性和稳定性。通过融合多个系统的信息,组合PPP可以更快地获取ZTD 的估值。为了更加直观体现测站的收敛时间,统计了所有测站不同时间段下的收敛时间,如图4所示。

图4 不同时间段各测站平均收敛时间Fig.4 Average convergence time of each station in different time periods

由图4可以看出,在4个不同时间段的测站中,基于最终精密星历获取的单系统G 和C双频B1C/B2a的平均收敛时间分别为32,33,34,33和31,32,34,32 min,双系统GC 的平均收敛时间分别为29,30,31和30 min,多系统固定解GREC 平均收敛时间分别为23,23,26 和24 min,多系统浮点解GREC平均收敛时间分别为26,27,28,27 min。根据研究数据表明:与G,C,GC 相比,多系统固定解GREC的收敛时间分别缩短了27%,25%和20%。此外,多系统固定解相比于多系统浮点解的收敛时间缩短11%。

这些收敛时间表明多系统固定解GREC 在不同时间段的收敛速度较快,相较于单系统和双系统组合,具有更高的效率和性能。由于多系统在单、双系统的基础之上获得了更多的卫星数量,观测几何图形也得到了改进。这有助于PPP的收敛,进而改善ZTD 的收敛时间和PWV 估值精度。而多系统浮点解也表现出一定的优势,虽不如多系统固定解GREC的收敛速度快,但相较于单系统G 和C固定解仍有改进。因此,多系统PPP-AR 技术在大气水汽探测中表现出显著的优势。

2.2.3 GNSS PPP反演的ZTD 精度分析

为分析不同时间段、不同星座组合下的ZTD 估值精度,选取了ALIC,BIK0两个测站,分别计算了不同方案下获取的ZTD 估值与对应参考值IGSZTD 之间的偏差。以平均偏差和均方根误差作为精度评价指标,ALIC,BIK0测站不同时间段、不同星座组合ZTD 精度对比见表2。

表2 ALIC、BIK0测站不同时间段、不同星座组合ZTD精度对比Tab.2 Comparison of ZTD accuracy of ALIC and BIK0 stations in different time periods and constellation combinations mm

表2 统计了ALIC 与BIK0 测站在不同时间段、不同星座组合下的ZTD 精度,GREC 组合获得的ZTD 估值性能比单、双系统有所提高,能够获得更高精度的ZTD 估值。可见,与单系统和双系统解算相比,多系统PPP解算ZTD 估值精度更高。

2.2.4 GNSS PPP反演PWV 的性能分析

以探空站PWV(radiosondes precipitable water vapor,RS-PWV)为参考值,对上述5种方案反演的PWV 精度进行分析。由于探空站每天只有两个时段(UTC时间00.00和12.00)的数据,因此在分析时只保留与探空站时段一致的PWV。另外,HOB2站与探空站(编号94975)相距25.1 km,WUH2与探空站(编号57494)相距5.2 km,且对流层延迟具有较强的时空差异性,故分析时只采用这两个站的PWV。5种方案获得HOB2 测站的PWV 与探空站94975获取的PWV 对比见图5(a);5种方案获得WUH2 测站的PWV 与探空站57494 获取的PWV 对比见图5(b)。

图5 HOB2和WUH2站与探空站(编号94975,57494)对比Fig.5 Comparison of HOB2 and WUH2 stations with sounding stations(No.94975,57494)

由图5可知,G,C,GC,GREC,float-GREC这5种方案得到PWV 与RS-PWV 之间的差值大部分在10 mm 之内,GREC组合获取的PWV 整体的变化更趋向于RS-PWV。以探空站RS-PWV 为参考值,5种方案计算的WUH2和HOB2站的PWV 与其相应的探空站PWV 的均方根误差见表3。

表3 WUH2与HOB2站PWV的精度对比Tab.3 Comparison of PWV accuracy between WUH2 and HOB2 stations

由表3可知,相较于单系统(G)固定解、单系统(C)固定解、双系统固定解及多系统浮点解,多系统固定解获取PWV 估值展示了更好的性能。WUH2和HOB2测站双系统GC解算均方根误差较单系统G,C双频B1C/B2a减少0.21 mm,0.29 mm 和0.10 mm,0.05 mm;相较于GC 组合,GREC 组合的PWV 精度进一步提升,WUH2 和HOB2 测站均方根误差减少0.02 mm 和0.10 mm;相较于float-GREC 组合,GREC组合精度也有所提升,WUH2 测站均方根误差减少0.02 mm,HOB2测站均方根误差减少0.08 mm。