李雪珍 章浙濤 劉歡 袁海軍

0 引言

全球?qū)Ш叫l(wèi)星系統(tǒng)(GNSS)包括全球的、區(qū)域的衛(wèi)星導(dǎo)航系統(tǒng)以及增強(qiáng)系統(tǒng)，如美國(guó)的全球定位系統(tǒng)(GPS)、中國(guó)的北斗衛(wèi)星導(dǎo)航系統(tǒng)(BDS)、俄羅斯的全球衛(wèi)星導(dǎo)航系統(tǒng)(GLONASS)、歐洲的伽利略衛(wèi)星定位系統(tǒng)(Galileo)、日本的準(zhǔn)天頂系統(tǒng)(QZSS)等.GNSS在全球或區(qū)域范圍內(nèi)提供定位、導(dǎo)航和授時(shí)(PNT)服務(wù)，被廣泛應(yīng)用于自然災(zāi)害監(jiān)測(cè)、變形監(jiān)測(cè)、航空運(yùn)輸、防震減災(zāi)、車輛駕駛等領(lǐng)域.BDS是中國(guó)自主研發(fā)、獨(dú)立運(yùn)行的衛(wèi)星導(dǎo)航系統(tǒng)，其發(fā)展分為3個(gè)階段：北斗衛(wèi)星導(dǎo)航試驗(yàn)系統(tǒng)(BDS-1)、區(qū)域北斗系統(tǒng)(BDS-2)和全球北斗系統(tǒng)(BDS-3)[1-2].如今BDS-3已經(jīng)建成，由3顆地球靜止軌道(GEO)衛(wèi)星、3顆傾斜地球同步軌道(IGSO)衛(wèi)星和24顆中圓地球軌道(MEO)衛(wèi)星組成，于2020年向全球用戶提供導(dǎo)航定位服務(wù)[3].BDS-3在繼承BDS-2原有B1I和B3I信號(hào)的同時(shí)，新增了B1C、B2a和B2b頻點(diǎn)[4].目前BDS能接收更多的可見(jiàn)衛(wèi)星信號(hào)，得到更優(yōu)的有效觀測(cè)數(shù)據(jù)，實(shí)現(xiàn)良好的衛(wèi)星幾何構(gòu)型，進(jìn)而提升定位精度.BDS全球服務(wù)水平和高程定位精度達(dá)到10 m，在亞太地區(qū)可達(dá)到5 m，系統(tǒng)服務(wù)可用性優(yōu)于95%[5].

GNSS實(shí)際應(yīng)用環(huán)境是復(fù)雜且惡劣的.在城市峽谷環(huán)境中，高大建筑物的表面及街道上的樹(shù)木和車輛等物體使衛(wèi)星信號(hào)產(chǎn)生一定的反射和折射現(xiàn)象[6].衛(wèi)星信號(hào)的衰減、中斷導(dǎo)致可見(jiàn)衛(wèi)星數(shù)量的減少和衛(wèi)星幾何分布不理想，出現(xiàn)較多較為嚴(yán)重的粗差.少量的粗差，能對(duì)參數(shù)估計(jì)的結(jié)果造成干擾，進(jìn)而嚴(yán)重影響GNSS定位的精度和可靠性.目前用戶對(duì)于定位的實(shí)時(shí)性、精度、可靠性的需求急速增加.在多頻多系統(tǒng)和城市峽谷環(huán)境等復(fù)雜條件下，如何有效抵御粗差的影響是GNSS數(shù)據(jù)處理中的棘手問(wèn)題，更是數(shù)據(jù)質(zhì)量控制的關(guān)鍵點(diǎn).為剔除粗差，提高定位精度，北斗/GNSS數(shù)據(jù)的粗差處理主要可概括為兩大類方法：基于均值漂移的粗差處理方法[7-8]和基于方差膨脹的粗差處理方法[9].

本文以最常用的數(shù)據(jù)探測(cè)及穩(wěn)健估計(jì)為例闡述兩類粗差處理方法的流程，并研究給出了適用于復(fù)雜條件的閾值方案，最后利用實(shí)測(cè)數(shù)據(jù)的GPS+BDS組合定位進(jìn)行驗(yàn)證.

1 基于均值漂移的粗差處理

GNSS函數(shù)模型中的粗差通常會(huì)導(dǎo)致最小二乘估計(jì)存在偏差.基于均值漂移的粗差處理即粗差納入函數(shù)模型，認(rèn)為含有粗差的觀測(cè)值，其期望發(fā)生平移[10].最早提出且目前使用最為廣泛的基于均值漂移的粗差處理方法是數(shù)據(jù)探測(cè)(data snooping)[7].

具體而言，數(shù)據(jù)探測(cè)依賴于在原假設(shè)模型和一組備選假設(shè)模型之間進(jìn)行假設(shè)檢驗(yàn)，要求檢驗(yàn)統(tǒng)計(jì)量必須遵循正態(tài)分布、τ分布、x2分布和F分布[11-12].此外，由于觀測(cè)模型的幾何形狀[11]、假設(shè)之間的可分離性[13]、選擇的統(tǒng)計(jì)量[14]、預(yù)先確定的臨界值[15]等原因，不可避免地會(huì)產(chǎn)生漏檢、誤報(bào)和錯(cuò)誤識(shí)別，因此并不能完全消除偏差并確保最終參數(shù)估計(jì)的無(wú)偏性.數(shù)據(jù)探測(cè)最初只是針對(duì)一個(gè)粗差研究的.在復(fù)雜條件下,一組觀測(cè)值包含多個(gè)粗差時(shí)，數(shù)據(jù)探測(cè)只能一個(gè)接一個(gè)地探測(cè)，難以滿足實(shí)際的應(yīng)用需求.而多維粗差的同時(shí)定位與定值(ELGE法)不僅能確定多個(gè)粗差的位置，還能求出各個(gè)粗差的數(shù)值大小[16].部分最小二乘平差把觀測(cè)值是否含有粗差分成兩組，在不含粗差的一組里實(shí)施最小二乘平差，同樣能求出多個(gè)粗差的位置與數(shù)值[17].與上述以殘差為對(duì)象的方法不同的是擬準(zhǔn)檢定法，它依據(jù)觀測(cè)值的真誤差判斷粗差的位置[18].在觀測(cè)值相互獨(dú)立且等權(quán)時(shí)，以上幾種方法可以認(rèn)為是等價(jià)的[19-20].

Baarda[7]最早提出測(cè)量系統(tǒng)的可靠性理論.可靠性分析為粗差探測(cè)提供了相應(yīng)的理論基礎(chǔ).可靠性分析的3項(xiàng)重要指標(biāo)分別是最小可檢測(cè)偏差(MDB)、內(nèi)部可靠性和外部可靠性[21].其中MDB指在一定正確檢驗(yàn)概率條件下可被檢測(cè)到的最小偏差的絕對(duì)值.內(nèi)部可靠性指的是發(fā)現(xiàn)粗差的能力，外部可靠性是指不可發(fā)現(xiàn)粗差對(duì)平差結(jié)果的影響程度[22].對(duì)于多個(gè)備選假設(shè)的情況，最小可分離偏差(MSB)和最小可識(shí)別偏差(MIB)被研究和定義，它們可衡量與另一個(gè)備選假設(shè)成功分離或在多個(gè)備選假設(shè)中以一定的概率正確識(shí)別的最小偏差的大小[13,23].將MDB代入?yún)?shù)解得到外部可靠性，以衡量未檢測(cè)到的粗差對(duì)估計(jì)結(jié)果的影響顯著性.其他相關(guān)的可靠性測(cè)量也被研究，如可靠性數(shù)據(jù)指標(biāo)[24]和可控性(區(qū)域三角)[25].目前，基于均值漂移的粗差處理方法已被廣泛應(yīng)用于各種領(lǐng)域，如大地測(cè)量網(wǎng)的質(zhì)量控制和完好性監(jiān)測(cè)等.

2 基于方差膨脹的粗差處理

另一類基于方差膨脹的粗差處理方法認(rèn)為含有粗差的觀測(cè)值的方差發(fā)生變化，但期望不變[5].基于方差膨脹的粗差處理方法的研究主要集中在穩(wěn)健估計(jì)(抗差估計(jì))[22,26].穩(wěn)健估計(jì)將粗差納入隨機(jī)模型，通過(guò)逐次迭代平差結(jié)果，不斷改變觀測(cè)值的權(quán)，最終使含有粗差的觀測(cè)值權(quán)為零.穩(wěn)健估計(jì)主要依賴等價(jià)權(quán)函數(shù)來(lái)檢測(cè)和排除GNSS觀測(cè)數(shù)據(jù)中的粗差[27].穩(wěn)健估計(jì)主要包括M估計(jì)(廣義極大似然估計(jì))[28]、R估計(jì)(秩檢驗(yàn)估計(jì))[29]、L1范數(shù)估計(jì)[30]、最小截平方和估計(jì)[31]、符號(hào)約束的穩(wěn)健估計(jì)[23].M估計(jì)的抗差性和有效性主要取決于所采用的參數(shù)的可靠性、等價(jià)權(quán)函數(shù)及其臨界值的合理性[32]，因此不同權(quán)函數(shù)消除或減弱粗差影響的能力不盡相同.已有許多不同的權(quán)函數(shù)能高效和準(zhǔn)確處理獨(dú)立觀測(cè)值和相關(guān)觀測(cè)值，如雙重法[33]、Huber法[27]、Hampel法[34]、Danish法[35]、IGGI法[36]和IGGIII法[9,37].其中IGGI作為國(guó)內(nèi)早期經(jīng)典的質(zhì)量控制算法，它依據(jù)殘差的假設(shè)檢驗(yàn)結(jié)果調(diào)整權(quán)陣來(lái)減弱含粗差的觀測(cè)值對(duì)定位的影響.在IGGI的基礎(chǔ)上，利用等價(jià)權(quán)原理，進(jìn)一步發(fā)展到相關(guān)觀測(cè)值的IGGIII方法.同樣基于等價(jià)權(quán)，Yang等[38]利用方差分量的穩(wěn)健估計(jì)來(lái)減弱粗差的影響.早期，在隨機(jī)模型中的觀測(cè)值的方差被認(rèn)為是相等的[39].后來(lái)采用異方差假設(shè)來(lái)代替GNSS觀測(cè)的不切實(shí)際的同方差.通常采用衛(wèi)星高度角和載噪比作為觀測(cè)值異方差假設(shè)的兩個(gè)主要指標(biāo).其中高度角和觀測(cè)值精度之間的關(guān)系可以用余弦函數(shù)或指數(shù)函數(shù)來(lái)充分描述[40].由于載噪比與GNSS觀測(cè)值記錄在同一個(gè)跟蹤環(huán)路中，因此載噪比與GNSS觀測(cè)值精度高度一致[41].

此外，雙因子方差膨脹模型和雙因子等價(jià)權(quán)模型同樣有效地控制粗差對(duì)參數(shù)估值的影響[42].可從污染誤差模型入手，研究最小均方差準(zhǔn)則的參數(shù)穩(wěn)健估計(jì)[43].楊元喜[44]根據(jù)多種抗差濾波及其性質(zhì)構(gòu)建抗差自適應(yīng)濾波理論體系.同樣地，基于方差膨脹的粗差處理方法經(jīng)常被用于大地測(cè)量數(shù)據(jù)處理，例如變形分析[29]、基準(zhǔn)變換[45]、GNSS完好性監(jiān)測(cè)[46]、卡爾曼濾波[47].但是它們?cè)谀Ｐ徒Y(jié)構(gòu)不強(qiáng)、粗差較多的情況下，并不能完全消除粗差的影響.因此在多系統(tǒng)多頻和城市峽谷環(huán)境等復(fù)雜條件下，對(duì)粗差處理方法的研究顯得尤為重要.

3 兩種粗差處理方法

3.1 數(shù)據(jù)探測(cè)

在復(fù)雜條件下，GNSS觀測(cè)值通常不可避免地包含粗差.采用數(shù)據(jù)探測(cè)的方法對(duì)粗差進(jìn)行后驗(yàn)處理，具體分為以下3個(gè)步驟[7-8]：

1)探測(cè)

采用整體檢驗(yàn)，構(gòu)造Tq統(tǒng)計(jì)量：

(1)

2)識(shí)別

采用標(biāo)準(zhǔn)化殘差構(gòu)造wi檢驗(yàn)統(tǒng)計(jì)量:

(2)

3)調(diào)節(jié)

最后利用剩余的觀測(cè)值重新進(jìn)行最小二乘估計(jì)得到可靠的結(jié)果，并將其輸出.

3.2 穩(wěn)健估計(jì)

選權(quán)迭代法作為M估計(jì)的重要方法，同時(shí)是穩(wěn)健估計(jì)中最為常用且計(jì)算簡(jiǎn)便的方法，其基本步驟如下：

1)建立數(shù)學(xué)模型

(3)

2)按最小二乘法求解參數(shù)估值及其殘差

(4)

(5)

3)求解觀測(cè)值的等價(jià)權(quán)矩陣，應(yīng)用抗差最小二乘迭代計(jì)算，即

(6)

選定某一值ω作為迭代閾值(后續(xù)實(shí)驗(yàn)所討論的閾值)，當(dāng)?shù)趉次與第k-1次迭代所得的估值之差的絕對(duì)值小于等于該值，則停止迭代.

(7)

4)最后求解

(8)

(9)

選權(quán)迭代法的關(guān)鍵在于選擇權(quán)函數(shù)，最經(jīng)典的是IGGⅢ方法.IGGⅢ方法的權(quán)函數(shù)pi確定如下:

(10)

(11)

(12)

式中，q為多余觀測(cè)數(shù).

本文中討論了面向復(fù)雜條件時(shí)ω、k0、k1的選取，由于過(guò)多地剔除觀測(cè)值，會(huì)導(dǎo)致定位精度降低，因此可以通過(guò)迭代閾值(ω)來(lái)控制迭代次數(shù)，使穩(wěn)健估計(jì)既具抗差性，又高效.

4 實(shí)際數(shù)據(jù)分析

為驗(yàn)證本文提出的粗差處理方案的有效性并評(píng)價(jià)定位精度，本文選取2021年年積日第288天24小時(shí)的5組西南地區(qū)滑坡監(jiān)測(cè)數(shù)據(jù)進(jìn)行實(shí)驗(yàn).5組數(shù)據(jù)分別命名為Test1、Test2、Test3、Test4和Test5，其基線長(zhǎng)度分別約92 m、143 m、53 m、133 m和75 m.實(shí)驗(yàn)主要利用GPS的L1+L2觀測(cè)值和BDS的B1+B2 觀測(cè)值進(jìn)行處理，電離層延遲和對(duì)流層延遲分別利用Klobuchar模型和Saastamoinen模型進(jìn)行改正，模糊度固定采用最小二乘模糊度降相關(guān)平差(LAMBDA)方法.具體的解算策略如表1所示.

表1 解算策略

表2 粗差處理方案

表3 觀測(cè)數(shù)據(jù)的模糊度固定

限于篇幅，本文僅給出Test1的詳細(xì)數(shù)據(jù)處理結(jié)果分析.在衛(wèi)星定位中，可視衛(wèi)星數(shù)越多，衛(wèi)星幾何構(gòu)型越穩(wěn)定.圖1為Test1觀測(cè)數(shù)據(jù)的GPS+BDS的可視衛(wèi)星數(shù)和定位精度因子(PDOP).可以看出，觀測(cè)數(shù)據(jù)的可視衛(wèi)星數(shù)平均值為21，但衛(wèi)星數(shù)存在明顯的波動(dòng)，接收衛(wèi)星信號(hào)不太穩(wěn)定，表明觀測(cè)質(zhì)量不佳.所有歷元的PDOP值都大于1，其PDOP均值約1.6，有些歷元的PDOP值超過(guò)3，因此整體看衛(wèi)星的空間幾何分布良好，但有時(shí)存在較差的空間結(jié)構(gòu).從整體看，隨著衛(wèi)星數(shù)的減少，PDOP的值隨之增加.圖2為Test1觀測(cè)數(shù)據(jù)的衛(wèi)星天空視圖.可以發(fā)現(xiàn)GPS和BDS存在部分低高度角衛(wèi)星，導(dǎo)致衛(wèi)星信號(hào)存在反射、衍射和被阻擋的現(xiàn)象.因此在數(shù)據(jù)處理過(guò)程中設(shè)置了截止高度角為15°，剔除了不符合要求的衛(wèi)星.

圖1 Test1數(shù)據(jù)的衛(wèi)星數(shù)和PDOP值Fig.1 Satellite number and PDOP value of Test1

圖2 Test 1數(shù)據(jù)的衛(wèi)星天空視圖Fig.2 Skyplot of Test 1

圖3和圖4分別給出了A方案的雙差偽距殘差和雙差載波殘差.由圖3可知，GPS和BDS的偽距殘差絕大多數(shù)集中在10 m之內(nèi)，其中GPS L2的偽距殘差較大，可以看到在某些歷元存有偽距粗差.由圖4可看出，GPS和BDS的載波殘差絕大多數(shù)集中在0.1 m之內(nèi).其中載波粗差沒(méi)有完全剔除，而且較大的載波粗差會(huì)影響歷元其他觀測(cè)值，對(duì)定位結(jié)果產(chǎn)生嚴(yán)重影響.總之，需要加入合適的粗差處理方法來(lái)提高定位精度.

圖3 A方案的雙差偽距殘差Fig.3 Double-differenced pseudorange residuals of scheme A

圖5、6、7分別描繪了各個(gè)方案的定位誤差結(jié)果.需要說(shuō)明的是，定位結(jié)果中包含了固定解和浮點(diǎn)解的結(jié)果，以進(jìn)一步整體評(píng)估其定位性能.由圖5可以看出，無(wú)任何粗差處理方法的A方案定位結(jié)果是不太理想的，尤其是U方向的波動(dòng)較為明顯.對(duì)比圖5，圖6的3個(gè)數(shù)據(jù)探測(cè)方案都能探測(cè)到觀測(cè)值的粗差，并進(jìn)行粗差剔除，其中B方案在U方向上絕大多數(shù)歷元的定位誤差都比其他兩個(gè)方案小.圖7中大迭代閾值(ω)的穩(wěn)健估計(jì)方案比小迭代閾值(ω)的穩(wěn)健估計(jì)方案在相同的閾值(k0、k1)下有更為明顯的收斂.其中在小閾值(k0、k1)時(shí)，E1方案在U方向上明顯比E方案的誤差小.同樣E1方案對(duì)比其他穩(wěn)健估計(jì)方案有更為明顯的粗差處理效果，在衛(wèi)星定位中表現(xiàn)了更好的定位性能.因此小閾值(k0、k1)結(jié)合大迭代閾值(ω)的穩(wěn)健估計(jì)方案更適用于復(fù)雜條件的粗差處理.

圖4 A方案的雙差載波殘差Fig.4 Double-differenced carrier residuals of scheme A

圖5 A方案的定位誤差Fig.5 Positioning error of scheme A

圖6 數(shù)據(jù)探測(cè)方案的定位誤差Fig.6 Positioning error of data snooping

表5 穩(wěn)健估計(jì)方案的定位RMSE統(tǒng)計(jì)

為驗(yàn)證各方案粗差處理的有效性，表6和表7分別列出數(shù)據(jù)探測(cè)與穩(wěn)健估計(jì)的定位可用性，可以看到數(shù)據(jù)探測(cè)方案的可用性在任何的水平分量都優(yōu)于A方案.同樣穩(wěn)健估計(jì)的E、F方案的定位可用性略差于A方案，但其他方案的定位可用性比A方案有所提升，其中E1方案在水平分量為2.0 m時(shí)，定位可用性達(dá)到96.50%.上述結(jié)果同樣證明大迭代閾值(ω)配合小閾值(k0、k1)的穩(wěn)健估計(jì)方案可以更好地識(shí)別粗差，對(duì)觀測(cè)值權(quán)重的調(diào)整更為合理.

圖7 穩(wěn)健估計(jì)方案的定位誤差Fig.7 Positioning error of robust estimation

表6 數(shù)據(jù)探測(cè)方案的定位可用性

表7 穩(wěn)健估計(jì)方案的定位可用性

5 結(jié)語(yǔ)

觀測(cè)數(shù)據(jù)中含有異常粗差時(shí)，會(huì)對(duì)定位的精度和模糊度收斂產(chǎn)生嚴(yán)重的影響.本文系統(tǒng)研究了面向復(fù)雜條件時(shí)粗差處理的方法，針對(duì)復(fù)雜條件下的北斗/GNSS數(shù)據(jù)提出多種粗差處理方案，并以五組實(shí)驗(yàn)數(shù)據(jù)分析各方案的模糊度收斂速率，著重以一組數(shù)據(jù)對(duì)比分析各方案的粗差處理效果和定位性能.

綜上，數(shù)據(jù)探測(cè)和穩(wěn)健估計(jì)的處理粗差效果相當(dāng)，但由于數(shù)據(jù)探測(cè)的粗差處理相對(duì)更精細(xì)化，因此在本文實(shí)驗(yàn)中數(shù)據(jù)探測(cè)可以表現(xiàn)出更好的適用性.