旋轉(zhuǎn)長基線干涉儀“維特比-卡爾曼”濾波定位算法

馬 爽，蔡 凱，劉元華

（北京遙感信息研究所，北京 100011）

0 引言

旋轉(zhuǎn)長基線干涉儀（RLBI）通過不斷旋轉(zhuǎn)形成時序上的多條虛擬基線來完成相位差解模糊和定位。與傳統(tǒng)的多通道干涉儀和陣列測向等方法相比，旋轉(zhuǎn)基線體制具有設(shè)備簡單、通道一致性要求低等優(yōu)點(diǎn)，但是，長基線干涉儀只能測量得到模糊相位差，如何利用模糊相位差高效準(zhǔn)確地定位，是旋轉(zhuǎn)干涉儀定位體制需要解決的首要問題。

要解決該問題，最直接的思路是首先將模糊相位差恢復(fù)成真實(shí)相位差，再進(jìn)行定位解算。文獻(xiàn)[1-4]分別研究了圓陣、GPS測向、導(dǎo)彈導(dǎo)引頭測向和近場測向等條件下相位差解模糊的方法，實(shí)現(xiàn)了針對具體應(yīng)用場景的定位解算。但是，該類方法對信號測量間隔和精度要求高，一般需要滿足相對無模糊條件才能實(shí)現(xiàn)。為了擴(kuò)展適用范圍，文獻(xiàn)[5-7]探索了利用模糊相位差直接定位的方法，在定位解空間利用各種優(yōu)化算法直接搜索最優(yōu)解。該類方法不再需要恢復(fù)真實(shí)相位差，但是，計算量隨著可能定位區(qū)域的擴(kuò)大迅速增加，且定位精度與搜索粒度密切相關(guān)，在具體應(yīng)用中受到一定限制。

本文在研究旋轉(zhuǎn)長基線干涉儀相位差和模糊數(shù)變化特點(diǎn)的基礎(chǔ)上，提出了相位差“條帶”和“路徑”的概念，將需要搜索的解空間限制在模糊數(shù)序列空間中，與直接在定位空間中搜索相比，大大降低了解空間的維度。用隱馬爾可夫模型（HMM）對“條帶”轉(zhuǎn)移過程和相位差觀測過程進(jìn)行建模，通過結(jié)合“卡爾曼”濾波過程和“維特比”狀態(tài)解碼過程，同步實(shí)現(xiàn)了相位差解模糊和濾波定位。該算法能夠?qū)?shù)據(jù)進(jìn)行序貫處理，實(shí)時完成解模糊和定位，

定位性能與使用真實(shí)相位差進(jìn)行定位時的“卡爾曼”濾波算法一致。

1 旋轉(zhuǎn)長基線干涉儀相位差觀測過程HMM模型

旋轉(zhuǎn)基線通過不斷旋轉(zhuǎn)來獲取目標(biāo)在不同基線轉(zhuǎn)角下的相位差測量值。在旋轉(zhuǎn)過程中，目標(biāo)與基線的相對位置不斷變化，對應(yīng)的模糊數(shù)和相位差測量值也隨之改變，這一過程可以用HMM[8]來建模。

1.1 模糊相位差、“條帶”和“路徑”

旋轉(zhuǎn)基線干涉儀采用單根基線測量信號相位差，真實(shí)相位差為：

式中，d為基線長度，λ是信號波長，θ為信號到達(dá)方向與基線法平面的夾角。

干涉儀輸出的測量相位差φ′∈(-π,π]，當(dāng)d＞λ時，存在相位差模糊，測量相位差φ′與真實(shí)相位差φ之間滿足如下關(guān)系：

式中，mod2π(·)表示按照2π取模，整數(shù)n稱為模糊數(shù)。

模糊數(shù)n的取值范圍為[-N,N]，N為不超過dλ的整數(shù)。例如，當(dāng)信號頻率為3 GHz、基線長度為1 m時，模糊數(shù)的取值范圍為[-10,10]。

假設(shè)干涉儀基線在定位區(qū)域上方水平放置，定位區(qū)域內(nèi)每個點(diǎn)對應(yīng)的真實(shí)相位差φ、模糊相位差φ′和模糊數(shù)n具有如下特點(diǎn)：

1）真實(shí)相位差φ連續(xù)分布，并沿基線矢量方向單調(diào)變化。φ的取值范圍由基線波長比d/λ決定，d/λ越大，取值范圍也越大；基線法平面與定位平面交線對應(yīng)的相位差為0°。

2）模糊相位差φ′是真實(shí)相位差φ向(-π,π]區(qū)間內(nèi)塌縮的結(jié)果。由φ＝0的區(qū)域向外，每當(dāng)φ超出(-π,π]范圍時，就加上或減去2π的整數(shù)倍，使其值保持在(-π,π]區(qū)間中。

3）定位區(qū)域被劃分為一系列具有相同模糊數(shù)n的區(qū)域，區(qū)域的邊界由φ′＝π對應(yīng)等相位差線確定，這些區(qū)域按照-N到N的順序沿基線矢量方向依次排列，隨著基線的旋轉(zhuǎn)不斷變化。

將定位區(qū)域上具有相同模糊數(shù)的區(qū)域稱為“條帶”。隨著基線的旋轉(zhuǎn)，目標(biāo)位置所在條帶不斷變化。將每個測量時刻目標(biāo)所在的條帶n連接起來，就構(gòu)成了一條“路徑”。圖1顯示一條可能的路徑在條帶空間的變化情況。

圖1 每個測量時刻目標(biāo)所在條帶構(gòu)成的路徑

1.2 觀測過程的HMM模型

目標(biāo)所在條帶與其位置和基線轉(zhuǎn)角有關(guān)，將條帶作為隱變量，旋轉(zhuǎn)基線相位差序列的觀測過程可用HMM來描述。

HMM是一種時序的概率模型，描述由一個隱藏的馬爾可夫鏈隨機(jī)生成不可觀測的狀態(tài)隨機(jī)序列，再由各個狀態(tài)生成一個觀測而產(chǎn)生觀測隨機(jī)序列的過程，其中隱藏的隨機(jī)序列稱為狀態(tài)序列,觀測產(chǎn)生的隨機(jī)序列稱為觀測序列。

設(shè)Q＝{q1,q2,q3,…,q K}代表所有狀態(tài)的集合，V＝{v1,v2,v3,…,v M}代表所有可能觀測的集合，HMM可由如下三元組來描述：

式 中，A＝[aij]K×K為 狀 態(tài) 轉(zhuǎn) 移 概 率 矩 陣，aij＝P(st＋1＝qj|st＝qi)表示從時刻t到時刻t＋1由狀態(tài)qi轉(zhuǎn)移到狀態(tài)qj的概率；B＝[bij]K×M為觀測概率矩陣，bij＝P(ot＝v j|st＝qi)表示時刻t在qi狀態(tài)下觀測到v j的概率；Π＝[πi]K×1為初始狀態(tài)概率向量，πi＝P(s1＝qi)代表t＝1時刻處于狀態(tài)qi的概率。

具體到旋轉(zhuǎn)基線干涉儀定位的情況，目標(biāo)所在條帶n是不可觀測的隱藏狀態(tài)，狀態(tài)集合Q＝{-N,-(N-1),…,N-1,N}。觀測集合對應(yīng)相位差測量值，是(-π,π]之間的連續(xù)變量。狀態(tài)轉(zhuǎn)移矩陣A用于描述基線旋轉(zhuǎn)導(dǎo)致的目標(biāo)在各“條帶”之間的轉(zhuǎn)移概率。觀測概率矩陣B表示目標(biāo)位于某一“條帶”時觀測到測量相位差的概率。Π代表目標(biāo)位置最初所在條帶的概率。

由上文可知，任意一條“路徑”對應(yīng)一個狀態(tài)序列s1s2s3…sT，如果能夠由測量相位差序列o1o2o3…o T推斷出目標(biāo)真實(shí)位置對應(yīng)的狀態(tài)序列，就實(shí)現(xiàn)了相位差解模糊。

接下來，通過將“卡爾曼”濾波(Kalman filtering)過程嵌入到“維特比”算法(Viterbi algorithm)中，同步實(shí)現(xiàn)相位差動態(tài)解模糊和濾波定位。

2 “維特比-卡爾曼”濾波定位算法

2.1 假設(shè)真實(shí)相位差已知情況下的“卡爾曼”濾波定位過程

“卡爾曼”濾波[9]是一種經(jīng)典的遞歸濾波器，典型應(yīng)用是從一組包含噪聲的觀測序列中預(yù)測物體的運(yùn)動狀態(tài)。在旋轉(zhuǎn)干涉儀定位中，假設(shè)真實(shí)相位差已知，卡爾曼濾波定位過程如下。

1）狀態(tài)方程和觀測方程

為了討論方便，假定目標(biāo)固定不動，對應(yīng)的狀態(tài)方程為：

狀態(tài)S＝rT，rT＝[x y z]T是目標(biāo)的直角坐標(biāo)矢量。

假設(shè)能夠測量到真實(shí)相位差φ，觀測方程為：φ＝

式中，d為干涉儀基線長度，b＝[x b yb zb]T為單位基線矢量，r1＝[x1y1z1]T是干涉儀所在位置的直角坐標(biāo)矢量。

2）卡爾曼濾波過程

根據(jù)以上的狀態(tài)方程和觀測方程，卡爾曼濾波的具體過程如下。

預(yù)測：

最小預(yù)測均方誤差(MSE)矩陣：

卡爾曼增益：

修正：

最小MSE矩陣：

式中，σ2[n]為真實(shí)相位差φ的測量噪聲方差，H[n]為觀測方程φ＝h(S)的雅克比矩陣

濾波過程中v(n)＝φ[n]-h(S?[n|n-1]稱為“新息”，表示觀測值與預(yù)測值的差。

令?[-1|-1]等于目標(biāo)的初始位置[x0y0z0]T，M[-1|-1]等于目標(biāo)初始位置方差矩陣為相位差φ的觀測方差，通過遞歸計算，就可以得到對系統(tǒng)狀態(tài)S的估計。

3）雅克比矩陣H[n]的推導(dǎo)

雅克比矩陣H[n]是真實(shí)相位差φ關(guān)于目標(biāo)位置rT＝[x y z]T的導(dǎo)數(shù)：

式中，

令：

得：

以上是假設(shè)真實(shí)相位差φ已知情況下的卡爾曼濾波過程，但實(shí)際中干涉儀只能輸出測量相位差φ′，要完成濾波定位過程，就需要實(shí)時解出目標(biāo)所在的條帶n，接下來，通過引入“維特比”算法來解決這一問題。

2.2 “維特比”算法解碼最大概率路徑

“維特比”算法是HMM模型中根據(jù)觀測序列解碼狀態(tài)序列的標(biāo)準(zhǔn)算法。給定觀測序列o1o2o3…o T，要找出最可能的狀態(tài)序列s1s2s3…sT，使式（16）中的條件概率最大：

根據(jù)馬爾可夫齊次性假設(shè)和觀測獨(dú)立性假設(shè)，等價于:

如果直接計算式（17），計算復(fù)雜度將是O(K T)，式中,K為狀態(tài)數(shù)，T為序列長度，計算量非常大?！熬S特比”算法通過在“籬笆網(wǎng)絡(luò)（Lattice）”中應(yīng)用動態(tài)規(guī)劃技術(shù)，將計算復(fù)雜度降低為O(TK2)，有效解決了狀態(tài)解碼問題。

圖2給出了籬笆網(wǎng)絡(luò)的結(jié)構(gòu)。路徑在網(wǎng)絡(luò)中從左向右傳播，由于所有路徑必然經(jīng)過網(wǎng)絡(luò)中的每一列，且如果最大概率路徑經(jīng)過節(jié)點(diǎn)X ij，那么該路徑上從起始時刻到第i個時刻的部分路徑也必然是該時段所有路徑中的最大概率路徑。因此，在任一時刻t，只要考慮K條最大概率路徑即可[10]。

圖2 籬笆網(wǎng)絡(luò)（Lattice）的結(jié)構(gòu)

基于以上性質(zhì)，籬笆網(wǎng)絡(luò)中最大概率的路徑的計算過程如下：

1）初始化。計算第一列各狀態(tài)節(jié)點(diǎn)的初始概率δ1，并記錄每個節(jié)點(diǎn)前一時刻的狀態(tài)ψ1：

式中，πk為狀態(tài)k的初始概率，bk(o1)為狀態(tài)k條件下觀測到o1的概率。由于第一列沒有前一時刻，因此ψ1(k)統(tǒng)一記為0。

2）遞推。對于t＝2,3,…,T，根據(jù)下式遞推計算籬笆網(wǎng)絡(luò)中每一列的K個節(jié)點(diǎn)的概率δt，并記錄每個節(jié)點(diǎn)前一個時刻的狀態(tài)ψt：

式中，ajk為狀態(tài)j轉(zhuǎn)移到狀態(tài)k的概率，bk(ot)為狀態(tài)k條件下觀測到ot的概率。

3）終止。當(dāng)遞推到達(dá)時刻T后，遞推終止，此時，按照下式計算最優(yōu)路徑概率P*和最優(yōu)路徑的最后一個狀態(tài)k*T：

4）回溯。找到最優(yōu)路徑的最后一個狀態(tài)k*T后，根據(jù)下式回溯出最大概率路徑：

最后求得最大概率路徑I*：

接下來，通過將“卡爾曼”濾波步驟嵌入到“維特比”算法過程中，擴(kuò)展上述遞歸過程，完成狀態(tài)轉(zhuǎn)移概率A＝[aij]K×K和觀測概率B t＝[bi(ot)]K×1的計算，同步實(shí)現(xiàn)旋轉(zhuǎn)基線干涉儀相位差解模糊和濾波定位。

2.3 “維特比-卡爾曼”濾波過程

2.3.1 初始化籬笆網(wǎng)絡(luò)和濾波狀態(tài)

根據(jù)式（18）、式（19），初始化籬笆網(wǎng)絡(luò)的狀態(tài)，由于目標(biāo)的位置未知，信號到達(dá)的初始時刻基線的旋轉(zhuǎn)角度也是隨機(jī)的，因此，可以簡單地認(rèn)為目標(biāo)位于各個條帶的概率是相等的，即πi＝1/k,i＝1,2,…,K。同樣，由于目標(biāo)在條帶中的具體位置未知，觀測概率的最初狀態(tài)概率也認(rèn)為是相等的，bk(o1)＝1/k,i＝1,2,…,K，因此：

在籬笆網(wǎng)絡(luò)中，初始時刻的每個狀態(tài)代表目標(biāo)所在的條帶，也就是可能的模糊數(shù)n∈[-N,N]。根據(jù)測量相位差φ′和模糊數(shù)n可以恢復(fù)出K＝2N＋1個可能的真實(shí)相位差φ(k)，利用φ(k)根據(jù)式（6）-（10）執(zhí)行一步濾波定位，得到初始時刻每個狀態(tài)下的一次定位位置S1(k)和最小均方誤差矩陣M1(k)。

2.3.2 遞推計算濾波位置和最大概率路徑

注意到“卡爾曼”濾波和“維特比”最大概率路徑計算過程都采用遞歸的方式計算，將2種過程融合在一起，用籬笆網(wǎng)絡(luò)的狀態(tài)k將測量相位差φ′轉(zhuǎn)換為真實(shí)相位差φ(k)，用其進(jìn)行一步濾波定位，同時，用該濾波定位結(jié)果對應(yīng)的估計相位差φ?(k)與真實(shí)相位差φ(k)之間的相似程度來計算觀測概率，完成一步遞歸計算，具體步驟如下：

1）執(zhí)行一步濾波定位

對于籬笆網(wǎng)絡(luò)中的節(jié)點(diǎn)X tk，利用該節(jié)點(diǎn)對應(yīng)的模糊數(shù)n將測量相位差φ′t轉(zhuǎn)換為真實(shí)相位差φt(k)：

分別以t-1時刻的K個節(jié)點(diǎn)對應(yīng)的濾波位置St-1(k)，k＝1,2,…,K和 最 小 均 方 誤 差 矩 陣M t-1(k),k＝1,2,…,K為狀態(tài)，φt(k)為觀測值，根據(jù)式（6）-（10）執(zhí)行一步濾波定位，得到節(jié)點(diǎn)X tk下的K個 濾 波 位 置Stk(j),j＝1,2,…,K和M tk(j),j＝1,2,…,K。

2）計算狀態(tài)轉(zhuǎn)移概率

狀態(tài)轉(zhuǎn)移概率矩陣A＝[aij]K×K，其中aij表示隨著基線的旋轉(zhuǎn)，在相鄰時刻目標(biāo)在從條帶i轉(zhuǎn)移到條帶j的概率。根據(jù)在時刻t對系統(tǒng)狀態(tài)了解程度的不同，按照以下幾種情況來計算狀態(tài)轉(zhuǎn)移概率。

①無任何先驗(yàn)知識的情況

在目標(biāo)位置無任何先驗(yàn)知識，且時刻t-1到t之間的時間間隔隨機(jī)產(chǎn)生的情況下，隨著基線的旋轉(zhuǎn)，目標(biāo)位置可能在任意2個條帶之間轉(zhuǎn)移，此時，可以簡單認(rèn)為轉(zhuǎn)移概率是平均分布的，即：

②目標(biāo)粗略位置已知的情況

假設(shè)已經(jīng)獲取了目標(biāo)的粗略位置，例如在濾波定位過程中不斷接近目標(biāo)真實(shí)位置，那么，在t-1時刻，根據(jù)位置利用式（2）和式（5）可以估計出目標(biāo)可能的條帶范圍，例如連續(xù)的3個條帶(i-1,i,i＋1)，基線旋轉(zhuǎn)到t時刻，同樣可以計算出另一個可能的條帶范圍，例如(j-1,j,j＋1)，這樣，從t-1時刻到t時刻，轉(zhuǎn)移概率可以表示為：

具體的條帶范圍需要根據(jù)估計位置與真實(shí)位置的接近程度來設(shè)置，越接近，則可能的條帶范圍越小，否則，需要給多個可能的條帶賦予非零概率。

③相對相位差無模糊的情況

當(dāng)t-1和t時刻的真實(shí)相位差滿足|φt-φt-1|≤π時，稱為相位差相對無模糊。此時，可以恢復(fù)出測量相位差φ′t＋1和φ′t之間真實(shí)的差值。文獻(xiàn)[1]中證明，當(dāng)|Δω|≤2arcsin(λ(4d))時，干涉儀覆蓋范圍內(nèi)任意位置的相位差相對無模糊。此時，轉(zhuǎn)移概率可以表示為：

當(dāng)φ′t-φ′t-1≤-π時，

當(dāng)φ′t-φ′t-1＞π時，

當(dāng)-π＜φ′t-φ′t-1≤π時，

3）計算觀測概率

觀測概率bk(ot)表示t時刻在狀態(tài)k下觀測到ot的概率。相位差是由位置決定的，經(jīng)過一次濾波后，在節(jié)點(diǎn)X tk下有K個濾波位置Stk(j)，j＝1,2,…,K，根據(jù)式（5）可計算出對應(yīng)的K個相位差：

將φtk(j)與真實(shí)相位差φt(i)之間的相似程度作為計算觀測概率bk(ot)的依據(jù)。首先利用sigmoid函數(shù)將φtk(j)-φt(i)變換到(-1,1)的范圍內(nèi)：

然后按照標(biāo)準(zhǔn)正態(tài)分布計算觀測概率：

4）保存最大概率狀態(tài)

對于進(jìn)入節(jié)點(diǎn)X tk的K條路徑，選取當(dāng)前最大的部分路徑概率對應(yīng)的一條作為X tk的狀態(tài)。根據(jù)下式計算籬笆網(wǎng)絡(luò)中時刻t的K個節(jié)點(diǎn)的概率δt，并記錄每個節(jié)點(diǎn)前一個時刻的狀態(tài)ψt：

記錄最大概率路徑對應(yīng)的濾波位置和最小均方誤差矩陣：

5）遞歸終止和路徑回溯

當(dāng)遞推到達(dá)時刻T后，遞推終止，此時，按照下式計算最優(yōu)路徑概率P*、最優(yōu)路徑的最后一個狀態(tài)k*T和濾波定位位置S*T：

找到最優(yōu)路徑的最后一個狀態(tài)k*T后，根據(jù)下式回溯出最大概率路徑：

最后求得最大概率路徑I*：

至此，完成了“維特比-卡爾曼”濾波定位的全過程。該算法利用“維特比”算法根據(jù)測量相位差序列解碼最優(yōu)路徑，相當(dāng)于解出每次測量時對應(yīng)的模糊數(shù)，從而利用真實(shí)相位差完成濾波定位。定位性能方面與“卡爾曼”濾波一致，此處不再贅述。

3 仿真分析

仿真場景設(shè)置如下：目標(biāo)位置的經(jīng)緯度為[123.117 5°,23.553 1°]，旋 轉(zhuǎn) 基 線 位 置 的 經(jīng) 緯 度 為[125°,28°]，在500 km高度水平放置，基線長度為2 m，按照360°/s的速度旋轉(zhuǎn)。輻射源信號脈沖重復(fù)間隔在[1 000μs,10 000μs]區(qū)間內(nèi)隨機(jī)選取，脈沖列長度為100。濾波初始位置經(jīng)緯度隨機(jī)選取為[124.132 7°,22.080 9°]。

采用“維特比-卡爾曼”濾波算法進(jìn)行處理，經(jīng)過100步濾波后得到定位位置的經(jīng)緯度為[123.114 7°，23.551 3°]，定位誤差為351 m，與利用真實(shí)相位差進(jìn)行卡爾曼濾波后的定位結(jié)果一致。濾波定位每一步的結(jié)果如圖3所示。其中，圖3(a)給出了濾波定位過程中估計位置不斷接近真實(shí)位置的情況；圖3(b)是濾波過程中最大概率路徑對應(yīng)的相位差“新息”的變化情況，可見，隨著濾波的進(jìn)行，相位差“新息”在不斷減小，說明估計相位差與真實(shí)相位差越來越接近；圖3(c)給出了籬笆網(wǎng)格中由最終狀態(tài)回溯的所有路徑（13條）的每一步的定位結(jié)果，可見，非真實(shí)路徑對應(yīng)的定位結(jié)果始終是發(fā)散的，而真實(shí)路徑對應(yīng)的結(jié)果逐漸收斂到真實(shí)位置（紅色星號）；圖3(d)給出了在維特比遞推過程中籬笆網(wǎng)絡(luò)中每一列最大概率的變化情況，可見，隨著遞歸的進(jìn)行，真實(shí)路徑上的概率比重不斷增加，而其它路徑上的概率由于采用了錯誤的相位差進(jìn)行濾波，無法收斂到真實(shí)位置，導(dǎo)致概率比重逐漸減小，通過選取最大概率的狀態(tài)并回溯后，找出真實(shí)路徑，從而實(shí)現(xiàn)相位差解模糊和定位。

圖3 “維特比-卡爾曼”濾波算法一次仿真中的定位過程

4 結(jié)束語

本文研究了旋轉(zhuǎn)長基線干涉儀利用模糊相位差直接定位問題。首先用HMM對旋轉(zhuǎn)基線定位過程建模，利用維特比算法根據(jù)測量相位差解碼模糊數(shù)序列，在籬笆網(wǎng)格計算過程中引入卡爾曼濾波步驟，同步計算目標(biāo)估計位置和最優(yōu)路徑概率，解碼最優(yōu)狀態(tài)序列，實(shí)現(xiàn)了利用模糊相位差直接濾波定位。最后用仿真實(shí)驗(yàn)驗(yàn)證了算法的有效性?！?/p>