城市复杂场景下GNSS定位的因子图优化方法及其抗差性能分析 Factor Graph Optimization for Urban Environment GNSS Positioning and Robust Performance Analysis 1 2 3 first-author 张小红，博士，教授，主要从事导航定位技术及其应用研究。xhzhang@sgg.whu.edu.cn 张小红，博士，教授，主要从事导航定位技术及其应用研究。xhzhang@sgg.whu.edu.cn xhzhang@sgg.whu.edu.cn 1 1 2 corresp 朱锋，博士，特聘副研究员。fzhu@whu.edu.cn 朱锋，博士，特聘副研究员。fzhu@whu.edu.cn fzhu@whu.edu.cn 武汉大学测绘学院，湖北 武汉，430079 武汉大学测绘学院，湖北 武汉，430079 School of Geodesy and Geomatics, Wuhan 430079, China School of Geodesy and Geomatics, Wuhan 430079, China 湖北珞珈实验室，湖北 武汉，430079 湖北珞珈实验室，湖北 武汉，430079 Hubei Luojia Laboratory, Wuhan University, Wuhan 430079, China Hubei Luojia Laboratory, Wuhan University, Wuhan 430079, China 武汉大学中国南极测绘研究中心，湖北 武汉，430079 武汉大学中国南极测绘研究中心，湖北 武汉，430079 Chinese Antarctic Center of Surveying and Mapping, Wuhan University, Wuhan 430079, China Chinese Antarctic Center of Surveying and Mapping, Wuhan University, Wuhan 430079, China 北斗/全球导航卫星系统（global navigation satellite system,GNSS）在开阔环境下可以提供连续可靠的高精度导航定位服务，但是在城市复杂场景下，GNSS多路径与非视距信号严重、粗差与周跳发生频繁，导航定位能力仍然存在不足。相较于扩展卡尔曼滤波(extended Kalman filter, EKF)方法，因子图优化能够充分利用历史观测，通过窗口内历元间约束与冗余观测信息共同抑制异常数据影响。构建了基于滑动窗口因子图优化的GNSS定位模型，通过验后残差迭代分析进行粗差探测，并从最小可探测误差、粗差探测成功率、定位精度提升等方面深入分析因子图优化与EKF的抗差性能。以城市复杂场景数据进行处理验证，结果表明，因子图优化的最小可探测误差减小了11.92%~32.56%，粗差探测成功率提升了3.84%~10.47%，GNSS定位精度提升了11.29%~25.99%。总体而言，对于城市复杂场景下的GNSS导航定位应用，因子图优化具备更好的抗差性能和定位精度，有望取代现有基于单历元观测值的EKF模型。 Objectives 2 BeiDou satellite navigation system （BeiDou）/global navigation satellite system（GNSS） can provide continuous and reliable high-precision navigation and positioning service in open-sky, yet in urban environment, suffering from outliers and cycle slips caused by severe multi-path and non-line-of-sight signals, the navigation and positioning capability remains inadequate. Compared with the extended Kalman filter (EKF), factor graph optimization (FGO) can make full use of historical measurements and restrain the influence of anomalous data through constraints and redundant measurements inside the window. Methods 2 A GNSS positioning model based on sliding-window FGO is constructed, adopting posteriori residuals test as outlier detection method, and analysis of robust performance between EKF and FGO is presented in terms of minimum detectable bias, correct detection rate and positioning accuracy. Results 2 The result of urban environment experiment shows that the minimum detectable bias is reduced by 11.92%-32.56%, the correct detection rate is improved by 3.84%-10.47%, and the GNSS positioning accuracy is improved by 11.29%-25.99%. Conclusions 2 Overall，for Application of GNSS navigation and positioning in urban environment，FGO has better robustness and positioning accuracy，could replace EKF models based on single epoch observations. 城市复杂场景 因子图优化 粗差探测 定位精度提升 urban environment factor graph optimization outlier detection positioning accuracy promotion 全球卫星导航系统（global navigation satellite system， GNSS）能够在全球范围内提供全时段连续可用的高精度导航定位信息 ［ 1 1 ］ ，在开阔场景下已广泛应用于海陆空等各种载体的时空信息确定。然而，城市复杂场景下的GNSS信号容易受到建筑物遮挡或反射，多路径与非视距信号严重、粗差与周跳发生频繁，GNSS定位精度与可靠性无法得到保障 ［ 2 2 - 3 3 2-3 ］ 。 一直以来，扩展卡尔曼滤波（extended Kalman filter， EKF）是卫星导航领域的常用算法 ［ 4 4 ］ ，EKF通过状态方程将历史观测转化为当前状态的先验信息进行传递，传递后的历史观测所对应的状态线性化点就会被固定，当历史观测存在未探测的粗差时，线性化点误差较大，容易导致先验信息受到污染，从而影响定位精度。相比而言，因子图优化（factor graph optimization， FGO）通过重复迭代更新线性化点，能够充分利用历史观测，挖掘观测信息在时间维度上的约束信息，抑制粗差的影响 ［ 5 5 ］ 。 Kschischang等 ［ 6 6 ］ 在纠错码与和积算法的基础上提出了因子图的概念，为后续研究提供了理论基础，Dellaert等 ［ 7 7 ］ 将其作为优化方法引入导航问题中，因子图优化是目前同步定位与建图领域中后端优化的常用估计器之一 ［ 8 8 ］ 。2021—2022年连续两届谷歌智能手机分米级定位比赛第一名的定位算法均采用因子图优化模型 ［ 9 9 - 10 10 9-10 ］ ，证明了因子图优化在GNSS定位问题中的可行性。Watson等 ［ 11 11 ］ 实现了基于因子图优化的PPP（precise point positioning）定位模型，收敛阶段的定位精度相较于EKF方法有大幅度提升。文伟松等 ［ 12 12 ］ 基于双差观测值构建了实时差分定位（real - time kinematic，RTK）因子图模型，通过城市峡谷实验对定位结果的抗差效果进行了分析，证明了因子图优化对粗差具有较好的抑制能力 ［ 13 13 ］ 。徐浪等 ［ 14 14 ］ 针对标准单点定位（standard point positioning， SPP）伪距解与RTK浮点解，比较了单历元最小二乘与因子图优化的定位性能，然而两种优化方法的观测信息实质上是不对等的。Bai等 ［ 15 15 ］ 和Jiang等 ［ 16 16 ］ 通过历元间差分观测值以增强因子图的历元间约束，并在此基础上研究了深度学习、分组解算法在因子图抗差上的应用。然而，上述研究多在定位结果层面对因子图优化方法的抗差效果进行评价，缺少对EKF与FGO方法在抗差机理上的深入分析。因此，本文基于滑动窗口因子图优化的GNSS定位模型，从最小可探测误差、粗差探测成功率与定位精度等方面对EKF与FGO方法的抗差机理和性能进行全面的研究与分析。 1　滑动窗口的因子图优化 1 在导航定位问题中，因子图以图的形式对最优估计问题建模，根据最大验后估计（maximum a posteriori estimate， MAP）准则对待估状态进行求解。在优化过程中，图规模必然随时间不断增长，导致图优化的实时性下降，因此引入滑动窗口平衡精度与效率，提高系统实时性。 1.1　因子图优化 2 因子图是概率图模型的一种，可以简明直观地表达一个优化问题。因子图在具体形式上是有两种节点类型的无向二分图，以因子节点编码局部函数 <math id="M1"><mi>ϕ</mi><mfenced separators="|"><mrow><mi>x</mi></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079729&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079721&type= 8.21266651 3.97933316 ，变量节点编码待估变量 x ，而边则表示了函数与变量间的相关关系，因子图表示了全局函数与局部函数间的因式分解关系。例如， 图1 图1 所对应的全局函数 <math id="M2"><mi>ϕ</mi><mfenced separators="|"><mrow><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">3</mn></mrow></msub></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079732&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079722&type= 19.55800056 3.97933316 即为： <math display="block" id="M3"><mtable><mtr><mtd><mi>ϕ</mi><mfenced separators="|"><mrow><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">3</mn></mrow></msub></mrow></mfenced><mo>=</mo><msub><mrow><mi>ϕ</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></mrow></mfenced><msub><mrow><mi>ϕ</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></mrow></mfenced><msub><mrow><mi>ϕ</mi></mrow><mrow><mn mathvariant="normal">3</mn></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></mrow></mfenced><mo>×</mo></mtd></mtr><mtr><mtd><msub><mrow><mi>ϕ</mi></mrow><mrow><mn mathvariant="normal">4</mn></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">3</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></mrow></mfenced><msub><mrow><mi>ϕ</mi></mrow><mrow><mn mathvariant="normal">5</mn></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">3</mn></mrow></msub></mrow></mfenced></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079735&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079741&type= 63.66933060 9.82133293 而在优化问题中，因子图中的因子节点对应了概率密度函数，则全局函数对应联合概率密度函数。对于任一因子图表示的优化问题，都可以根据最大验后估计准则进行求解。最大验后估计可以归结为最大化所有因子节点的乘积： <math display="block" id="M4"><msup><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi mathvariant="normal">M</mi><mi mathvariant="normal">A</mi><mi mathvariant="normal">P</mi></mrow></msup><mo>=</mo><munder><mrow><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow><mrow><mi mathvariant="bold-italic">X</mi></mrow></munder><mi>ϕ</mi><mfenced separators="|"><mrow><mi mathvariant="bold-italic">X</mi></mrow></mfenced><mo>=</mo><munder><mrow><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">x</mi></mrow><mrow><mi mathvariant="bold-italic">X</mi></mrow></munder><munder><mo movablelimits="false">∏</mo><mrow><mi>i</mi></mrow></munder><mrow/><msub><mrow><mi>ϕ</mi></mrow><mrow><mi>i</mi></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079749&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079756&type= 68.15666962 6.60400009 在实际问题中，变量的真实概率密度函数复杂且难以获得，通常假设为高斯分布，因此所有因子都具备以下形式： <math display="block" id="M5"><msub><mrow><mi>ϕ</mi></mrow><mrow><mi>i</mi></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced><mo>∝</mo><mi mathvariant="normal">e</mi><mi mathvariant="normal">x</mi><mi mathvariant="normal">p</mi><mfenced open="{" close="}" separators="|"><mrow><mo>-</mo><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mn mathvariant="normal">2</mn></mrow></mfrac><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">z</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>-</mo><msub><mrow><mi>h</mi></mrow><mrow><mi>i</mi></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi mathvariant="bold-italic">i</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079768&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079751&type= 53.84799957 9.05933285 式中， <math id="M6"><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><mo>·</mo></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi mathvariant="bold-italic">i</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079761&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079770&type= 6.68866634 5.67266655 表示马氏范数。对 式（3） 式（3） 取负对数并乘以2，就可以转化为一个最小二乘估计问题： <math display="block" id="M7"><msup><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi mathvariant="normal">M</mi><mi mathvariant="normal">A</mi><mi mathvariant="normal">P</mi></mrow></msup><mo>=</mo><munder><mrow><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">n</mi></mrow><mrow><mi mathvariant="bold-italic">X</mi></mrow></munder><munder><mo movablelimits="false">∑</mo><mrow><mi>i</mi></mrow></munder><mrow/><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">z</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>-</mo><msub><mrow><mi>h</mi></mrow><mrow><mi>i</mi></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi mathvariant="bold-italic">i</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079775&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079763&type= 54.10200119 7.53533268 由上述推导可得，在高斯噪声分布的假设下，因子图实质上是以图形式表达的最小二乘问题，因子即为优化问题中的代价函数。 1.2　滑动窗口优化 2 在导航定位问题中，观测信息以时间序列先后到达，图规模也随之增长。为了降低图优化的计算成本以兼顾实时性，通常引入滑动窗口来限制图规模。当图中变量节点数量超过窗口长度时，则在新节点加入后，通过边缘化删去窗口中最早的变量节点，并向图中加入边缘化后获得的先验因子以保留历史观测信息 ［ 17 17 ］ 。 如 图2 图2 所示，假设窗口中待边缘化变量为 <math id="M8"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>m</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079872&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079878&type= 3.89466691 3.72533321 ，与其相连接的因子集合为 Z m ，与 <math id="M9"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>m</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079872&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079878&type= 3.89466691 3.72533321 具有相关性的变量为 <math id="M10"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079863&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079853&type= 3.21733332 3.72533321 ，则边缘化仅在 Z m 形成的子图上进行，即： <math display="block" id="M11"><msup><mrow><mfenced open="{" close="}" separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>m</mi></mrow></msub><mo>,</mo><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub></mrow></mfenced></mrow><mrow><mi mathvariant="normal">M</mi><mi mathvariant="normal">A</mi><mi mathvariant="normal">P</mi></mrow></msup><mo>=</mo><munder><mrow><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">n</mi></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>m</mi></mrow></msub><mo>,</mo><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub></mrow></munder><munder><mo movablelimits="false">∑</mo><mrow><msubsup><mrow><mi mathvariant="bold-italic">z</mi></mrow><mrow><mi>m</mi></mrow><mrow><mi>i</mi></mrow></msubsup><mo>∈</mo><msub><mrow><mi mathvariant="bold-italic">Z</mi></mrow><mrow><mi>m</mi></mrow></msub></mrow></munder><mrow/><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msubsup><mrow><mi mathvariant="bold-italic">z</mi></mrow><mrow><mi>m</mi></mrow><mrow><mi>i</mi></mrow></msubsup><mo>-</mo><msub><mrow><mi>h</mi></mrow><mrow><mi>i</mi></mrow></msub><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi mathvariant="bold-italic">i</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079808&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079806&type= 65.44733429 8.12800026 将 式（5） 式（5） 对应的法方程记为： <math display="block" id="M12"><mfenced open="[" close="]" separators="|"><mrow><mtable><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>m</mi><mi>m</mi></mrow></msub></mtd><mtd><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>m</mi><mi>r</mi></mrow></msub></mtd></mtr><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>r</mi><mi>m</mi></mrow></msub></mtd><mtd><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>r</mi><mi>r</mi></mrow></msub></mtd></mtr></mtable></mrow></mfenced><mfenced open="[" close="]" separators="|"><mrow><mtable><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>m</mi></mrow></msub></mtd></mtr><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub></mtd></mtr></mtable></mrow></mfenced><mo>=</mo><mfenced open="[" close="]" separators="|"><mrow><mtable><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>m</mi></mrow></msub></mtd></mtr><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>r</mi></mrow></msub></mtd></mtr></mtable></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079818&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079816&type= 36.66066742 10.24466610 式中， <math id="M13"><mi mathvariant="bold-italic">Λ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079828&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079819&type= 2.96333337 2.79399991 为参数对应的信息矩阵，其逆矩阵即为协方差阵 <math id="M14"><mi mathvariant="bold-italic">D</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079823&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079829&type= 3.04800010 2.79399991 。利用舒尔补消元可以得到 <math id="M15"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079863&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079853&type= 3.21733332 3.72533321 的分布为 <math id="M16"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub><mo>~</mo><mi>N</mi><mfenced separators="|"><mrow><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>t</mi></mrow></msub><mo>,</mo><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079845&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079827&type= 30.47999954 4.74133301 ，其中， <math display="block" id="M17"><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>t</mi></mrow></msub><mo>=</mo><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>r</mi></mrow></msub><mo>-</mo><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>r</mi><mi>m</mi></mrow></msub><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>m</mi><mi>m</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>m</mi></mrow></msub></mtd></mtr><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow></msub><mo>=</mo><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>r</mi><mi>r</mi></mrow></msub><mo>-</mo><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>r</mi><mi>m</mi></mrow></msub><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>m</mi><mi>m</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>m</mi><mi>r</mi></mrow></msub></mtd></mtr></mtable></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079849&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079847&type= 36.91466904 9.90600014 注意该分布对应边缘化时的线性化点为 <math id="M18"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><msub><mrow><mi>r</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079859&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079842&type= 3.89466691 3.64066648 ，在下一次优化时， <math id="M19"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079863&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079853&type= 3.21733332 3.72533321 的值已经发生改变，需要进行更新，更新后的先验分布为： <math display="block" id="M20"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub><mo>~</mo><mi>N</mi><mfenced separators="|"><mrow><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>t</mi></mrow></msub><mo>+</mo><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub><mo>-</mo><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><msub><mrow><mi>r</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub></mrow></msub></mrow></mfenced><mo>,</mo><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079867&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079855&type= 50.20733261 5.41866684 将其作为先验约束加入图中，代价函数为： <math display="block" id="M21"><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msup><mrow><mi mathvariant="bold-italic">e</mi></mrow><mrow><mi>M</mi></mrow></msup></mrow></mfenced></mrow><mrow><msup><mrow><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow></msub></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msup></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo>=</mo><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msup><mrow><msub><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow></msub></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msup><msub><mrow><mi mathvariant="bold-italic">g</mi></mrow><mrow><mi>t</mi></mrow></msub><mo>+</mo><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>r</mi></mrow></msub><mo>-</mo><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><msub><mrow><mi>r</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow></msub></mrow></msub></mrow></mfenced><mo>-</mo><msup><mrow><mi>h</mi></mrow><mrow><mi>M</mi></mrow></msup><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>k</mi></mrow></msub></mrow></mfenced></mrow></mfenced></mrow><mrow><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079876&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079874&type= 69.08799744 7.02733374 被边缘化的状态 <math id="M22"><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>m</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079872&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079878&type= 3.89466691 3.72533321 线性化点在其对应的观测因子 Z m 中固定，会造成一定的精度损失，因此滑动窗口是平衡效率与精度的策略。 2　GNSS定位的因子图优化模型 1 基于滑动窗口的因子图优化，本文实现的GNSS伪距定位算法如 图3 图3 所示，融合了边缘化先验信息、多频多系统的伪距观测值与GNSS测速结果，并对窗口内观测值进行验后粗差探测以提升定位结果的可靠性。 2.1　GNSS定位因子 2 伪距观测值定位以距离交会的方式测定接收机的绝对位置，其观测方程为： <math display="block" id="M23"><msubsup><mrow><mi>P</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>i</mi><mo>,</mo><mi>f</mi></mrow></msubsup><mo>=</mo><msubsup><mrow><mi>ρ</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>i</mi></mrow></msubsup><mo>+</mo><mi>c</mi><mfenced separators="|"><mrow><mi mathvariant="normal">δ</mi><msub><mrow><mi>t</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>-</mo><mi mathvariant="normal">δ</mi><msup><mrow><mi>t</mi></mrow><mrow><mi>i</mi></mrow></msup></mrow></mfenced><mo>+</mo><msubsup><mrow><mi>I</mi></mrow><mrow><mi>k</mi></mrow><mrow/></msubsup><mo>+</mo><msubsup><mrow><mi>T</mi></mrow><mrow><mi>k</mi></mrow><mrow/></msubsup><mo>+</mo><mi>ε</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079886&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079885&type= 57.57332993 4.48733330 式中， <math id="M24"><msubsup><mrow><mi>P</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>i</mi><mo>,</mo><mi>f</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079904&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079903&type= 5.07999992 3.97933316 表示第 k 历元上第 i 颗卫星在 f 频点上的伪距观测值； <math id="M25"><msubsup><mrow><mi>ρ</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>i</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079909&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079905&type= 3.13266683 4.74133301 表示卫星到接收机的几何距离； <math id="M26"><mi mathvariant="normal">δ</mi><msub><mrow><mi>t</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079921&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079911&type= 3.80999994 3.72533321 表示接收机钟差； <math id="M27"><mi mathvariant="normal">δ</mi><msup><mrow><mi>t</mi></mrow><mrow><mi>i</mi></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079924&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079916&type= 3.64066648 3.04800010 、 <math id="M28"><msubsup><mrow><mi>I</mi></mrow><mrow><mi>k</mi></mrow><mrow/></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079936&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079925&type= 2.70933342 3.97933316 、 <math id="M29"><msubsup><mrow><mi>T</mi></mrow><mrow><mi>k</mi></mrow><mrow/></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079940&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079937&type= 3.89466691 3.97933316 分别表示卫星钟差、电离层延迟与对流层延迟； <math id="M30"><mi>ε</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079944&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079942&type= 1.43933344 2.79399991 表示未建模的其他误差，包括大气残余误差、多路径误差、接收机噪声等。随机模型则采用高度角定权模型： <math display="block" id="M31"><msup><mrow><msub><mrow><mi>σ</mi></mrow><mrow><mi>P</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>=</mo><msup><mrow><mi>a</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>+</mo><msup><mrow><mi>b</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>/</mo><mi mathvariant="normal">s</mi><mi mathvariant="normal">i</mi><msup><mrow><mi mathvariant="normal">n</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mi>E</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079950&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079946&type= 31.32666588 4.40266657 式中， a 和 b 为伪距噪声参数； E 为观测值对应的卫星高度角。伪距因子为二元因子，与对应历元的位置变量和接收机钟差变量连接，其代价函数为： <math display="block" id="M32"><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msubsup><mrow><mi mathvariant="bold-italic">e</mi></mrow><mrow><mi>k</mi></mrow><mrow><msub><mrow><mi>P</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></msubsup></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>P</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo>=</mo><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msubsup><mrow><mi>P</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>i</mi><mo>,</mo><mi>f</mi></mrow></msubsup><mo>-</mo><msup><mrow><mi>h</mi></mrow><mrow><mi>P</mi></mrow></msup><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>,</mo><mi mathvariant="normal">δ</mi><msub><mrow><mi>t</mi></mrow><mrow><mi>k</mi></mrow></msub></mrow></mfenced></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>P</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079968&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079951&type= 46.22800064 6.34999990 为了在因子图中建立GNSS观测值在时间上的相关性，向图中引入GNSS测速结果对相邻历元的位置变量进行约束，观测方程为： <math display="block" id="M33"><msub><mrow><mi mathvariant="bold-italic">v</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>=</mo><mrow><mrow><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>-</mo><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>k</mi><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msub></mrow></mfenced></mrow><mo>/</mo><mrow><mi mathvariant="normal">Δ</mi><mi>t</mi></mrow></mrow></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079972&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079970&type= 30.90333366 4.48733330 式中， <math id="M34"><msub><mrow><mi mathvariant="bold-italic">v</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079976&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079974&type= 3.04800010 3.72533321 为GNSS测速结果； <math id="M35"><mi mathvariant="normal">Δ</mi><mi>t</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079980&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079978&type= 3.72533321 2.79399991 为相邻历元时间差。速度约束因子的代价函数为： <math display="block" id="M36"><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msubsup><mrow><mi mathvariant="bold-italic">e</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>V</mi></mrow></msubsup></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>V</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo>=</mo><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">v</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>-</mo><msup><mrow><mi>h</mi></mrow><mrow><mi>V</mi></mrow></msup><mfenced separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>,</mo><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>k</mi><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msub></mrow></mfenced></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>V</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079954&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079952&type= 47.49800110 6.43466663 2.2　GNSS因子图模型 2 基于上述因子编码形式，实现了如 图4 图4 所示的GNSS因子图定位模型。模型的状态空间 <math id="M37"><msup><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079987&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079956&type= 4.23333359 3.04800010 包括接收机位置与钟差： <math display="block" id="M38"><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>=</mo><msup><mrow><mfenced open="[" close="]" separators="|"><mrow><mtable><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">x</mi></mrow><mrow><mi>k</mi></mrow></msub></mtd><mtd><mi mathvariant="normal">δ</mi><msub><mrow><mi mathvariant="bold-italic">t</mi></mrow><mrow><mi>k</mi></mrow></msub></mtd></mtr></mtable></mrow></mfenced></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup></mtd></mtr><mtr><mtd><msup><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup><mo>=</mo><mfenced open="{" close="}" separators="|"><mrow><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo>,</mo><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo>,</mo><mo>⋯</mo><mo>,</mo><msub><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi>n</mi></mrow></msub></mrow></mfenced></mtd></mtr></mtable></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079991&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079989&type= 35.56000137 11.26066685 式中， <math id="M39"><mi mathvariant="normal">δ</mi><msub><mrow><mi mathvariant="bold-italic">t</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079961&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079992&type= 3.89466691 3.72533321 表示窗口内第 k 历元的接收机钟差向量，包括该历元接收到所有卫星系统的钟差； n 表示窗口长度，设置为10。需要说明的是，为了因子图的简洁直观，一个历元上仅展示了两个伪距因子与两种系统的接收机钟差，实际情况下与接收到的数量一致。 结合 式（9） 式（9） 、 式（12） 式（12） 、 式（14） 式（14） 三种因子，GNSS因子图模型的优化目标函数为： <math display="block" id="M40"><msup><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup><mo>=</mo><munder><mrow><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">n</mi></mrow><mrow><msup><mrow><mi mathvariant="bold-italic">X</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup></mrow></munder><munder><mo movablelimits="false">∑</mo><mrow><mi>i</mi><mo>,</mo><mi>k</mi></mrow></munder><mrow/><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msubsup><mrow><mi mathvariant="bold-italic">e</mi></mrow><mrow><mi>k</mi></mrow><mrow><msub><mrow><mi>P</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></msubsup></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>P</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo>+</mo><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msubsup><mrow><mi mathvariant="bold-italic">e</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>V</mi></mrow></msubsup></mrow></mfenced></mrow><mrow><msub><mrow><mi mathvariant="bold-italic">Σ</mi></mrow><mrow><mi>V</mi></mrow></msub></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup><mo>+</mo><msubsup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><msup><mrow><mi mathvariant="bold-italic">e</mi></mrow><mrow><mi>M</mi></mrow></msup></mrow></mfenced></mrow><mrow><msubsup><mrow><mi mathvariant="bold-italic">Λ</mi></mrow><mrow><mi>t</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080001&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45079999&type= 67.22533417 7.61999989 2.3　滑动窗口粗差探测 2 在城市复杂场景中，多路径效应与非视距信号发生频繁，严重影响定位精度与可靠性。因此，在完成 式（16） 式（16） 的优化后，对窗口内的所有观测进行粗差探测与剔除。优化后获得窗口内GNSS观测值的验后残差为： <math display="block" id="M41"><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><msup><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup><mo>=</mo><msup><mrow><mfenced open="[" close="]" separators="|"><mrow><mtable><mtr><mtd><mtable><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub></mtd><mtd><msub><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub></mtd><mtd><mo>⋯</mo></mtd></mtr></mtable></mtd><mtd><msub><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mi>n</mi></mrow></msub></mtd></mtr></mtable></mrow></mfenced></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup></mtd></mtr><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mi>k</mi></mrow></msub><mo>=</mo><msup><mrow><mfenced open="[" close="]" separators="|"><mrow><mtable><mtr><mtd><mtable><mtr><mtd><msubsup><mrow><mi>r</mi></mrow><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msubsup></mtd><mtd><msubsup><mrow><mi>r</mi></mrow><mrow><mi>k</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msubsup></mtd><mtd><mo>⋯</mo></mtd></mtr></mtable></mtd><mtd><msubsup><mrow><mi>r</mi></mrow><mrow><mi>k</mi></mrow><mrow><mi>m</mi></mrow></msubsup></mtd></mtr></mtable></mrow></mfenced></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup></mtd></mtr></mtable></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080005&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080003&type= 39.11600113 11.93799973 式中， <math id="M42"><msup><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080012&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080011&type= 2.87866688 3.04800010 表示窗口内所有历元观测值的残差向量； <math id="M43"><msub><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mi>k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080014&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080013&type= 2.62466669 3.72533321 表示窗口内第 k 个历元的残差向量。对残差向量构建卡方检验量： <math display="block" id="M44"><mi>T</mi><mo>=</mo><msup><mrow><msup><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup><msubsup><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">Z</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup><msup><mrow><mi mathvariant="bold-italic">r</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080016&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080015&type= 21.84399986 3.97933316 式中， <math id="M45"><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">Z</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080020&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080018&type= 4.48733330 3.72533321 表示观测值的协方差阵，并根据显著性水平为 <math id="M46"><mi>α</mi><mo>=</mo><mn mathvariant="normal">0.1</mn><mi mathvariant="normal">%</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080038&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080022&type= 15.32466698 3.47133350 设置阈值 <math id="M47"><msub><mrow><mi>T</mi></mrow><mrow><mi>α</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080042&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080040&type= 3.80999994 3.72533321 。若检验通过，则认为窗口内观测值不存在粗差；若检验不通过，则基于均值漂移模型假设对窗口内最大的残差进行剔除并重新进行优化，通过重新线性化以充分利用窗口内观测值的相关性完成抗差。 同时，为了检验抗差性能，本文还对方法的最小可探测误差（minimal detectable bias，MDB） <math id="M48"><msub><mrow><mo>∇</mo></mrow><mrow><mi mathvariant="normal">M</mi><mi mathvariant="normal">D</mi><mi mathvariant="normal">B</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080046&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080044&type= 7.61999989 4.06400013 进行了评估，计算公式如下： <math display="block" id="M49"><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><msub><mrow><mo>∇</mo></mrow><mrow><mi mathvariant="normal">M</mi><mi mathvariant="normal">D</mi><mi mathvariant="normal">B</mi></mrow></msub><mo>=</mo><mfrac><mrow><mi>δ</mi></mrow><mrow><mroot><mrow><msubsup><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">Z</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">r</mi></mrow></msub><msubsup><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">Z</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup></mrow><mrow/></mroot></mrow></mfrac></mtd></mtr><mtr><mtd><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">r</mi></mrow></msub><mo>=</mo><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">Z</mi></mrow></msub><mo>-</mo><mi mathvariant="bold-italic">H</mi><msup><mrow><mfenced separators="|"><mrow><msup><mrow><mi mathvariant="bold-italic">H</mi></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup><msubsup><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">Z</mi></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msubsup><mi mathvariant="bold-italic">H</mi></mrow></mfenced></mrow><mrow><mo>-</mo><mn mathvariant="normal">1</mn></mrow></msup><msup><mrow><mi mathvariant="bold-italic">H</mi></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup></mtd></mtr></mtable></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080048&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080025&type= 51.22333527 17.61066628 式中， <math id="M50"><msub><mrow><mi mathvariant="bold-italic">D</mi></mrow><mrow><mi mathvariant="bold-italic">r</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080027&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080050&type= 4.06400013 3.72533321 表示残差的协方差阵； <math id="M51"><mi mathvariant="bold-italic">H</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080029&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080052&type= 3.47133350 2.79399991 表示观测雅各比矩阵； <math id="M52"><mi>δ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080032&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080054&type= 1.94733346 2.79399991 表示非中心化参数，通常取值为4.13。 3　实验测试与分析 1 为了验证本文GNSS伪距定位的因子图优化方法性能，使用两组城市复杂场景的车载动态数据，并与EKF的伪距单点定位进行对比分析。解算策略如 表 1 表 1 所示，为了保证分析结果的可靠性，两种方法在输入数据、解算策略上保持一致。 项目 策略 观测数据 多频多系统数据 星历 广播星历 截止高度角/（°） 15 截止信噪比/dBHz 20 速度约束 GNSS测速结果 随机模型 高度角模型 粗差探测策略 卡方检验+最大化残差 两组实验均在城市复杂场景下进行，实验一采用测量型接收机，实验路线与观测环境如 图5 图5 所示，包含城市公路、高楼遮挡、高架遮挡等典型复杂场景；实验二采用低成本u-blox接收机，仅能接收单频数据，实验路段包含城市峡谷等场景，两组数据的参考真值均由高精度组合导航设备经过紧耦合后处理得到。实验数据信息如 表2 表2 所示。 注： G、C、E、R、J分别代表GPS、BDS、Galileo、GLONASS、QZSS卫星。 统计项 测量型设备实验 低成本设备实验 采集地点 上海市 武汉市 采集设备 NovAtel SPAN u-blox 采集时间 2022-03-09 2019-08-01 数据时长/min 20 50 采样间隔/s 1 1 实验场景 城市车载 城市车载 解算频点 G：L1/L2 C：B1/B2 E：E1/E5a R：L1/L2 J：L1/L2 G：L1 C：B1 E：E1 3.1　测量型设备实验 2 图6 图6 给出了测量型设备实验中，FGO、EKF与单历元最小二乘解（single epoch， SEP）的平面定位误差累积分布，可知在90%处误差分别为3.64 m、4.29 m、7.97 m，SEP定位结果精度与FGO、EKF两种方法存在明显差距。受限于模型结构强度，SEP方法无法通过边缘化或状态方程传递历史观测信息，即在输入数据上SEP与FGO、EKF存在实质差别，定位精度显然较差。因此，后续主要对FGO与EKF两种方法进行对比分析。 图7 图7 给出了该组实验的定位误差序列图，可以看出相较于EKF方法，FGO在东（east，E）、北（north，N）、天（up，U）3个方向上的定位精度分别提升了51%、16%、2%，三维定位精度提升了11.29%。在误差序列出现异常值的情况下，FGO充分利用了历元间观测的相关性，具备更好的抗差性能，能够获得更为平滑的定位结果。进一步对两种方法的抗差结果进行分析，各卫星各频点上的最小可探测误差平均值分布情况如 图8 图8 所示，EKF与FGO的MDB平均值分别为2.340 m与2.061 m，由于在观测时间相关性上的优势，FGO的MDB指标整体均优于EKF，表明具有更好的内部可靠性。 此外，本文将位置参考真值代入观测方程 式（10） 式（10） 中，并利用基站数据进行站间单差，消除卫星端和大气残余误差，从而获得真实的残差作为参考，其中大于10 m的残差标记为粗差。 图7 图7 中700~800 s时间段定位误差较大，表明粗差出现较为频繁，从中选取某一历元并绘制其验后残差。实验一验后残差序列如 图9 图9 所示，L1频点的真实残差按大到小依次排序，其中前7颗卫星残差超过10 m，即为粗差，而剩余卫星残差小于10 m，即为正常。以此为参考，FGO中仅有C13卫星被误判为粗差，EKF存在6颗卫星被误判为粗差，且EKF残差与真实残差分布差异较大，表明EKF发生了比较严重的粗差转移现象，导致粗差探测失效，定位结果较差，而FGO能够抑制粗差转移，具有较好的粗差识别能力。从 图10 图10 中可以看到，EKF的定位误差序列突刺较多，而FGO利用窗口内冗余信息能够对异常观测起到较好的抑制作用，相较于EKF在三维定位精度上提升了25.99%，且定位误差序列更为平稳。 另外，本文还统计了粗差探测成功率（correct detection rate， CDR），计算式为： <math display="block" id="M53"><mi mathvariant="normal">粗差</mi><mi mathvariant="normal">探测</mi><mi mathvariant="normal">成功</mi><mi mathvariant="normal">率</mi><mo>=</mo><mfrac><mrow><mi mathvariant="normal">成功</mi><mi mathvariant="normal">探测</mi><mi mathvariant="normal">粗差</mi><mi mathvariant="normal">数</mi></mrow><mrow><mi mathvariant="normal">粗差</mi><mi mathvariant="normal">总数</mi></mrow></mfrac><mo>×</mo><mn mathvariant="normal">100</mn><mi mathvariant="normal">%</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080105&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080114&type= 71.54333496 8.80533314 表3 表3 为两种方法的定位精度、MDB和CDR统计值。其中，FGO的CDR为76.60%，较EKF提升了10.47%，具有更好的粗差正确识别能力。 2 模型 3 定位误差RMS/m 2 <math id="M54"><mover><mrow><mi mathvariant="normal">M</mi><mi mathvariant="normal">D</mi><mi mathvariant="normal">B</mi></mrow><mo>¯</mo></mover></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080111&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080109&type= 6.60400009 2.70933342 /m 2 CDR/% E N U EKF 1.908 2.951 4.033 2.340 66.13 FGO 0.934 2.473 3.941 2.061 76.60 3.2　低成本设备实验 2 低成本设备在复杂环境下的伪距噪声以及多路径和非视距信号相对严重，定位精度影响更为明显，整体结果统计如 表4 表4 所示。 2 模型 3 定位误差RMS/m 2 <math id="M55"><mover><mrow><mi mathvariant="normal">M</mi><mi mathvariant="normal">D</mi><mi mathvariant="normal">B</mi></mrow><mo>¯</mo></mover></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080111&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=45080109&type= 6.60400009 2.70933342 /m 2 CDR/% E N U EKF 0.859 1.098 2.207 5.565 66.43 FGO 0.825 0.925 1.482 3.753 70.27 图11 图11 给出了各卫星的MDB平均值，其中EKF与FGO分别为5.565 m、3.753 m，由于FGO的MDB值小于EKF，说明FGO的内部可靠性较好。进一步的，FGO的CDR为70.27%，相较于EKF的66.43%提升了3.84%。在2 200~2 900 s时间段内，实验车辆经过城市峡谷场景，粗差频繁发生，整体数据质量较差，FGO通过窗口内冗余观测信息和历元间约束，在一定程度上可以提升复杂场景下的抗差性能和定位精度。 4　结　语 1 针对城市复杂场景下GNSS多路径与非视距信号严重、粗差发生频繁、定位精度差等问题，本文构建了基于滑动窗口因子图优化的GNSS定位模型，融合了边缘化先验信息、多频多系统的伪距观测值与GNSS测速结果，通过窗口内历元间约束与冗余观测信息抑制粗差影响，并结合验后残差迭代分析进行粗差探测与剔除。在输入信息与解算策略一致的前提下，通过城市复杂场景实验对EKF与FGO两种方法进行对比分析。实验结果表明，相较于EKF方法，FGO的MDB指标减小了11.92%~32.56%，说明FGO具有更好的内部可靠性，且验后残差分布与真实残差符合程度更好，可以有效识别粗差并进行正确剔除，因此，FGO的CDR提升了3.84%~10.47%，三维定位精度提升了11.29%~25.99%。总体而言，对于城市复杂场景下的GNSS导航定位应用，FGO具备更好的抗差性能和定位精度，证明了方法的实用价值，有望取代现有基于单历元观测值的EKF模型。未来将在此基础上结合自适应因子、噪声估计等方法精化模型，提升FGO在城市复杂场景中的定位性能。 参 考 文 献 Pan L ， Zhang X ， Liu J ， et al . Performance Evaluation of Single-frequency Precise Point Positioning with GPS， GLONASS， BeiDou and Galileo ［J］. The Journal of Navigation ， 2017 ， 70 （ 3 ）： 465 - 482 . Performance Evaluation of Single-frequency Precise Point Positioning with GPS， GLONASS， BeiDou and Galileo Geng Jianghui ， Chang Hua ， Guo Jiang ， et al . Three Multi-frequency and Multi-system GNSS High-Precision Point Positioning Methods and Their Performance in Complex Urban Environment ［J］. Acta Geodaetica et Cartographica Sinica ， 2020 ， 49 （ 1 ）： 1 - 13 . Three Multi-frequency and Multi-system GNSS High-Precision Point Positioning Methods and Their Performance in Complex Urban Environment 耿江辉 ， 常华 ， 郭将 ， 等 . 面向城市复杂环境的3种多频多系统GNSS单点高精度定位方法及性能分析 ［J］. 测绘学报 ， 2020 ， 49 （ 1 ）： 1 - 13 . 面向城市复杂环境的3种多频多系统GNSS单点高精度定位方法及性能分析 Zhang Xiaohong ， Tao Xianlu ， Wang Yingzhe ， et al . MEMS-Enhanced Smartphone GNSS High-Precision Positioning for Vehicular Navigation in Urban Conditions ［J］. Geomatics and Information Science of Wuhan University ， 2022 ， 47 （ 10 ）： 1740 - 1749 . MEMS-Enhanced Smartphone GNSS High-Precision Positioning for Vehicular Navigation in Urban Conditions 张小红 ， 陶贤露 ， 王颖喆 ， 等 . 城市场景智能手机GNSS/MEMS融合车载高精度定位 ［J］. 武汉大学学报 （信息科学版） ， 2022 ， 47 （ 10 ）： 1740 - 1749 . 城市场景智能手机GNSS/MEMS融合车载高精度定位 Takasu T ， Yasuda A . Development of the Low-Cost RTK-GPS Receiver with an Open Source Program Package RTKLIB ［C］// International Symposium on GPS/GNSS ， International Convention Center Jeju Korea Seogwiposi， Korea ， 2009 . Development of the Low-Cost RTK-GPS Receiver with an Open Source Program Package RTKLIB Watson R M ， Gross J N . Robust Navigation in GNSS Degraded Environment Using Graph Optimization ［C］// The 30th International Technical Meeting of the Satellite Division of the Institute of Navigation （ION GNSS+ 2017），Portland， Oregon，USA， 2017 . Robust Navigation in GNSS Degraded Environment Using Graph Optimization Kschischang F R ， Frey B J ， Loeliger H A . Factor Graphs and the Sum-Product Algorithm ［J］. IEEE Transactions on Information Theory ， 2001 ， 47 （ 2 ）： 498 - 519 . Factor Graphs and the Sum-Product Algorithm Dellaert F ， Kaess M . Square Root SAM： Simultaneous Localization and Mapping via Square Root Information Smoothing ［J］. The International Journal of Robotics Research ， 2006 ， 25 （ 12 ）： 1181 - 1203 . Square Root SAM： Simultaneous Localization and Mapping via Square Root Information Smoothing Wang Jinke ， Zuo Xingxing ， Zhao Xiangrui ， et al . Review of Multi-source Fusion SLAM： Current Status and Challenges ［J］. Journal of Image and Graphics ， 2022 ， 27 （ 2 ）： 368 - 389 . Review of Multi-source Fusion SLAM： Current Status and Challenges 王金科 ， 左星星 ， 赵祥瑞 ， 等 . 多源融合SLAM的现状与挑战 ［J］. 中国图象图形学报 ， 2022 ， 27 （ 2 ）： 368 - 389 . 多源融合SLAM的现状与挑战 Suzuki T . First Place Award Winner of the Smartphone Decimeter Challenge： Global Optimization of Position and Velocity by Factor Graph Optimization ［C］// The 34th International Technical Meeting of the Satellite Division of the Institute of Navigation （ION GNSS+ 2021），LouisSt.， Missouri， 2021 . First Place Award Winner of the Smartphone Decimeter Challenge： Global Optimization of Position and Velocity by Factor Graph Optimization Suzuki T . 1st Place Winner of the Smartphone Decimeter Challenge： Two-Step Optimization of Velocity and Position Using Smartphone’s Carrier Phase Observations ［C］// The 35th International Technical Meeting of the Satellite Division of the Institute of Navigation ， Denver， USA ， 2022 . 1st Place Winner of the Smartphone Decimeter Challenge： Two-Step Optimization of Velocity and Position Using Smartphone’s Carrier Phase Observations Watson R M ， Gross J N . Evaluation of Kinematic Precise Point Positioning Convergence with an Incremental Graph Optimizer ［C］// 2018 IEEE/ION Position ， Location and Navigation Symposium （PLANS）， Monterey， CA， USA ， 2018 . Evaluation of Kinematic Precise Point Positioning Convergence with an Incremental Graph Optimizer Wen W ， Hsu L T . Towards Robust GNSS Positioning and Real-time Kinematic Using Factor Graph Optimization ［C］// 2021 IEEE International Conference on Robotics and Automation （ICRA） ， Xi’an， China ， 2021 . Towards Robust GNSS Positioning and Real-time Kinematic Using Factor Graph Optimization Gao H ， Li H ， Huo H ， et al . Robust GNSS Real-Time Kinematic with Ambiguity Resolution in Factor Graph Optimization ［C］// 2022 International Technical Meeting of the Institute of Navigation ， Long Beach， California ， 2022 . Robust GNSS Real-Time Kinematic with Ambiguity Resolution in Factor Graph Optimization Xu Lang ， Chen Xiaodong ， Huang Pan . An Improved GNSS Positioning Algorithm Based on Factor Graph Optimization ［C］// The 13rd China Satellite Navigation Conference： Autonomous Navigation and Intelligent Operation ， Beijing，China ， 2022 . An Improved GNSS Positioning Algorithm Based on Factor Graph Optimization 徐浪 ， 陈晓东 ， 黄攀 . 基于因子图优化的GNSS定位算法研究 ［C］// 第十三届中国卫星导航年会： 自主导航与智能运行 ， 中国北京 ， 2022 . 基于因子图优化的GNSS定位算法研究 Bai X ， Wen W ， Hsu L T . Time-Correlated Window-Carrier-Phase-Aided GNSS Positioning Using Factor Graph Optimization for Urban Positioning ［J］. IEEE Transactions on Aerospace and Electronic Systems ， 2022 ， 58 （ 4 ）： 3370 - 3384 . Time-Correlated Window-Carrier-Phase-Aided GNSS Positioning Using Factor Graph Optimization for Urban Positioning Jiang Y ， Gao Y ， Ding W ， et al . GNSS Precise Positioning for Smartphones Based on the Integration of Factor Graph Optimization and Solution Separation ［J］. Measurement ， 2022 ， 203 ： 111924 . GNSS Precise Positioning for Smartphones Based on the Integration of Factor Graph Optimization and Solution Separation Sibley G ， Matthies L ， Sukhatme G . Sliding Window Filter with Application to Planetary Landing ［J］. Journal of Field Robotics ， 2010 ， 27 （ 5 ）： 587 - 608 . Sliding Window Filter with Application to Planetary Landing 10.13203/j.whugis20230203 P228 A 1671-8860（2023）07-1050-08 张小红,张元泰,朱锋.城市复杂场景下GNSS定位的因子图优化方法及其抗差性能分析[J].武汉大学学报（信息科学版）,2023,48(7):1050-1057.DOI: 10.13203/j.whugis20230203 http://dx.doi.org/10.13203/j.whugis20230203 张小红,张元泰,朱锋.城市复杂场景下GNSS定位的因子图优化方法及其抗差性能分析[J].武汉大学学报（信息科学版）,2023,48(7):1050-1057.DOI:10.13203/j.whugis20230203 ZHANG Xiaohong,ZHANG Yuantai,ZHU Feng.Factor Graph Optimization for Urban Environment GNSS Positioning and Robust Performance Analysis[J].Geomatics and Information Science of Wuhan University,2023,48(7):1050-1057.DOI: 10.13203/j.whugis20230203 http://dx.doi.org/10.13203/j.whugis20230203 ZHANG Xiaohong,ZHANG Yuantai,ZHU Feng.Factor Graph Optimization for Urban Environment GNSS Positioning and Robust Performance Analysis[J].Geomatics and Information Science of Wuhan University,2023,48(7):1050-1057.DOI:10.13203/j.whugis20230203 国家重点研发计划（2020YFB0505803）；国家自然科学基金（42104021）；湖北省科技重大项目（2021AAA010）；湖北珞珈实验室专项基金（2201000038）。 2023-06-04 2023-07-05 2023-06-29 39072592 武汉大学学报（信息科学版） 7 48