一種基于粗差檢驗(yàn)整體平差的周跳探測與修復(fù)方法

趙祥偉，陳正宇

(中國能源建設(shè)集團(tuán)江蘇省電力設(shè)計(jì)院有限公司，江蘇 南京 211102)

0 引言

GNSS動態(tài)測量中，載波相位觀測值對于定位結(jié)果的影響較大，同時(shí)是系統(tǒng)性與長時(shí)間的影響[1-3]。為了得到高精度的定位結(jié)果，必須要準(zhǔn)確地修復(fù)觀測值中可能存在的周跳。但是，因?yàn)镚NSS工作的機(jī)動性較強(qiáng)，觀測環(huán)境變化較大，受到對流程、電離層、多路徑等因素的影響也各不相同，這對載波觀測值的質(zhì)量造成了多方面的影響[4]。因而，周跳的現(xiàn)象經(jīng)常發(fā)生。如何快速高效地探測修復(fù)周跳成為了測繪研究者們急需解決的問題。

隨著研究的深入，對于周跳探測修復(fù)工作的研究成果也越來越多。方法更是層出不窮，目前已有的檢測周跳的方法主要有多項(xiàng)式擬合法[5-6]、差分法[7]、多普勒頻移法[3,8]、相位/偽距組合法[9]、電離層殘差法[2]、濾波法[10]和站間或站星間差分法[11]。

本文選取現(xiàn)有的相位/偽距組合法，結(jié)合粗差探測[12-14]方法，對周跳探測方法進(jìn)行適當(dāng)?shù)母倪M(jìn)，通過實(shí)例分析得出該方法能夠探測出 2周以上的周跳值，且探測方法簡單，探測速度較快，探測結(jié)果可靠。

1 周跳檢驗(yàn)量構(gòu)成

1.1 偽距/載波組合法

在周跳探測的諸多方法中，一般常利用模糊度參數(shù)進(jìn)行周跳的探測修復(fù)。在信號觀測方程中，偽距觀測方程為[15]：

式(1)中，偽距觀測方程中不包括模糊度項(xiàng)，因而認(rèn)為偽距觀測值不會受到周跳的影響。相應(yīng)的載波相位觀測方程可表示為：

綜合式(1)和(2)，為分離出模糊度參數(shù)值，通常需要對其進(jìn)行組合，將相位觀測值直接減去偽距觀測值，消除接收機(jī)鐘差和衛(wèi)星鐘差等公共誤差，忽略次要因素的影響，如多路徑效應(yīng)、噪聲項(xiàng)等。相位偽距觀測方程作差，表示為：

式中：λ表示載波的波長；?表示載波相位觀測值；P表示偽距觀測值；I表示電離層延遲量；a表示整周模糊度；m?和mP分別表示載波相位觀測值和偽距觀測值的觀測值噪聲；δ?和δP分別表示載波相位觀測值和偽距觀測值殘差。

1.2 構(gòu)成周跳檢驗(yàn)量

如式(3)所示，經(jīng)過相位偽距觀測值作差后，方程里剩余有電離層殘差、多路徑殘差以及觀測噪聲等因素影響。當(dāng)歷元采樣率較小時(shí)，電離層變化值不會太大，而多路徑殘差以及觀測噪聲等影響同樣變化不大，理論上通過歷元間作差可以消除絕大部分的值，假設(shè)觀測值中含有周跳，便可以通過殘差序列反映出來。此時(shí)，周跳檢驗(yàn)量D表示為：

式中：t表示對應(yīng)的檢驗(yàn)時(shí)間；Δt為數(shù)據(jù)的采樣間隔。

如果信號采樣間隔較短，因而電離層殘差、多路徑效應(yīng)在時(shí)間上相關(guān)性較大，歷元間求差后，其對周跳檢驗(yàn)量的影響較小，可以忽略。由于一般民用用戶不能獲得雙頻P碼觀測值，通常是可獲得L1頻率的C/A碼偽距和L2頻率的P2碼數(shù)據(jù)，這時(shí)可將P2碼偽距代替C/A碼偽距計(jì)算檢驗(yàn)量D，以便獲得相當(dāng)精度的D檢驗(yàn)量。

此時(shí)，周跳檢驗(yàn)量D的精度主要受偽距測量誤差的影響。偽距測量的測距精度對于P碼約為29 cm；而載波相位相應(yīng)的波長為λ1=19.03 cm，λ2=24.42 cm。對 式(4)運(yùn) 用 誤 差傳播定律得：①對于L1：mD≈1.8(周)；②對于L2：mD≈2.2(周)。

因此在能獲得L2的P碼觀測值的情形下，以三倍檢驗(yàn)量方差為限，D檢驗(yàn)量能探測出大于7～8周的周跳。

2 粗差檢驗(yàn)與整體平差

對于上文構(gòu)成的周跳探測檢驗(yàn)量D，其基本服從正態(tài)分布。對其中出現(xiàn)的異常值，通常作以下處理：

令Li=λ?(ti)-P(ti),i=1, 2, …,n。n為歷元個(gè)數(shù)，則相應(yīng)的檢驗(yàn)量殘差陣可表示為：

根據(jù)式(7)可知，系數(shù)矩陣B為(n-1)×n的矩陣，無法根據(jù)最小二乘方法獲得方程的唯一解或最優(yōu)解。

2.1 粗差檢驗(yàn)

對式(7)進(jìn)行粗差定位，結(jié)合檢驗(yàn)量構(gòu)成的誤差方程，利用粗差探測方法，對發(fā)生周跳的歷元進(jìn)行定位，具體方案流程如下。

通 過式(7)中的常 數(shù)項(xiàng)l可知，E(l)=E(L)-L0。觀測值的數(shù)學(xué)期望與近似觀測值的差是一個(gè)很小的數(shù)，一般規(guī)定小于三倍的中誤差，即為觀測值的閉合差。通過對常數(shù)項(xiàng)l進(jìn)行檢驗(yàn)可以確定發(fā)生周跳的具體位置，具體做法為：①首先結(jié)合誤差方程計(jì)算出常數(shù)項(xiàng)l，給出每個(gè)常數(shù)項(xiàng)l的具體初值。②根據(jù)誤差傳播定律，求出l的中誤差，其中L0的中誤差初值根據(jù)先驗(yàn)經(jīng)驗(yàn)值給出，以三倍中誤差m為限，給出具體常數(shù)項(xiàng)的分類范圍，分為大于m和小于m兩個(gè)數(shù)據(jù)集合。③用小于m的li并集的補(bǔ)和大于m的li的集合分別求交集，假設(shè)結(jié)果是單個(gè)元素，那么可以確定該元素值為粗差，此刻標(biāo)記上粗差的具體位置i。④假設(shè)集合中出現(xiàn)兩個(gè)或兩個(gè)以上的元素，那么需要重新建立誤差方程式進(jìn)行甄別，循環(huán)求解。

2.2 整體平差修復(fù)周跳

此時(shí)，假定在觀測的歷元i處發(fā)生了周跳，那么從1～(i-1)歷元處的觀測方程系數(shù)值均為0，此時(shí)，方程呈現(xiàn)出嚴(yán)重病態(tài)，可以不予考慮，即無周跳的歷元不予處理。在i～n之間的所有歷元均受到周跳的影響，此時(shí)需要進(jìn)行處理得到發(fā)生周跳值的大小。B矩陣的劃分如下：

結(jié)合式(7)和式(8)可得，最小二乘表達(dá)式轉(zhuǎn)化為：

此時(shí)，根據(jù)最小二乘原理即可以求得周跳改正的最優(yōu)解，循環(huán)計(jì)算直至i=n，此時(shí)便可修復(fù)觀測值中所有的周跳。

根據(jù)上述方法便可對任意觀測文件進(jìn)行周跳的探測與修復(fù)工作，該方法簡單易行，易于編程實(shí)現(xiàn)。

3 實(shí)例分析

為驗(yàn)證上述方法的可靠性與實(shí)用性，選取2019年7月1日上海佘山站當(dāng)天的觀測數(shù)據(jù)。截取2019年7月1日1時(shí)06分0秒～2019年7月1日2時(shí)13分0秒的22號衛(wèi)星觀測數(shù)據(jù)，選取觀測數(shù)據(jù)，根據(jù)式(4)組成檢驗(yàn)量，其殘差陣序列如圖 1所示。

圖 1 無周跳下檢驗(yàn)量D殘差序列

圖 1中，由偽距/載波構(gòu)成的殘差期望值分布在±0.5周左右，故而取三倍中誤差作為檢驗(yàn)范圍，閾值線如圖1所示。故可以判定該方法在不加入任何其他的改正手段后，其周跳探測范圍為±2周左右。

為驗(yàn)證本文方法對周跳探測的效果，人為在觀測值中加入周跳。加入方式：在10歷元處加入1周周跳；40歷元處加入5周周跳；80歷元處加入8周周跳，加入周跳后的檢驗(yàn)量殘差序列圖如圖2所示。

圖 2 加入周跳后檢驗(yàn)量殘差分布圖

圖 2中，加在10歷元處的1周的周跳無法通過閾值線反映出來，周跳的值與噪聲值相當(dāng)，淹沒在噪聲中，故而會出現(xiàn)探測的誤差；而加在40歷元處的5周的周跳與80歷元處的8周的周跳在圖中很明顯可以辨別出來，周跳的值在相應(yīng)的閾值線之上，此時(shí)通過閾值可以判別出具體哪個(gè)歷元發(fā)生了周跳。

在本文方法中，根據(jù)殘差陣進(jìn)行粗差探測，探測出殘差陣中異常值發(fā)生的具體位置，探測流程按照2.1節(jié)中的步驟進(jìn)行。通過粗差探測，可以得到只有在i=40和i=80時(shí)，用于求交集的兩個(gè)誤差方程常數(shù)項(xiàng)集合的交集才為單個(gè)元素，表示在該地方存在粗差；而i=1處的粗差因?yàn)閭尉?載波組合法的周跳探測能力限制而無法進(jìn)行粗差辨別，出現(xiàn)探測 盲區(qū)。

當(dāng)i=40時(shí)，進(jìn)行第一次循環(huán)平差計(jì)算時(shí)，結(jié)果見表1所列。

表1 第一次循環(huán)平差周跳修復(fù)結(jié)果

由表1可以看出，在40歷元處，經(jīng)過第一次整體平差，修正周跳后的載波觀測值與真值相比相差了1周的誤差，這是因?yàn)榧釉?0歷元處的1周周跳沒有正確修復(fù)的緣故，故而存在 1周的修復(fù)誤差，后面的歷元類似；在80歷元處修正后的值與真值差值為9周，其原因是人為加在上面的8周周跳以及1周周跳粗差。第一次平差計(jì)算結(jié)束后，i指向80歷元處，進(jìn)行第二次平差計(jì)算。

當(dāng)i=80時(shí)，進(jìn)行第二次循環(huán)平差計(jì)算時(shí)，結(jié)果見表2所列。

表2 第二次循環(huán)平差周跳修復(fù)結(jié)果

由表2的平差結(jié)果可以看出，經(jīng)過第二次平差后，加入周跳后觀測值的周跳值得到了修復(fù)，40和80歷元處的較大周跳得到了正確修復(fù)，而加在10歷元處的周跳因?yàn)樘綔y方法的局限性，無法得到正確辨別，因而存在一定的探測盲區(qū)。

4 結(jié)論

本文對經(jīng)典方法偽距/載波組合法的原理進(jìn)行了詳細(xì)闡述，并結(jié)合粗差探測的方法對周跳進(jìn)行探測修復(fù)。從原理上來說，該方法可行而且易于實(shí)現(xiàn)編程。結(jié)合實(shí)例，驗(yàn)證了本文方法的優(yōu)缺點(diǎn)。

通過仿真算例實(shí)例驗(yàn)證得出，本文方法對于2周以上的周跳具有較強(qiáng)的探測修復(fù)能力，且探測速度較快，不需要復(fù)雜的判斷迭代過程，探測準(zhǔn)確，修復(fù)精度高；但是對于一些小周跳，由于初始方法的周跳探測能力限制，使其對小周跳的探測無能為力。因此下一步的工作主要是如何利用更為精確的周跳探測方法結(jié)合更加優(yōu)化的粗差探測方法，同時(shí)可以考慮多頻觀測條件下的周跳探測方法，進(jìn)一步優(yōu)化本文提出的方法，提升該方法的周跳探測適用 范圍。