虢 盛，徐進(jìn)立，張紹成，茍鴻飛

（中國(guó)地質(zhì)大學(xué)（武漢） 地理與信息工程學(xué)院，武漢 430078）

0 引言

全球定位系統(tǒng)（global positioning system, GPS）及北斗衛(wèi)星導(dǎo)航系統(tǒng)（BeiDou navigation satellite system，BDS）進(jìn)行組合導(dǎo)航定位時(shí)，不僅能通過增加衛(wèi)星數(shù)量來提高導(dǎo)航定位精度，還能通過不同系統(tǒng)觀測(cè)數(shù)據(jù)相互校驗(yàn)來保證用戶導(dǎo)航定位的可靠性。

20世紀(jì)60年代，荷蘭巴爾達(dá)（Baarda）教授提出了可靠性研究的基礎(chǔ)理論，并在假設(shè)單一粗差和驗(yàn)前單位權(quán)中誤差已知的情況下，推導(dǎo)得出了著名的“數(shù)據(jù)探測(cè)法”[1]，文獻(xiàn)[2]將 Baarda 的數(shù)據(jù)探測(cè)理論應(yīng)用于GPS基線解算中，提出了基于三種 GPS基線解算模型下的內(nèi)部可靠性檢驗(yàn)方法。文獻(xiàn)[3]提出將粗差視為函數(shù)模型一部分的重要思想，進(jìn)一步完善了可靠性研究理論。文獻(xiàn)[4]探究了觀測(cè)值相關(guān)條件下的內(nèi)部可靠性原理，提出一種新的可靠性指標(biāo)，解決了觀測(cè)多余度值域超過[0,1]區(qū)間后的可靠性度量問題。

在GNSS接收機(jī)自主完好性評(píng)估方面，文獻(xiàn)[5]研究了 GPS單系統(tǒng)及多系統(tǒng)組合增強(qiáng)型接接收機(jī)自主完好性監(jiān)測(cè)（advanced receiver autonomous integrity monitoring, ARAIM）可用性性能。文獻(xiàn)[6]則利用GPS/BDS兼容接收機(jī)雙頻實(shí)測(cè)數(shù)據(jù)及BDS布滿星座假設(shè)下的仿真數(shù)據(jù)，通過改進(jìn)不同衛(wèi)星頻點(diǎn)特征向量選擇方法，提高了基于奇偶矢量法的接收機(jī)自主完好性性能，文獻(xiàn)[7]將GPS/BDS衛(wèi)星按星座分為兩組，驗(yàn)證了最優(yōu)加權(quán)平均解（optimal weighted average solution，OWAS）算法對(duì)30 m以上誤差的雙星多故障檢測(cè)的可行性。

本文基于GNSS接收機(jī)中GPS和BDS觀測(cè)數(shù)據(jù)的相互校驗(yàn)原理，來實(shí)現(xiàn)單一粗差探測(cè)。首先研究基于偽距殘差的內(nèi)部可靠性模型及經(jīng)典接收機(jī)自主完好性監(jiān)測(cè)（receiver autonomous integrity monitoring, RAIM）在 GPS/BDS雙模導(dǎo)航定位中的基本理論；然后根據(jù)GPS/BDS偽距觀測(cè)殘差方差特性，基于均值漂移模型提出一種系統(tǒng)間均值校驗(yàn)的接收機(jī)自主完好性粗差探測(cè)方法；最后通過實(shí)測(cè)數(shù)據(jù)的異常模擬，與經(jīng)典 RAIM 法進(jìn)行對(duì)比，來驗(yàn)證了本文算法的可行性和優(yōu)越性。

1 內(nèi)部可靠性數(shù)學(xué)模型

1.1 內(nèi)部可靠性模型

可靠性檢驗(yàn)理論作為用戶導(dǎo)航定位完好性監(jiān)測(cè)的重要研究?jī)?nèi)容，主要分為外部可靠性和內(nèi)部可靠性兩方面。外部可靠性主要研究未剔除粗差對(duì)結(jié)果影響的大小，即抵抗粗差的能力[8]。內(nèi)部可靠性是指在平差過程中發(fā)現(xiàn)并剔除粗差的能力[9]。

本文重點(diǎn)討論單一粗差假設(shè)下的GPS/BDS偽距殘差內(nèi)部可靠性檢驗(yàn)方法。當(dāng)偽距觀測(cè)值中不存在粗差時(shí)，H0假設(shè)為

當(dāng)偽距觀測(cè)值中存在粗差時(shí)，H1假設(shè)為

基于上述假設(shè)，構(gòu)建內(nèi)部可靠性檢驗(yàn)?zāi)Ｐ蛯ふ铱赡艹霈F(xiàn)的衛(wèi)星導(dǎo)航定位服務(wù)異常事件。

當(dāng)服從H0假設(shè)時(shí)，線性化誤差模型為

式中：l為n ×1維觀測(cè)值向量；v為n×1維改正數(shù)向量；B為n×m維觀測(cè)值矩陣；x︿為m×1維待估向量。觀測(cè)值權(quán)陣P為n×n階矩陣。

根據(jù)最小二乘準(zhǔn)則vTPv=min可知

參數(shù)估值為

對(duì)應(yīng)改正數(shù)為

改正數(shù)的協(xié)因數(shù)陣為

根據(jù)可靠性理論可得

通常將矩陣R對(duì)角線上的值作為可靠性量度（ri=Rii），其數(shù)值越小，表明對(duì)應(yīng)的觀測(cè)值可靠性越差；反之可靠性越好。當(dāng)觀測(cè)值之間相互獨(dú)立時(shí)，可靠性量度ri?[0,1]；而當(dāng)觀測(cè)值相關(guān)時(shí)，可靠性量度ri取值可能大于1。在觀測(cè)值相互獨(dú)立時(shí)，ri=1表明對(duì)應(yīng)觀測(cè)值為必要觀測(cè)；ri=0則表明對(duì)應(yīng)觀測(cè)為多余觀測(cè)[10]。

若觀測(cè)值向量無粗差，則?x的估計(jì)值是無偏的，用間接平差求解即可；若觀測(cè)值向量中包含粗差，那么上述的估計(jì)值是有偏的，即拒絕H0假設(shè)，接受H1假設(shè)，線性化誤差模型為

式中：H為n ×k維粗差向量；s︿為k×1維未知粗差[11]。

根據(jù)均值漂移函數(shù)模型，可得對(duì)應(yīng)法方程為

因此可知未知粗差的估值為

1.2 經(jīng)典粗差探測(cè)法

經(jīng)典 RAIM 算法主要有快照式算法、偽距殘差比較法和奇偶矢量法，且已有學(xué)者證明了三者的等價(jià)性[12]。這里介紹的經(jīng)典粗差探測(cè)法是在偽距殘差比較法基礎(chǔ)上，將 GPS/BDS觀測(cè)值按 2.1節(jié)總結(jié)的基于最小二乘的間接平差法處理，從而構(gòu)建卡方分布檢驗(yàn)統(tǒng)計(jì)量探測(cè)完成粗差探測(cè)的過程。統(tǒng)計(jì)量F定義

式中：n為GPS和BDS衛(wèi)星數(shù)之和；v為衛(wèi)星觀測(cè)值殘差；P為相應(yīng)的權(quán)陣。

觀測(cè)方程含測(cè)站三維坐標(biāo)及對(duì)應(yīng)GPS/BDS不同系統(tǒng)接收機(jī)鐘差 5個(gè)待估參數(shù)。若偽距觀測(cè)值無粗差，則理論上殘差向量滿足零均值的正態(tài)分布，進(jìn)而F滿足自由度為n-5的卡方分布。

根據(jù)用戶需求定義完好性誤警概率PFA即可確定對(duì)應(yīng)的檢驗(yàn)限值F，具體公式[13]為

式中：fχ(n?5)(x)為卡方分布的概率密度函數(shù)；為GPS與BDS偽距觀測(cè)方差。

1.3 基于均值漂移的粗差探測(cè)法

經(jīng)典粗差探測(cè)法是將粗差視為隨機(jī)模型一部分，有粗差的觀測(cè)值相比其他觀測(cè)值具有相同的期望、較大的方差，因而通常稱為“方差擴(kuò)大模型”[14]；此外，還可將粗差視為函數(shù)模型的一部分，含粗差的觀測(cè)值相比其他觀測(cè)值具有不同的期望，在方差相同、期望的不同條件下，觀測(cè)量會(huì)呈現(xiàn)非中心化分布，因此也稱為“均值漂移模型”[15]。

無粗差時(shí)，GPS/BDS偽距改正數(shù)vFull滿足零均值的正態(tài)分布，即

但在單個(gè)歷元有限可視衛(wèi)星條件下，偽距殘差均值μ通常不絕對(duì)等于零，其統(tǒng)計(jì)結(jié)果可表達(dá)為

式中：μ1、μ2分別為 GPS和 BDS偽距殘差期望，對(duì)應(yīng)殘差方差分別為和。

以武漢地區(qū)的國(guó)際 GNSS服務(wù)（International GNSS Service, IGS）九峰站（JFNG）2020年12月1日24小時(shí)GPS/BDS觀測(cè)數(shù)據(jù)為例，GPS和BDS L1和B2頻點(diǎn)殘差序列統(tǒng)計(jì)值如圖1所示。

圖1 JFNG 2020年11月30日GPS/BDS的殘差分布

由圖1可見，無粗差條件下GPS/BDS偽距殘差趨近于零均值；但GPS殘差方差σ12=2.41 m2而BDS方差σ22=4.44 m2，存在一定差異。

按均值漂移模型對(duì)GPS/BDS偽距觀測(cè)值分別平差，并構(gòu)建一個(gè)新的內(nèi)部可靠性檢驗(yàn)統(tǒng)計(jì)量，即偏離量K。K定義為BDS和GPS獨(dú)立偽距殘差分布期望值之差，計(jì)算方法為

進(jìn)一步構(gòu)建T統(tǒng)計(jì)量，檢驗(yàn)偏移量K值差異，并按式（15）計(jì)算對(duì)應(yīng)檢驗(yàn)門限值。通過比較檢驗(yàn)統(tǒng)計(jì)量T與檢驗(yàn)門限值大小即可完成故障檢測(cè)。

假定兩獨(dú)立樣本方差未知且不等情況下，采用薩特思韋特（Satterthwaite）近似法構(gòu)造T統(tǒng)計(jì)量[16]為

對(duì)應(yīng)自由度N表達(dá)為

式中，ni（i=1,2）為GPS/BDS衛(wèi)星數(shù)量。

采用均值漂移定律解算含粗差偽距值，雖然能大幅度減弱殘差樣本方差的放大，但在 GPS/BDS可視衛(wèi)星20余顆的條件下，殘余方差失真效應(yīng)不可忽略。本文基于GPS/BDS偽距殘差先驗(yàn)信息取代 GPS/BDS偽距殘差方差，根據(jù)圖1統(tǒng)計(jì)結(jié)果，設(shè)定GPS和BDS先驗(yàn)方差分別為2.5 m和4.5 m。

本文改進(jìn)內(nèi)部可靠性檢驗(yàn)方法概括如下：

1）建立假設(shè)檢驗(yàn)：H0認(rèn)為該歷元存在粗差；H1認(rèn)為該歷元不存在粗差。并假設(shè)第i顆衛(wèi)星存在粗差，此時(shí)式(9)中粗差向量H(i)=1,其余元素均為0。

2）根據(jù)式（9）至式（13）計(jì)算對(duì)應(yīng)的粗差估值s︿、偽距殘差v，并按 GPS/BDS方差先驗(yàn)信息，重構(gòu)GPS/BDS觀測(cè)值方差。最后根據(jù)式（21）至（22）式計(jì)算檢驗(yàn)統(tǒng)計(jì)量T。

3）循環(huán)n次（n表示 GPS/BDS衛(wèi)星數(shù)量之和），并選擇最大值Tmax與按式(15)計(jì)算的檢驗(yàn)門限值比較，Tmax較大則認(rèn)為該歷元存在粗差，并求得對(duì)應(yīng)的偽距殘差v 和粗差估值s︿ 。本文式(15)中誤警概率PFA≤0.01。

1.4 最不可探測(cè)衛(wèi)星

在內(nèi)部可靠性檢驗(yàn)中，最小可探測(cè)誤差（minimal detectable bias ,MDB）是一個(gè)很重要的指標(biāo)，其表示在一定的檢驗(yàn)功效以及一定的顯著水平下，可發(fā)現(xiàn)粗差的下界。

若只有一個(gè)粗差的單個(gè)備選假設(shè)下，其最小可探測(cè)誤差為

式中：Mi為第i顆衛(wèi)星對(duì)應(yīng)的MDB值；λ0為非中心化參數(shù)；σ0為驗(yàn)前單位權(quán)中誤差；P為權(quán)；Qvv為改正數(shù)的協(xié)因數(shù)陣。根據(jù)式（8）可靠性矩陣R，可將式（22）改寫為

從式（23）中可以看出，當(dāng)非中心化參數(shù)、驗(yàn)前單位權(quán)中誤差以及對(duì)應(yīng)觀測(cè)值的權(quán)一定時(shí)，M和可靠性量度值成反比。當(dāng)可靠性度量參數(shù)ir越大時(shí)，觀測(cè)值可靠性越高，所對(duì)應(yīng)的M值越小，進(jìn)而認(rèn)為MDB值最大的衛(wèi)星為最不可探測(cè)衛(wèi)星。

2 算例分析

選取2020年12月1日J(rèn)FNG站的GPS和BDS雙模觀測(cè)數(shù)據(jù)，設(shè)定衛(wèi)星截止高度角 10°，L1/C2信噪比大于30 dB，并摒棄軌道誤差較大的BDS 地球靜止軌道（geostationary Earth orbit,GEO）衛(wèi)星后，畫出當(dāng)天可見衛(wèi)星數(shù)和位置精度衰減因子（position dilution of precision, PDOP）值如圖2所示，該觀測(cè)時(shí)段內(nèi)可見GPS/BDS衛(wèi)星數(shù)分別在 6~13顆和 12~19顆之間，對(duì)應(yīng) PDOP值在 1~2之間， GPS/BDS可視衛(wèi)星具有較好的空間幾何分布。

圖2 JFNG站點(diǎn)可見衛(wèi)星數(shù)/PDOP值

基于改進(jìn)粗差探測(cè)方法計(jì)算無粗差情況下檢驗(yàn)效果如圖3所示，當(dāng)GPS/BDS偽距觀測(cè)值中未包含粗差時(shí)，T統(tǒng)計(jì)量整體小于檢驗(yàn)門限。表明當(dāng)觀測(cè)值中無粗差時(shí)，BDS和GPS系統(tǒng)殘差期望無明顯差異。

圖3 無粗差下T統(tǒng)計(jì)量

當(dāng)對(duì) BDS最不可探測(cè)衛(wèi)星 C2頻點(diǎn)上分別添加10、20、30 m粗差后，其檢驗(yàn)效果如圖4所示。

圖4 本文RAIM法添加10，20和30 m粗差后T統(tǒng)計(jì)量示意圖

統(tǒng)計(jì)分析結(jié)果表明，當(dāng)人為在最不可探測(cè)衛(wèi)星上添加大小為10 m的粗差后，共有1 381個(gè)歷元的檢驗(yàn)量T大于其對(duì)應(yīng)歷元下門限值；添加粗差量增至20 m后，觀測(cè)時(shí)間段內(nèi)有2 794個(gè)歷元的T值大于門限值，表明有95.5%的概率可以探測(cè)并剔除大小為20 m的粗差；添加粗差量再次增加至30 m，共有2 863個(gè)歷元無法通過檢驗(yàn)，對(duì)應(yīng)粗差探測(cè)成功百分比為99.4%。

如圖5所示，對(duì)比經(jīng)典 RAIM探測(cè)效率圖并將改進(jìn)方法與經(jīng)典 RAIM 方法探測(cè)成功百分比差異統(tǒng)計(jì)如表1。經(jīng)典探測(cè)法對(duì)20 m以下的微小粗差較不敏感，人為添加10、20 m粗差的探測(cè)成功百分比僅為 0和 3.2%。而本文的改進(jìn)方法對(duì) 10、20 m粗差探測(cè)成功率分別提高了 15.5%、36.7%，表明該方法較經(jīng)典 RAIM 法而言，能夠有效剔除量級(jí)為20 m的粗差。

圖5 經(jīng)典RAIM法添加10，20和30 m粗差后F統(tǒng)計(jì)量示意圖

表1 經(jīng)典RAIM算法不同粗差探測(cè)率統(tǒng)計(jì)表

3 結(jié)束語

本文通過分析 GPS/BDS偽距殘差期望和方差特性，提出了一種基于GPS/BDS偽距殘差均值差異的內(nèi)部可靠性檢驗(yàn)方法。利用偽距觀測(cè)值中存在粗差時(shí)，GPS/BDS獨(dú)立偽距殘差分布之間存在非中心化這一特點(diǎn)，構(gòu)造統(tǒng)計(jì)量T，在漏檢率PFA＜1%條件下確定對(duì)應(yīng)檢驗(yàn)門限值，判斷系統(tǒng)是否存在硬件級(jí)粗差。模擬粗差探測(cè)實(shí)驗(yàn)是在GPS/BDS實(shí)測(cè)數(shù)據(jù)上人為添加粗差，驗(yàn)證了本文所提出的可靠性檢驗(yàn)方法能夠準(zhǔn)確、高效地實(shí)現(xiàn)粗差探測(cè)。GPS/BDS雙系統(tǒng)不僅能夠提高定位精度，而且還可以進(jìn)一步保障導(dǎo)航數(shù)據(jù)的完好性。