基于值域的多約束GNSS單頻單歷元定姿新算法

吳洪濤，趙修斌，龐春雷，張良，馮波

空軍工程大學(xué) 信息與導(dǎo)航學(xué)院，西安 710077

單頻單歷元姿態(tài)測量由于所具有實(shí)時(shí)性高、低成本、不受整周跳變影響的優(yōu)勢,已成為當(dāng)前航空工程領(lǐng)域的研究熱點(diǎn),但單頻單歷元解算模型強(qiáng)度太弱，整周模糊度的難度較大。利用多天線可以顯著增強(qiáng)單頻單歷元解算模型的強(qiáng)度[1]，實(shí)現(xiàn)載體航向角、俯仰角和橫滾角的全姿態(tài)測量[2]。Giorgi和Teunissen利用坐標(biāo)旋轉(zhuǎn)矩陣建立多天線定姿模型，將LAMBDA算法擴(kuò)展為多基線多約束LAMBDA (MC-LAMBDA)算法[3-5]。龔昂將利用解析方法求解最優(yōu)旋轉(zhuǎn)矩陣，進(jìn)一步簡化了MC-LAMBDA算法，提高了算法效率[6]。李青松將坐標(biāo)旋轉(zhuǎn)矩陣的單位正交約束條件融入了模糊度浮點(diǎn)解和協(xié)方差矩陣的求解過程，提高了單歷元固定的成功率和正確率[7]。這些方法都是基于LAMBDA算法從模糊度域確定搜索空間，解算效率隨模糊度維數(shù)增加快速下降，計(jì)算量隨天線數(shù)量增加顯著增長。

除了上述算法外，模糊度解算算法還有部分模糊度法[8-11]，長短基線法[12-13]、模糊度函數(shù)法[14-15]等。部分模糊度法基于分組的思想，將模糊度分為基本組和剩余組，直接對(duì)基本組模糊度進(jìn)行求解，但分步求解降低了理論成功率，且一般很難解析給出模糊度分組算法。長短基線法利用無模糊度短基線輔助長基線解算模糊度，計(jì)算效率高，但成倍地增加了硬件成本。模糊度函數(shù)法通過將多維模糊度的搜索轉(zhuǎn)化為三維坐標(biāo)域的搜索，計(jì)算效率與模糊度維數(shù)無關(guān)，但在多天線條件下建立坐標(biāo)域搜索涉及到多基線之間的相對(duì)位置轉(zhuǎn)換，實(shí)現(xiàn)較為復(fù)雜。

文獻(xiàn)[16]提出的基于值域的模糊度搜索方法，在基線長度已知的條件下，將對(duì)模糊度域的搜索轉(zhuǎn)化為對(duì)姿態(tài)域的搜索，方法簡單，易于實(shí)現(xiàn)。首先，對(duì)多天線姿態(tài)測量算法進(jìn)行了分析；然后，將姿態(tài)約束融入解算模型中，建立了約束條件下值域搜索模型，提高了模糊度單歷元搜索效率，推導(dǎo)了基于最優(yōu)姿態(tài)解近似估計(jì)的模糊度固定算法。該算法較標(biāo)準(zhǔn)迭代算法計(jì)算效率大大增加，且具有與標(biāo)準(zhǔn)算法相當(dāng)?shù)男阅?，能夠?qū)崿F(xiàn)模糊度的單頻單歷元可靠固定。

1 多天線定姿算法及分析

1.1 多天線定姿模型

基線矢量在載體坐標(biāo)系和東北天(ENU)坐標(biāo)系之間的轉(zhuǎn)換關(guān)系為

a′=Ra

(1)

式中：a為載體坐標(biāo)系下的基線矢量；a′為ENU坐標(biāo)系下的基線矢量；R為由載體坐標(biāo)系到ENU坐標(biāo)系的坐標(biāo)旋轉(zhuǎn)矩陣。R可由姿態(tài)角表示為

R=

(2)

式中：Cξ=cosξ；Sξ=sinξ；θ、β和γ分別為載體的航向角、俯仰角和橫滾角；θ∈[0,2π]、β∈[-π,π]、γ∈[-π,π]。設(shè)

(3)

R中的元素滿足單位正交條件,即

(4)

理論上只要得到R前兩行元素,即可解得全部姿態(tài)角。

(5)

β=-arcsinr31

(6)

(7)

當(dāng)存在多條基線時(shí)，式(1)可寫為矩陣形式，n+1個(gè)天線對(duì)應(yīng)n條相互獨(dú)立的基線。設(shè)天線平面與載體坐標(biāo)系重合，則

(8)

A′=RA

(9)

式中:A為載體坐標(biāo)系下的基線矢量矩陣；A′為ENU坐標(biāo)系下的基線矢量矩陣。多天線雙差載波相位觀測方程可寫為[17]

[Φ1Φ2…Φn]=

(10)

式中：G為雙差視線向量矩陣；Φ1、Φ2、…、Φn分別為基線1、基線2、…、基線n的雙差載波相位觀測向量；N1、N2、…、Nn分別為基線對(duì)應(yīng)的雙差整周模糊度向量；λ為載波波長。

對(duì)式(10)進(jìn)行向量化運(yùn)算可得

(11)

式中：?為Kronecker積；vec(·)為向量化運(yùn)算。同理可得到多天線雙差偽距觀測方程為

(12)

式中：P1、P2、…、Pn分別為基線1、基線2、…、基線n的雙差偽距觀測向量。

令

則式(11)和式(12)有[17]

(13)

P=BM

(14)

當(dāng)?shù)玫秸_的模糊度空間之后，利用最小二乘原理固定模糊度，代價(jià)函數(shù)為

(15)

式中：QΦ和QP分別為雙差載波相位觀測量和雙差偽距觀測量的協(xié)方差矩陣。

當(dāng)成功固定整周模糊度N之后，可根據(jù)式(13)解出M,由式(5)～式(7)可得到待求解的姿態(tài)角。

1.2 算法分析

多天線增強(qiáng)了單頻單歷元解算模型的強(qiáng)度，但模糊度維數(shù)隨天線數(shù)量成倍增加，從模糊度域確定搜索空間的計(jì)算量將顯著增長?；谥涤虻哪：人阉魉惴ǎ瑢?duì)模糊度域的搜索轉(zhuǎn)化為對(duì)姿態(tài)域的搜索，搜索過程與模糊度維數(shù)無關(guān)，適用于多天線姿態(tài)測量。但該算法沒有將實(shí)際應(yīng)用中的姿態(tài)約束信息融入解算模型，單歷元模糊度搜索效率不高。一般通過對(duì)低成本慣性測量單元(Inertial Measurement Unit，IMU)輸出數(shù)據(jù)進(jìn)行解算，或水平儀、微機(jī)電系統(tǒng)傳感器(Micro Electro Mechanical Systems，MEMS)等其他低成本傳感器可以得到俯仰角、橫滾角的強(qiáng)約束[18-19]。

當(dāng)?shù)玫秸_的模糊度空間之后，利用總體最小二乘原理固定模糊度?？紤]俯仰角、橫滾角的約束信息，代價(jià)函數(shù)可寫為[20]

(16)

利用純載波相位方程求解N條件下M的最小二乘解，往往不能滿足單位正交約束條件。理論上需要在單位正交約束條件下，迭代求解N條件下坐標(biāo)旋轉(zhuǎn)矩陣M的最優(yōu)解，但迭代求解會(huì)大大增加計(jì)算量。為了減小計(jì)算量，可以通過求解M最優(yōu)解的近似解[7]，固定模糊度。近似解一方面需要具有與迭代解相當(dāng)?shù)男阅?，另一方面需要盡量減小計(jì)算量。

2 本文所提算法

2.1 多約束多天線模糊度搜索模型的建立

由式(2)M可以由θ、β、γ表示，即

M=M(θ,β,γ)

代入式(13)整理得

(17)

式中：下標(biāo)i為第i條基線;j為第i條基線對(duì)應(yīng)的第j個(gè)雙差觀測方程。Φij的誤差通常在1/4周以內(nèi)，若已知θ、β、γ,通過直接四舍五入取整即可確定整周模糊度:

(18)

式中：round(·)為四舍五入取整。整周模糊度可以看做姿態(tài)角的函數(shù)?；诖?，通過直接對(duì)姿態(tài)角進(jìn)行遍歷搜索，解算出相應(yīng)的模糊度浮點(diǎn)解，取整即可得到模糊度的搜索空間。

為了提高搜索效率，合理選擇搜索步長成為關(guān)鍵。為了確保整周模糊度不被漏搜，并考慮到由于取整時(shí)四舍五入導(dǎo)致的模糊度變化，要求搜索步長使模糊度的每次變化不能超過0.5周。以航向角為例，設(shè)Δθ為航向角的搜索步長，即

(19)

當(dāng)Δθ小于10°有sin(Δθ)≈Δθ,cos(Δθ)≈1代入式(19)整理得

i=1,2,…,n;j=1,2,…,t

(20)

式中：n為相互獨(dú)立的基線數(shù)；t為雙差觀測方程個(gè)數(shù)。

文獻(xiàn)[16]通過式(20)確定搜索步長時(shí)不考慮俯仰角、橫滾角約束，而將約束信息融入式(20)中能進(jìn)一步增大搜索步長，提高搜索效率。

2.2 姿態(tài)角約束條件下搜索步長的確定

設(shè)

(21)

(22)

(23)

由范數(shù)理論

(24)

設(shè)

則航向角搜索步長Δθ為

(25)

同理可得俯仰角、橫滾角的搜索步長Δβ、Δγ為

(26)

(27)

其中類似式(24)推導(dǎo)過程

(28)

(29)

(30)

(31)

(32)

(33)

2.3 姿態(tài)解近似估計(jì)的模糊度固定

利用歐拉角對(duì)式(16)重新參數(shù)化后,代價(jià)函數(shù)為

(34)

式中：

L=F(N,η)

(35)

(36)

利用最小二乘原理有

Δη≈(KT(η0)QK(η0))-1KT(η0)QΔL

(37)

式中：ΔL=L-F(N,η0)。

(38)

所提算法流程圖如圖1所示。

圖1 所提算法流程圖Fig.1 Flow chart of proposed algorithm

3 實(shí)驗(yàn)結(jié)果及分析

為驗(yàn)證算法效果，分別進(jìn)行了靜態(tài)驗(yàn)證實(shí)驗(yàn)和動(dòng)態(tài)實(shí)驗(yàn)。

3.1 靜態(tài)實(shí)驗(yàn)條件及結(jié)果分析

分別進(jìn)行了兩組靜態(tài)實(shí)驗(yàn)，利用GPS-703-GGG 天線和NovAtel OEM628板卡采集GPS系統(tǒng)L1頻點(diǎn)數(shù)據(jù)，采樣頻率為1 Hz，實(shí)驗(yàn)地點(diǎn)為空軍工程大學(xué)信息與導(dǎo)航學(xué)院科研樓頂，數(shù)據(jù)采用MATLAB2014A進(jìn)行事后處理，電腦硬件配置為Intel(R) Core(TM) i3-2350M雙核CPU，主頻2.3 GHz，內(nèi)存4 G。

第1組靜態(tài)實(shí)驗(yàn)天線配置如圖2所示，共采集63 488個(gè)歷元的數(shù)據(jù)。基線在載體系下的基線矢量矩陣為

(39)

根據(jù)事先測定，載體在ENU坐標(biāo)系中的航向角為277.424°，俯仰角為2.459°，橫滾角為1.814°。實(shí)驗(yàn)中，給俯仰角、橫滾角加入均方差為1°的零均值高斯噪聲，模擬由其他傳感器得到的俯仰角、橫滾角信息。俯仰角、橫滾角的搜索范圍為±5°，航向角的搜索范圍為±180°。

圖2 第1組靜態(tài)實(shí)驗(yàn)天線配置Fig.2 Configuration of antennas of first static experiment

第2組靜態(tài)實(shí)驗(yàn)天線配置如圖3所示，共采集1 716個(gè)歷元的數(shù)據(jù)，基線在載體系下的基線矢量矩陣為

圖3 第2組靜態(tài)實(shí)驗(yàn)天線配置Fig.3 Configuration of antennas of second static experiment

(40)

實(shí)驗(yàn)中，利用俯仰角、橫角滾先驗(yàn)信息由NovAtel SPAN-CPT慣性導(dǎo)航設(shè)備粗對(duì)準(zhǔn)獲得(標(biāo)準(zhǔn)差設(shè)為0.6°)。俯仰角、橫滾角的搜索范圍為±3°，航向角的搜索范圍仍為±180°

3.1.1 模糊度搜索算法比較

為了對(duì)本文模糊度搜索算法與原算法進(jìn)行比較，采用以下兩種算法對(duì)數(shù)據(jù)進(jìn)行解算，單歷元確定模糊度搜索空間。

算法1原算法，即不考慮姿態(tài)角約束推導(dǎo)搜索步長確定模糊度搜索空間。

算法2本文所提算法，即利用姿態(tài)角約束推導(dǎo)步長確定模糊度搜索空間。

利用兩種算法對(duì)兩組靜態(tài)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行解算，逐歷元搜索步長和模糊度空間大小如圖4所示，圖4(a)為第1組實(shí)驗(yàn)結(jié)果，圖4(b)為第2組實(shí)驗(yàn)結(jié)果。虛線代表算法1，實(shí)線代表算法2，平均搜索步長和模糊度空間平均大小如表1所示。由于第2組實(shí)驗(yàn)基線長度大于第1組實(shí)驗(yàn)，因此第2組實(shí)驗(yàn)姿態(tài)角搜索步長明顯減小，模糊度空間也顯著增加。但兩組實(shí)驗(yàn)結(jié)果中算法2均優(yōu)于算法1。

第1組實(shí)驗(yàn)中算法2與算法1相比，航向角平均搜索步長增大了34.5%，俯仰角平均搜索步長增大了150.6%，橫滾角平均搜索步長增大了126.9%，模糊度空間減小了65.8%。

第2組實(shí)驗(yàn)中算法2與算法1相比，航向角平均搜索步長增大了41.3%，俯仰角平均搜索步長增大了59.1%，橫滾角平均搜索步長增大了293.4%，模糊度空間減小了67.7%。

實(shí)驗(yàn)結(jié)果表明，本文所提算法充分利用了姿態(tài)角約束條件，單歷元模糊度搜索效率明顯提高。

圖4 兩種模糊度搜索算法實(shí)驗(yàn)結(jié)果Fig.4 Experimental results of two ambiguity search algorithms

表1 兩種模糊度搜索算法實(shí)驗(yàn)結(jié)果對(duì)比Table 1 Comparison of experimental results of two ambiguity search algorithms

3.1.2 模糊度固定算法比較

根據(jù)事先精確測量的姿態(tài)角，逐歷元反解得到正確的模糊度。在利用約束條件下的值域搜索模型得到正確的模糊度空間，采用以下4種方案對(duì)兩組靜態(tài)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行實(shí)驗(yàn)。

方案3本文所提算法,即基于最優(yōu)條件姿態(tài)解近似估計(jì)，計(jì)算代價(jià)函數(shù)式(17)。

方案4標(biāo)準(zhǔn)迭代算法,即利用牛頓迭代法求解最優(yōu)姿態(tài)解，計(jì)算代價(jià)函數(shù)式(17)。

4種方案計(jì)算的正確模糊度對(duì)應(yīng)的代價(jià)函數(shù)值如圖5所示，圖5(a)為第1組前500歷元實(shí)驗(yàn)結(jié)果,圖5(b)為第2組前500歷元實(shí)驗(yàn)結(jié)果。對(duì)4種方案計(jì)算的正確模糊度的平均代價(jià)函數(shù)值和固定模糊度的平均耗時(shí)進(jìn)行統(tǒng)計(jì)，結(jié)果如表2所示。

圖5 4種方案計(jì)算的前500歷元正確模糊度對(duì)應(yīng)的代價(jià)函數(shù)值Fig.5 Cost function value calculated from correct ambiguity for the first 500 epoch data by four schemes

表24種方案計(jì)算的正確模糊度平均代價(jià)函數(shù)值和固定模糊度的平均耗時(shí)
Table2Meanofcostfunctionvaluecalculatedfromcorrectambiguityandcomputation-timetofixambiguityoffourschemes

實(shí)驗(yàn)結(jié)果方案1方案2方案3方案4第1組實(shí)驗(yàn)正確模糊度的代價(jià)函數(shù)值0．0532．0170．2090．196固定模糊度平均耗時(shí)/s0．0070．0060．0450．226第2組實(shí)驗(yàn)正確模糊度的代價(jià)函數(shù)值0．0914．9332．6532．078固定模糊度平均耗時(shí)/s0．0470．0640．3141．301

方案1不考慮角度約束信息(式(16)中不包含角度約束信息加權(quán)部分)，因此得到正確模糊度的代價(jià)函數(shù)值最小。方案2沒有充分利用角度約束信息，計(jì)算得到正確模糊度的代價(jià)函數(shù)值遠(yuǎn)大于方案4。方案3充分考慮了角度約束信息，且利用最優(yōu)條件姿態(tài)角近似估計(jì)簡化了迭代過程。因此，方案3與方案4相比，正確模糊度的代價(jià)函數(shù)值幾乎一致，且計(jì)算時(shí)間大大減少。第2組實(shí)驗(yàn)中，方案4的模糊度固定時(shí)間大于1 s，幾乎不能實(shí)時(shí)解算，而方案3的計(jì)算時(shí)間遠(yuǎn)遠(yuǎn)小于1 s。

4種方案單歷元固定模糊度正確率，如表3所示。方案1僅考慮的基線約束，沒有考慮角度約束信息，固定模糊度正確率最低。方案2計(jì)算最優(yōu)條件解時(shí)利用角度約束條件不充分，固定模糊度正確率低于方案3。方案3與方案4由于充分利用了約束條件，固定模糊度正確率最高。

利用第1組實(shí)驗(yàn)數(shù)據(jù)，對(duì)4種方案在不同可視衛(wèi)星條件下單歷元固定模糊度的正確率進(jìn)行了統(tǒng)計(jì)，如表4所示。

方案1、方案2固定模糊度的正確率隨著可視衛(wèi)星數(shù)量的減少而降低，方案3由于充分利用了約束條件，固定模糊度的正確率不受影響，與方案4相同。

實(shí)驗(yàn)結(jié)果表明，利用俯仰角、橫滾角約束條件能夠進(jìn)一步提高模糊度固定的成功率和可靠性。本文所提算法相比直接法犧牲了一定的效率，但卻有與標(biāo)準(zhǔn)算法相當(dāng)?shù)牧己眯阅?，且?jì)算效率遠(yuǎn)遠(yuǎn)高于標(biāo)準(zhǔn)算法。

表3 4種方案固定模糊度正確率Table 3 Fixing ambiguity success rate of four schemes

表4不同可視衛(wèi)星數(shù)量下4種方案固定模糊度正確率

Table4Fixingambiguitysuccessrateoffourschemesonconditionofdifferentvisiblesatellites

3.2 動(dòng)態(tài)實(shí)驗(yàn)條件及結(jié)果分析

動(dòng)態(tài)實(shí)驗(yàn)時(shí)，將3個(gè)天線固定在實(shí)驗(yàn)車頂，采集GPS系統(tǒng)L1頻點(diǎn)的數(shù)據(jù)，采樣頻率為1 Hz。將NovAtel SPAN-CPT慣性導(dǎo)航設(shè)備水平固定在實(shí)驗(yàn)車的中軸線上，以初始化后的組合輸出作為姿態(tài)角參考值(測姿精度俯仰角為0.02°，橫滾角為0.02°，航向角為0.06°)，為了與靜態(tài)實(shí)驗(yàn)保持一致。實(shí)驗(yàn)地點(diǎn)為西安市東郊，實(shí)驗(yàn)場景如圖6所示，共采集896個(gè)歷元的數(shù)據(jù)，其中前347 s實(shí)驗(yàn)車處于靜止?fàn)顟B(tài)，用于標(biāo)定3個(gè)天線的姿態(tài)測量系統(tǒng)與組合導(dǎo)航系統(tǒng)的安裝誤差?；€在載體坐標(biāo)系下的基線矢量矩陣為

(41)

利用方案3對(duì)數(shù)據(jù)進(jìn)行解算，為了與靜態(tài)實(shí)驗(yàn)一致，俯仰角、橫滾角先驗(yàn)信息的標(biāo)準(zhǔn)差仍設(shè)為0.6°，搜索范圍為±3°，航向角的搜索范圍為±180°。

采用方案3和方案4對(duì)數(shù)據(jù)進(jìn)行單歷元解算，確定姿態(tài)信息。

方案4不考慮姿態(tài)角約束推導(dǎo)搜索步長確定模糊度搜索空間，利用牛頓迭代法計(jì)算代價(jià)函數(shù)式(17)，單歷元固定模糊度。

方案3考慮姿態(tài)角約束推導(dǎo)搜索步長確定模糊度搜索空間，利用最優(yōu)條件姿態(tài)解近似估計(jì)計(jì)算代價(jià)函數(shù)式(17)，單歷元固定模糊度。

圖6 動(dòng)態(tài)實(shí)驗(yàn)天線配置Fig.6 Configuration of antennas of dynamic experiment

兩種方案得到模糊度空間大小和平均耗時(shí)如表5所示。方案4幾乎不能實(shí)時(shí)輸出。與方案4相比，方案3具有更小的搜索模糊度空間和更高的計(jì)算效率，能夠做到實(shí)時(shí)輸出。

以姿態(tài)角參考值作為理論真值，將解算得到俯仰角、方位角與真值對(duì)比，得到誤差曲線如圖7所示，姿態(tài)角均方誤差(RMSE)如表6所示。

航向角的RMSE為0.170 2°，俯仰角的RMSE為0.390 2°，橫滾角的RMSE為0.368 9°。

實(shí)驗(yàn)結(jié)果表明，利用本文所提算法解算的姿態(tài)結(jié)果外符合精度良好，姿態(tài)信息測量準(zhǔn)確，模糊度得到了正確的固定。

表5 兩種方案計(jì)算效率對(duì)比Table 5 Efficiency of two schemes

圖7 動(dòng)態(tài)實(shí)驗(yàn)測量結(jié)果Fig.7 Attitude results of dynamic experiment

表6 航向角、俯仰角和橫滾角的 RMSETable 6 RMSE of heading, pitch and roll angles

4 結(jié) 論

1) 融入姿態(tài)約束的值域搜索模型相對(duì)不考慮姿態(tài)角約束，能夠進(jìn)一步增大搜索步長，提高搜索效率，減小模糊度空間。

2) 利用最優(yōu)姿態(tài)解近似估計(jì)固定模糊度，具有與標(biāo)準(zhǔn)迭代算法相當(dāng)?shù)恼_率，且計(jì)算效率遠(yuǎn)遠(yuǎn)高于標(biāo)準(zhǔn)迭代算法。

雖然利用俯仰角、橫滾角約束能夠提高值域搜索算法的效率，但當(dāng)基線較長時(shí)，利用本文算法確定的模糊度空間仍較大，計(jì)算效率不高。從靜態(tài)實(shí)驗(yàn)結(jié)果來看，在基線長度達(dá)到3.678 m時(shí)(第2組靜態(tài)實(shí)驗(yàn))，本文算法的實(shí)時(shí)性明顯下降。因此使用本文算法進(jìn)行單頻單歷元姿態(tài)測量時(shí)，基線長度一般不應(yīng)大于4 m。下一步重點(diǎn)考慮如何有效利用航向角約束和偽距觀測信息，進(jìn)一步提高值域類算法在較長基線條件下的適用性。