基于MATLAB的GPS周跳探測與修復(fù)

曹小鴻

(三和數(shù)碼測繪地理信息技術(shù)有限公司，甘肅 天水 741000)

GPS導(dǎo)航定位可以通過偽距和載波相位實(shí)現(xiàn)，而偽距的測量精度低于載波相位的精度，所以在高精度定位中一般用載波相位測量實(shí)現(xiàn)。而載波相位是由整周部分和不足一周的小數(shù)部分構(gòu)成，整數(shù)部分的連續(xù)變化值可以通過多普勒計(jì)數(shù)得到，載波相位只能測量相位中小于一周的值，但如果在測量過程中信號(hào)受到遮擋、信噪比低以及接收機(jī)出現(xiàn)故障等情況都會(huì)導(dǎo)致整周計(jì)數(shù)部分發(fā)生突變，即發(fā)生整周跳變的現(xiàn)象[1-4]，簡稱為周跳(cycle slip)。由于周跳在GPS高精度定位中影響特別大，所以必須進(jìn)行修復(fù)。

1 周跳探測原理

GPS觀測量是一種時(shí)間序列，周跳是GPS載波相位觀測數(shù)據(jù)序列中的一種粗差。探測GPS載波相位觀測量是否存在周跳，就是發(fā)現(xiàn)其是否發(fā)生突變，即通過比較觀測量在相鄰歷元間的變化是否連續(xù)，從而判斷周跳的位置和大小。周跳探測的方法比較多，各有優(yōu)缺點(diǎn)[5-7]。本文主要介紹電離層殘差法和Melbourne-Wübbena組合法。

1.1 電離層殘差法

電離層殘差法由Goad在1986年提出，利用雙頻載波相位測量的電離層殘差進(jìn)行探測和修復(fù)周跳[8]。它主要利用不同歷元電離層誤差的變化。在不考慮多路徑效應(yīng)和其他測量噪聲的條件下，同一歷元的雙頻載波相位測量值差為:

φ(t)=λ1φ1(t)-λ2φ2(t)

(1)

式中，t為觀測時(shí)刻；φ1(t)和φ2(t)為L1和L2兩個(gè)頻率的載波相位；λ1和λ2為L1和L2兩個(gè)頻率對(duì)應(yīng)的波長；N1和N2為L1和L2兩個(gè)頻率的整周未知數(shù)；f1和f2為L1和L2兩個(gè)頻率對(duì)應(yīng)的頻率;A(T)為其他項(xiàng)。將式(1)兩端同除以λ1，則有:

(2)

=Δion(t2)-Δion(t1)+ε

(3)

當(dāng)電離層變化不大時(shí)，檢測量D1的值會(huì)在0附近上下波動(dòng)。若存在周跳，令ΔN1、ΔN2分別為L1、L2的周跳值，如式(4)所示：

(4)

1.2 Melbourne-Wübbena組合法

Melbourne-Wübbena組合法是由Melbourne和Wübbena提出的，簡稱M-W組合法[9-10]。由于M-W組合的波長約為86 cm，因?yàn)椴ㄩL較長對(duì)探測周跳特別有利。M-W組合法是載波相位寬巷組合減去測碼偽距窄巷組合，因此，其組合觀測值如下：

(5)

式中，P1和P2為偽距觀測值；L1和L2為載波相位觀測值；f1和f2為對(duì)應(yīng)的頻率；λ1和λ2為對(duì)應(yīng)的波長。則M-W組合的寬巷整周模糊度為：

(6)

則M-W組合法周跳檢測量D2為：

D2=N(t2)-N(t1)

(7)

1.3 電離層殘差法和M-W組合法組合探測周跳

由于當(dāng)L1和L2周跳比值為77/60時(shí)，電離層殘差法無法探測出周跳，但M-W組合法可以探測出周跳；當(dāng)L1和L2周跳比值為1時(shí)，M-W組合法無法探測出周跳，但電離層殘差法可以探測出周跳。則可以采用電離層殘差法和M-W組合法組合起來探測周跳，聯(lián)立式(4)和式(6)得方程組：

(8)

則周跳值N為：

(9)

至此，GPS周跳發(fā)生的位置和大小可以求得。

2 算法實(shí)現(xiàn)

2.1 數(shù)據(jù)預(yù)處理

本文探測周跳的數(shù)據(jù)選擇國內(nèi)IGS站的SHAO站，由于MATLAB軟件在數(shù)據(jù)處理方面有著特別強(qiáng)大的功能，編寫了讀取O文件的程序，然后，根據(jù)探測周跳的方法完成了程序的編寫，電離層殘差法和M-W組合法組合探測周跳方法核心代碼如下：

f1=1575.42*10^6; %f1的頻率

f2=1227.60*10^6; %f2的頻率

C0=299792458; %光速

b1=C0/f1; %f1對(duì)應(yīng)的波長λ1

b2=C0/f2; %f2對(duì)應(yīng)的波長λ2

b3=C0/(f1-f2); %寬巷組合對(duì)應(yīng)的波長λ3

for k=1:1:t

MW(k,1)=Obs_L1(k)-Obs_L2(k)-(f1*Obs_P1(k)+f2*f1*Obs_P2(k))/(b3*(f1+f2));

end

for j=1:1:(t-1)

GFMW(j,1)=Obs_L1(j+1)-Obs_L1(j)-(f1/f2)*(Obs_L2(j+1)-Obs_L2(j));

GFMW(j,2)=MW(j+1,1)-MW(j,1);

end

A1=[1 -f1/f2;1 -1];

D1=GFMW(1:t-1,1:2);

CycleJump(1:t-1,1:2)=((inv(A1))*D1')';

CycleJumpFix=fix(CycleJump (1:t-1,1:2));

實(shí)驗(yàn)選取了該站10號(hào)GPS衛(wèi)星200個(gè)觀測歷元的連續(xù)觀測數(shù)據(jù)，先對(duì)該200個(gè)歷元的觀測數(shù)據(jù)進(jìn)行周跳探測與修復(fù)，形成一組無周跳的觀測數(shù)據(jù)，然后，根據(jù)各種方法的優(yōu)缺點(diǎn)人為的加入周跳，接著進(jìn)行周跳探測，處理過的周跳圖如圖1所示。

圖1 無周跳的原始數(shù)據(jù)

2.2 周跳探測

根據(jù)不同方法的優(yōu)缺點(diǎn)，人為的對(duì)L1和L2加入不同的周跳值，然后，分別用每種方法進(jìn)行周跳探測，加入的周跳數(shù)如表1所示。

表1 周跳探測與修復(fù)結(jié)果

分析圖2數(shù)據(jù)和表1內(nèi)容，將人為加入的周跳值通過兩種組合方法進(jìn)行探測。因?yàn)樾≈芴綔y比較困難，所以在第20和第50歷元處對(duì)L1和L2分別加入了周跳值為1，探測結(jié)果也比較理想，探測結(jié)果殘差值均較小(小于0.05周)，取整后的結(jié)果完全一致，完全實(shí)現(xiàn)了小周跳的探測。在第80和第100歷元處對(duì)L1和L2都加入較小的周跳值，再次驗(yàn)證組合方法對(duì)較小周跳的適用程度，殘差值為0.01周，取整后的結(jié)果完全一致，結(jié)果表明組合方法對(duì)較小周跳完全適用。由于當(dāng)L1和L2周跳比值為77/60時(shí)，電離層殘差法無法探測出周跳，但M-W組合法可以探測出周跳，所以在第140歷元的時(shí)候，在L1上加入了77周，L2上加入了66周，進(jìn)行組合法探測，探測結(jié)果是正確的。當(dāng)L1和L2周跳比值為1時(shí)，M-W組合法無法探測出周跳，但電離層殘差法可以探測出周跳，在第170歷元處對(duì)L1和L2都加入50周，探測結(jié)果殘差值均小于0.05周，取整后的結(jié)果完全一致，周跳探測結(jié)果理想，證明組合法既對(duì)大周跳探測準(zhǔn)確，又對(duì)L1和L2周跳比值為1的時(shí)候周跳探測是有效的。本文研究重點(diǎn)為通過人為模擬實(shí)驗(yàn)進(jìn)行周跳探測，該方法對(duì)定位結(jié)果影響程度分析將在后續(xù)具體的實(shí)測數(shù)據(jù)中進(jìn)行分析，本文將不再進(jìn)行討論。

圖2 周跳探測結(jié)果

2.3 周跳修復(fù)

根據(jù)探測出周跳值和所在的位置，計(jì)算發(fā)生周跳的歷元處的累計(jì)周跳值，然后，對(duì)兩個(gè)頻率的載波相位觀測值進(jìn)行修復(fù)，修復(fù)結(jié)果如圖3所示。根據(jù)修復(fù)結(jié)果，再進(jìn)行周跳探測，分析表明，未探測出明顯周跳，說明修復(fù)效果較為理想，基本修復(fù)了人為加入的周跳。

圖3 修復(fù)完成后的周跳探測圖

3 結(jié) 論

由于周跳在GPS精密定位中對(duì)定位精度影響極大，因此，在數(shù)據(jù)處理前必須先探測和修復(fù)周跳，由電離層殘差法和M-W組合法組合成的方法可以避開這兩種方法的缺點(diǎn)，利用其優(yōu)點(diǎn)相互補(bǔ)充，可以實(shí)現(xiàn)載波相位發(fā)生的各種大小和任意比例周跳的探測，從而提高GPS精密定位的定位精度。