大规模GNSS基准站网坐标解算方法研究
赵毅 吴晓莉 姜英明 冯彦同
连续运行基准站网及其系统(Continuously Operating Reference System,CORS)是全球导航卫星系统(Global Navigation Satellite System,GNSS)实现高精度导航与位置服务的支撑平台,也是建立坐标框架、监测地壳运动等科学和工程应用的重要基础设施,几乎可以应用于与位置和时间相关的所有领域[1-2].为了满足不同用户需求,国内外纷纷建立了不同尺度的GNSS基准站网,包括国际GNSS服务组织(International GNSS Service,IGS)建立的全球基准站和跟踪站网、美国连续运行参考站网系统(CORS)、德国卫星定位与导航服务系统(SAPOS)、日本GPS连续应变监测系统(COSMOS)、中国地壳运动观测网络、北斗地基增强系统网络和省市级的CORS网以及大坝等重要设施周边的工程基站网等[2-4].随着CORS站网的建设与发展,基于CORS网的相关科学研究主要分为两个方向:一个是基于CORS的高精度定位服务方法与应用研究[5-10],包括基于CORS站网的网络RTK(Real Time Kinematic)和PPP (Precise Point Positioning)定位方法、差分网定位性能,以及形变监测的应用等;另一个方向是基于CORS组网的坐标基准和坐标框架相关的研究[11-14].文献[11]将相邻省的区域CORS基准网(几十个站点)组网后经整体平差,建立统一的坐标基准;文献[12]提出融合全国的CORS成为国家级CORS,并据此升级、更新、维持2000国家大地坐标系统(CGCS2000);文献[13]分析了适用于大规模GNSS基准站网数据子网划分策略;文献[14]研究了基于多星多技术的毫米级坐标框架非线性建模方法.上述研究要么是基于区域CORS站网(几十个站),要么是基于IGS全球站网数据进行分析研究,对于国内上千个站网的坐标基准和坐标框架相关的研究还比较少.目前北斗地基增强系统“北斗全国一张网”已经有超过2 600 多个增强站,毫米级的基准坐标是提供实时高精度地基和星基服务的前提,也可以为北斗的坐标框架维持更新提供参考,因此需要研究基于GNSS大网的基准坐标解算方法.
本文首先介绍了GNSS坐标系统与主要的坐标框架,随后介绍了大规模GNSS基准站网的坐标解算方法,最后给出了基于实测数据的解算结果以及与IGS各分析中心的解算结果对比.
全球导航卫星系统的基本任务是为全球用户提供可用性、完好性、连续性和精度符合规定要求的导航定位服务,时间系统和坐标框架是GNSS赖以提供服务的基准,也是不同导航系统兼容和互操作必须要考虑的重要因素[15-18].GPS、GLONASS、Galileo、北斗(BDS)4大GNSS的坐标系统虽然定义各不相同,但都符合国际地球自转服务(IERS)规范,可以溯源到国际地球参考框架(International Terrestrial Reference Frame,ITRF),其中北斗坐标系(BeiDou Coordinate System,BDCS)与ITRF的差异在4 cm左右,并根据站速度至少每年对框架基准站点坐标进行更新[17].ITRF是由IERS建立的全球地面观测台站,通过GNSS、甚长基线干涉测量(Very Long Baseline Interferometry,VLBI)、卫星激光测距(Satellite Laser Ranging,SLR)、激光测月(Lunar Laser Ranging,LLR)、卫星多普勒定轨定位(Doppler Orbitography by Radiopositioning Integrated on Satellite,DORIS)等空间大地测量技术,对所有观测数据进行综合分析、处理,得到地面观测站的坐标和速度场,以及相应的地球定位定向参数(Earth Orientation Parameters,EOP)[19].ITRF已成为国际公认的应用最广泛、精度最高的地心坐标参考框架,经过多年发展与改进,现常用版本是ITRF2014,该版本框架考虑了由大气造成的非潮汐负载效应,并且新的观测站增加到975个站址、1 499 个测站,同时将GNSS观测网络与VLBI、SLR和DORIS大地测量技术相结合,对计算模型进行了优化,解算精度更高[19].ITRF的最新实现ITRF2020已于2022年上半年公布,但目前IGS公布的精密星历、钟差等产品仍基于IGb14(ITRF2014的IGS实现),为了将解算成果与IGS产品保持一致,本文后续的大规模GNSS基准站网(简称GNSS大网)基准坐标解算仍然基于ITRF2014框架.
从获取GNSS观测数据到基准坐标解算,需要经历全球框架站点选取、基准站网子网划分、基站单日数据处理和全球整网平差计算4个阶段.
2.1 全球框架站点选取
全球框架站点在基准站坐标解算过程中起到骨架和基石作用.利用公开的全球框架站点时序资料,按照国际通用标准,并基于如下原则选取ITRF站点:连续观测至少3年;远离板块边缘及变形区域;速度精度优于3 mm/年;至少3个不同分析中心解的速度残差小于3 mm/年.
本文共选取了231个ITRF站点,如图1所示.
图1 选择的全球ITRF站点分布Fig.1 Distribution of selected global ITRF sites
2.2 基准站子网划分
GNSS大网通常需要同步处理上千个站的GNSS观测数据,对计算和存储能力都提出了更高的要求,而且大多数IGS分析中心的高精度数据处理软件同时解算的测站数均小于100个,如美国麻省理工学院(MIT)和SCRIPPS海洋研究所共同开发的GAMIT软件目前支持同时解算的测站小于100个,瑞士伯尔尼大学研制的Bernese软件虽然对同时解算的测站数量没有软件限制,也受制于(云)服务器性能,通常同时解算测站数不超过几百个.因此对于GNSS大网需要进行分区处理,分区处理通常有间距分区法和区域分区法[20].间距分区法是根据各个测站之间的距离,将分布密集的测站抽稀,使距离较近的测站均匀地分布在不同的分区,即每个分区的测站分布较均匀,基线长度比较均匀;区域分区法就是按照地理区域进行子网划分.文献[14,20]表明间距分区法得到的基线相对精度明显优于区域分区法,因此本文采用间距分区法对全国2 600多个基站进行抽稀划分,共分为6个子网,每个子网的测站基本均匀分布在全国,如图2所示(不同分区用不同颜色点显示).
图2 基准站子网划分Fig.2 Division of base station subnet
2.3 基站单日数据处理
基站单日数据处理采用双差数据处理模型,即通过站间和星间两次作差消除接收机和卫星钟差、初始相位和硬件延迟等,大幅削弱对流层延迟和电离层延迟,获得高精度的定位定轨结果.本文通过Bernese软件实现基站的单日数据处理.将北斗全球增强站、IGS站、ITRF2014框架站共计3 000多个站点按照区域划分为6个子网分开解算,IGS核心站包含在每个子网中;将子网观测数据、IGS精密轨道、钟差、EOP等产品输入至Bernese软件,输出各子网的SINEX文件,作为后续整网平差的输入文件.
2.4 整网平差计算
整网平差计算采用自研软件以法方程叠加的方式实现网平差.整个解算过程共包括如下5个步骤:
Step1:对于第i个子网的SINEX文件,首先恢复法方程.
(1)
(2)
(3)
Step2：参数先验初值统一.
(4)
也就是说，先验初值的改变不影响法方程矩阵，只需对法方程向量进行改正即可，其改正量为-Nidx.
Step3:法方程叠加.
(5)
(6)
式(6)即为基于观测的综合解方法的观测方程，可以基于最小二乘准则求解.记每个解的权阵为Pi，则综合解的法方程为
(7)
其中，
(8)
(9)
需要注意的是，每个解的单位权方差可能不一致，每个解使用的数据处理软件、误差改正模型也可能是不同的.因此，某些情况下需要对输入解的权重进行适量调整.记每个解的权重调整因子为ai，则合并解的法方程变为
(10)
Step4:引入参数约束条件和基准条件，解算法方程，得到最终估值.
为了将解综合的成果约束到某一特定框架下，需要引入框架点的坐标.记框架点在同历元下的坐标为Xr，记约束平差的估值为Xc，根据Helmert转换关系：
Xr=Xc+Aθ,
(11)
式(11)中，θ=(Tx,Ty,Tz,D,Rx,Ry,Rz)T为框架转换的七参数，A为Helmert转换对应的设计矩阵.由式(11)可知，
(12)
E(Xr-Xc)=0.
(13)
式(13)对应的法方程为
(ETPθE)(Xr-Xc)=0,
(14)
这里Pθ为七参数的权阵，一般认为各参数相互独立.这里仅对七参数进行约束，因此可认为是一种最小约束解，施加最小约束后的法方程为
(N+ETPθE)x=b+(ETPθE)(Xr-Xc).
(15)
求解式 (15) 即可得到最小约束解.整网平差流程如图3所示.
图3 整网平差流程Fig.3 Whole network adjustment process
基于以上流程,若将输入的各子网SINEX文件替换为SINEX单天解文件,则可以生成多天解.
采用上述解算策略,选取GPS周2101～2120周(2020年4月12日—8月29日)共计20周的北斗全球增强站、IGS站、ITRF2014参考框架站等共计3 000 多个站的观测数据,对站坐标进行求解.同时给出IGS各分析中心的坐标解算与IGS周解对比.各机构的内符合精度使用其一周内基站坐标的重复性衡量,其表征的是解算固定基站时不同时段的坐标稳定性.
为说明本文采用解算策略解算站坐标内符合精度与IGS分析中心产品的差异,以南北方向3 mm、东西方向3 mm、高程6 mm为阈值,统计各IGS机构坐标重复性满足该阈值的百分比,如图4所示,从左到右分别给出了COD、COF、EMR、ESA、GFZ、GRG、JPL、MIG、MIT、NGS、SIO和QX(本文算法)的结果.图4中,COD、JPL和MIG是3天解产品,其他分析中心因为没有提供3天解产品,给出了单天解产品,本文算法(QX)分别给出单天解和3天解结果.
从图4中可以看出,COD的3天解重复性最好,本文算法计算的3天解重复性明显优于单天解,因为3天解实际上是由前后3天的单天解综合生成,实际上对基站坐标起了一定的平滑效果.ESA、GFZ单天解的水平重复性要优于JPL和MIG 3天解的结果,不同中心选择测站不同和处理策略的差异导致了内符合精度不同.因为高程方向的阈值�人椒较蛞�,所以高程方向的合格率更高.相比IGS各分析中心的产品,本文单天解的重复性与MIT相当,处于中等水平;3天解的重复性仅略低于COD,处于中上水平.
图4 本文算法和各IGS分析中心基站坐标重复性满足阈值百分比Fig.4 Percentages of base station coordinates repeatability meeting threshold for IGS analysis centers and the proposed method
以IGS周解结果作为真值,对本文解算周解的外符合精度进行评估,得到的外符合结果如表1所示.其中RMSE给出了以毫米为单位的整体统计结果,后两列给出了以IGS周解作为真值,分别以±1 cm 和±2 cm作为阈值,本文解算的坐标周解与IGS周解差值点位在南北方向、东西方向和高程方向的百分比.
表1 本文解算周解与IGS周解差值比较
从表1可以看出,本文所采用解算方法得到的基站坐标外符合精度在毫米级,南北和东西方向RMSE优于3 mm,高程RMSE优于6 mm.与IGS周解相比,水平方向差异较小,99.5%的水平点位差异在±1 cm 以内,100%的水平点位差异在±2 cm以内;高程方向差异较大,只有92.4%的高程点位差异在±1 cm 以内,仍有0.7%的高程点位差异超过±2 cm.
稳定的毫米级精度的基准站坐标是连续运行基准站网及其系统提供高精度地/星基增强服务的前提.本文采用自研的软件,基于双差数据处理和整网平差对国家北斗地基增强系统建设的“北斗全国一张网”2 600多个基准站坐标与IGS监测站和ITRF框架站联合解算,得到解算结果外符合精度南北和东西方向RMSE优于3 mm,高程RMSE优于6 mm,解算精度相比IGS各分析中心属于中上水平.
本文提出的解算方法可以支持3 000多个CORS站整网解算,后续可以用于国家GPS A、B级网以及中国地壳运动观测网络与监测网的基站坐标解算与框架维持.解算过程中需注意不同坐标框架的转换,如CGCS2000国家大地坐标系与ITRF框架之间的转换,也希望本文的研究方法和测站数据能为我国北斗坐标系框架的更新和维持做出一定的贡献.
猜你喜欢站网子网测站一种简单子网划分方法及教学案例＊计算机时代(2023年1期)2023-01-30GNSS钟差估计中的两种测站选取策略分析地理空间信息(2022年11期)2022-11-26考虑荷电状态的交直流微电网多模式协调控制策略电力建设(2022年10期)2022-09-28闪电定位网定位效率评估电子技术与软件工程(2022年6期)2022-07-07WiFi室内定位测站布设优化的DOP数值分析测绘地理信息(2022年2期)2022-04-02福海水文站气象要素对比分析陕西水利(2021年5期)2021-06-21子网划分问题研究及应用中国新通信(2019年21期)2019-03-30涪江桥流域雨量站网分布与面雨量误差关系研究中国农村水利水电(2017年10期)2017-03-22河南省水文站网现状分析评价治淮(2011年6期)2011-01-31">VLSM技术应用——以贺州学院行政办公楼网络为例科技传播(2010年19期)2010-04-12