徐海源，蘇成曉，黨衛(wèi)，高鳳輝

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

0 引言

星載無源定位系統(tǒng)具有覆蓋范圍大、不受地域限制等優(yōu)點，在輻射源監(jiān)測和定位中得到了廣泛應(yīng)用。其中，單星無源定位系統(tǒng)具有成本低、技術(shù)成熟等優(yōu)點［1-3］，利于靈活配置和大規(guī)模組網(wǎng)應(yīng)用。針對空中輻射源，單星定位系統(tǒng)需要假設(shè)高程進行定位，將帶來較大的定位誤差［4］。采用三星定位體制，通過增加測向等信息，并利用跟蹤濾波方法，可實現(xiàn)空中目標三維位置和速度的高精度估計［5-6］。

針對雙星體制對空中目標的定位，文獻［7］采用主星高精度測向和雙星時差聯(lián)合實現(xiàn)空中目標高精度定位，但需要較高的測向精度，傳統(tǒng)的干涉儀測向系統(tǒng)難以滿足要求。針對單星測向定位體制的多顆衛(wèi)星對同一目標存在共視的情況，文獻［4］提出了基于雙星測向結(jié)果求解空中目標位置的方法，典型場景下定位誤差仿真結(jié)果為幾千米，可提高空中目標的定位精度，但基于定位結(jié)果難以實現(xiàn)目標高程和運動速度的高精度估計。

本文針對2顆安裝二維正交干涉儀的單星定位系統(tǒng)對空中目標共視的應(yīng)用場景，從輻射源信號參數(shù)測量模型和處理流程出發(fā)，提出了基于雙星干涉儀相位差和時差融合的空中動目標定位跟蹤方法，以期提高空中目標的位置和速度估計精度。

1 雙星干涉儀相位差和時差融合跟蹤處理流程

1.1 雙星定位系統(tǒng)輻射源信號參數(shù)測量模型

2顆衛(wèi)星均安裝二維正交干涉儀的測向定位系統(tǒng)，在共視區(qū)域內(nèi)接收到空中目標輻射的信號，每顆衛(wèi)星可測得輻射源信號到達接收天線的相位差信息，基于2顆衛(wèi)星獲取的信號參數(shù)還可以進一步得到時差信息。雙星定位系統(tǒng)示意如圖1所示。

圖1 雙星定位系統(tǒng)示意

在單星二維正交干涉儀測向系統(tǒng)中，為簡化描述，假設(shè)星體坐標系與干涉儀坐標系重合，天線A、B分別安裝于星體坐標系Y、X軸，Z軸對地，天線O位于星體坐標系原點。原始觀測量為輻射源信號到達天線A和天線B與天線O之間的相位差，分別記為φA、φB。假設(shè)星體坐標系下的輻射源位置坐標為(xb，yb，zb)，則φA、φB可表示為［2］：

式中，n=，其中d為干涉儀基線長度，λ為信號波長。

對于2顆衛(wèi)星收到的輻射源同一脈沖信號，其時差τ可表示為：

式中，c為電磁波傳播速度，r0、r1分別為輻射源與衛(wèi)星0、衛(wèi)星1的距離。

1.2 定位跟蹤處理流程

本文針對運動目標輻射源發(fā)射脈沖信號的情況，考慮2種定位跟蹤處理流程：一種是基于單脈沖相位差測量數(shù)據(jù)和時差測量數(shù)據(jù)直接進行跟蹤濾波；第二種是針對多個脈沖序列按一定時間節(jié)拍進行相位差和時差參數(shù)估計后，再按時間節(jié)拍進行跟蹤濾波。

基于參數(shù)估計理論，直接采用原始測量數(shù)據(jù)進行跟蹤濾波能夠得到更高的精度，但在實際應(yīng)用中，考慮到處理效率，經(jīng)常采用第二種處理方式。2種處理流程如圖2所示。

圖2 2種定位跟蹤處理流程

2種處理流程的主要差別在于是否增加按時間節(jié)拍進行參數(shù)估計的環(huán)節(jié)，這里采用最小二乘法進行相位差和時差估計。

1.3 基于最小二乘法的相位差和時差估計

由于空中目標與2顆衛(wèi)星之間存在較大的相對運動速度，因此，在一定時間段內(nèi)對相位差和時差進行估計時，需要考慮切向相對運動速度帶來的相位差變化率和徑向相對運動速度帶來的時差變化率帶來的影響［3，8］。

1）相位差估計

利用觀測時間段內(nèi)的N個脈沖的相位差測量數(shù)據(jù)估計初始時刻的相位差φ0及相位差變化率φ?，最小二乘估計法為：

式中，φ=[φ0，φ1，…，φN-1]T為相位差測量數(shù)據(jù)構(gòu)成的向量，H為N×2矩陣：

2）時差估計

設(shè)觀測時間內(nèi)2顆衛(wèi)星均接收到同一輻射源的N個脈沖，脈沖配對后第k個脈沖對的到達時間（TOA）為t0k和t1k，τk=t1k-t0k表示第k個脈沖到達2個觀測站的時間差。

記τ=[τ0，…，τN-1]T為時差序列測量數(shù)據(jù)構(gòu)成的向量，則采用最小二乘估計法對初始時刻時差τ0及時差變化率τ?的估計公式為：

式中，H為與式（4）相同的N×2矩陣。

1.4 相位差及時差估計精度

最小二乘估計的方差為：

對于相位差估計，R=，其中為單脈沖的相位差測量方差。因此，相位差估計的方差可以進一步簡化為P=(HTH)-1。

假設(shè)相鄰2個脈沖的時間間隔Δtk=tk-tk-1相等，且等于Δt，通過進一步數(shù)學變換可得：

由此得到初始相位差φ0的估計誤差方差為=2(2N-1)/(N(N+1))。

2 雙星干涉儀相位差和時差融合跟蹤模型

2.1 狀態(tài)方程

考慮空中目標運動模型為勻速模型（CV），定義k時刻的狀態(tài)量為WGS-84坐標系下的目標位置及速度，即Xk=[xkykzkvxkvykvzk]T，則狀態(tài)方程為：

式中，A為狀態(tài)轉(zhuǎn)移矩陣，定義為：

式中，Δt為觀測時間間隔，wk-1為過程噪聲，其均值為0，協(xié)方差矩陣為Qk-1。

對于CV模型，狀態(tài)噪聲協(xié)方差Qk-1表示為［9］：

式中，I3為三維單位陣，為速度在各方向的隨機誤差的方差。

2.2 觀測方程

為了在WGS-84坐標系下進行目標狀態(tài)跟蹤濾波，首先需要將式（1）的相位差觀測量轉(zhuǎn)換到WGS-84坐標系。記星體坐標系相對WGS-84系的旋轉(zhuǎn)矩陣為M，其推導(dǎo)過程參見文獻［6］。記WGS-84坐標系下衛(wèi)星位置為（xs，ys，zs），輻射源位置為(x，y，z)，由于

式中，mij為坐標旋轉(zhuǎn)矩陣M中的元素。于是，式（1）可表示為：

記2顆衛(wèi)星獲得的輻射源信號相位差測量數(shù)據(jù)（或估計數(shù)據(jù)）序列為φ0，Ak、φ0，Bk、φ1，Ak、φ1，Bk，時差測量數(shù)據(jù)（或估計數(shù)據(jù)）序列為τk，則觀測方程為：

式中，(xi，k，yi，k，zi，k)為WGS-84坐標系下的衛(wèi)星位置，i=0，1；(xk，yk，zk)為WGS-84坐標系下輻射源的位置，ri，k=((xk-xi，k)2+(yk-yi，k)2+(zk-zi，k)2))1/2。mij為衛(wèi)星0星體坐標系相對于WGS-84坐標系的旋轉(zhuǎn)矩陣M的元素；m'ij為衛(wèi)星1星體坐標系相對于WGS-84坐標系的旋轉(zhuǎn)矩陣M'的元素；vk為均值為0，協(xié)方差矩陣為Rk的高斯白噪聲。

在基于單脈沖參數(shù)進行跟蹤濾波過程中，噪聲協(xié)方差Rk=diag分別為單脈沖時差、相位差測量誤差的方差。

若按照一定時間節(jié)拍采用脈沖積累先估計時差、相位差再進行跟蹤濾波，則需要根據(jù)積累脈沖數(shù)量N計算時差、相位差的誤差，根據(jù)1.4節(jié)分析可知，Rk=2(2N-1)/(N(N+1))·diag()。

3 EKF跟蹤算法

由于觀測方程為非線性方程，需要采用非線性濾波方法提高對運動目標位置和速度的估計精度。本文采用擴展卡爾曼濾波（EKF）進行空中運動目標跟蹤。

3.1 濾波初始化

1）目標狀態(tài)初始值

可以通過假定高程，采用其中一顆衛(wèi)星單星測向定位的方法計算得到初始位置(x0，y0，z0)［2］。目標初始速度可以設(shè)置為(0，0，0)，即目標狀態(tài)的初始值為X0/0=[x0y0z0000]T。

2）初始預(yù)測協(xié)方差矩陣

基于單星干涉儀相位差測量誤差、高程假設(shè)誤差、衛(wèi)星位置誤差等信息，可以計算得到初始點定位的協(xié)方差矩陣P1［2］。假定目標可能的最大航速為vmax，記P2=diag。設(shè)置初始預(yù)測協(xié)方差矩陣P0/0=diag(P1，P2)。

3.2 EKF跟蹤濾波步驟

基于上文建立的目標狀態(tài)方程和觀測方程，對目標進行EKF跟蹤過程如下：

1）目標狀態(tài)預(yù)測

2）協(xié)方差矩陣預(yù)測

3）計算Zk在預(yù)測點處的雅克比矩陣

4）卡爾曼增益矩陣計算

5）更新目標狀態(tài)

6）更新協(xié)方差矩陣

需要說明的是，最簡單的二維正交干涉儀測量的相位差會出現(xiàn)模糊現(xiàn)象，在實際應(yīng)用中可以采用多基線或長短基線結(jié)合等方法解模糊。在本文算法和步驟的描述中，均假定相位差觀測量為無模糊值。

4 仿真分析

為驗證提出的基于雙星干涉儀相位差和時差融合的空中運動目標跟蹤方法的有效性，建立典型場景進行了仿真分析。2顆衛(wèi)星軌道高度均為700 km，星間距80 km，每顆衛(wèi)星上均安裝二維正交干涉儀，基線長度3 m；目標高程為10 km，速度為300 m/s的勻速直線運動；目標輻射源頻率為3 GHz，輻射源脈沖重復(fù)頻率（PRF）為100 Hz。衛(wèi)星姿態(tài)測量誤差0.05°（3軸，1σ），干涉儀相位差測量誤差為20°（1σ），時差測量誤差為50 ns（1σ）。在雙星共同覆蓋區(qū)域內(nèi)選取一個目標位置，對基于脈沖測量參數(shù)直接進行EKF濾波的處理方法和基于時間節(jié)拍參數(shù)估計后進行EKF濾波的處理方法分別進行100次Monte-carlo仿真，并進行比較。同時，針對文獻［4］的雙星測向定位和文獻［7］的單星測向+雙星時差定位體制下的跟蹤濾波情況也進行了仿真比較。

1）單脈沖測量數(shù)據(jù)直接EKF仿真結(jié)果

基于單脈沖測量數(shù)據(jù)直接EKF濾波方法，針對觀測量分別采用雙星干涉儀相位差+時差、雙星干涉儀相位差、單星干涉儀相位差+時差3種情況，對目標位置（含高程）、航速和航向估計誤差進行了仿真比較，如圖3—5所示。

圖3 單脈沖直接EKF定位誤差

從圖3可以看出，本文提出的基于雙星干涉儀相位差+時差觀測量的EKF濾波方法相比另外2種方法具有更快的收斂速度和更高的位置估計精度，在仿真場景下跟蹤10 s后定位誤差小于1 km（含高程誤差），跟蹤100 s后定位誤差接近0.2 km（含高程誤差）。

從圖4—5可以看出，觀測量采用雙星干涉儀相位差+時差、單星干涉儀相位差+時差與采用雙星干涉儀相位差的EKF濾波方法相比具有更短的速度和航向估計收斂時間，原因在于衛(wèi)星和目標之間存在較大的相對速度，原始觀測量中的時差參數(shù)隱含了時差變化率信息，體現(xiàn)出了目標運動速度相關(guān)信息。

圖4 單脈沖直接EKF速度誤差

圖5 單脈沖直接EKF航向誤差

2）基于時間節(jié)拍參數(shù)估計的EKF仿真結(jié)果

基于上述仿真參數(shù)，采用1 s時間節(jié)拍進行參數(shù)估計（即每個節(jié)拍100個脈沖）后再進行EKF的方法，針對觀測量分別采用雙星干涉儀相位差+時差、雙星干涉儀相位差、單星干涉儀相位差+時差3種情況，對目標位置（含高程）、航速和航向估計誤差進行了仿真比較，如圖6—8所示。

圖6 基于時間節(jié)拍EKF定位誤差

圖7 基于時間節(jié)拍EKF速度誤差

圖8 基于時間節(jié)拍EKF航向誤差

3種觀測量的仿真結(jié)果比較與基于單脈沖測量數(shù)據(jù)直接EKF的結(jié)論基本一致，采用雙星干涉儀相位差+時差作為觀測量的濾波方法相比另外2種觀測量方法具有更高的目標位置估計精度和更短的速度、航向收斂時間，但收斂后的位置、速度和航向誤差均高于基于單脈沖測量數(shù)據(jù)直接EKF的方法。

3）2種處理流程的仿真結(jié)果比較

針對觀測量采用雙星干涉儀相位差+時差的EKF濾波方法，對單脈沖測量數(shù)據(jù)直接EKF和基于時間節(jié)拍參數(shù)估計EKF 2種處理流程的定位誤差仿真結(jié)果進行比較，如圖9所示，其中時間節(jié)拍設(shè)置了0.1 s、1 s 2種情況。

圖9 單脈沖EKF及不同時間節(jié)拍EKF定位誤差比較

從仿真結(jié)果可以看出，基于時間節(jié)拍參數(shù)估計的EKF處理方法在目標位置估計精度和收斂速度方面均低于基于單脈沖測量數(shù)據(jù)直接EKF的方法，且隨著時間節(jié)拍增大，估計精度下降明顯。原因在于，一是按時間節(jié)拍進行相位差和時差參數(shù)估計后，損失了相位差變化率和時差變化率等信息；二是用于時間節(jié)拍內(nèi)參數(shù)估計的最小二乘估計方法并非最優(yōu)估計方法，導(dǎo)致濾波過程中引入新的誤差。

雖然基于時間節(jié)拍參數(shù)估計的EKF處理方法在目標位置、速度等狀態(tài)估計精度上偏低，但其在計算量上相對于基于單脈沖測量數(shù)據(jù)直接EKF處理方法有明顯優(yōu)勢，適合應(yīng)用于處理資源受限的場景。

5 結(jié)束語

本文針對均具備二維正交干涉儀測向的2顆衛(wèi)星對空中輻射源定位跟蹤問題，提出了一種基于干涉儀相位差和時差融合的跟蹤濾波方法。首先，從輻射源信號參數(shù)原始觀測量出發(fā)，提出了2種定位跟蹤處理流程；其次，在WGS-84坐標系下建立了由時差和相位差構(gòu)成的觀測方程，并采用EKF方法進行跟蹤濾波。仿真結(jié)果表明，本文提出的目標跟蹤方法能夠?qū)崿F(xiàn)未知高程空中運動目標位置、航速和航向的高精度估計，相對于已有方法在收斂速度和估計精度上有明顯提高?；诒疚牡姆抡娣治鼋Y(jié)論，在處理資源充分的情況下，采用基于輻射源單脈沖測量數(shù)據(jù)直接EKF濾波方法能夠得到更高的空中目標定位跟蹤精度。