基于MCIS技術的相位差變化率單站無源定位

邢懷璽,張宇暉,陳游,*,周一鵬,何文波

1. 空軍工程大學 航空工程學院,西安 710038 2. 空軍裝備部 駐成都第三軍事代表室,成都 610000

單站無源定位技術由于避免了多站定位的時間同步和數據融合等問題,具有良好的機動性和靈活性,對于提高飛機的突防能力、生存能力和隱蔽攻擊能力都有重要的現實意義,一直是無源定位領域的研究熱點[1-5]。最早的單站無源定位利用輻射源來波方向(DOA)實現交叉定位[6-8],但這一方法受限于來波方向的測量精度,利用長基線干涉儀(LBI)可以精確獲得非合作輻射源信號到達陣元間的相位差及相位差變化率。許耀偉和孫仲康[9]提出了利用干涉儀測量相位差變化率實現對目標的測距,研究表明,采用相位差變化率定位較只測向定位大幅提高了定位精度和速度。相位差變化率是相位差對時間的導數,因此不需要對相位差的解模糊運算,同時可以消除天線間通道幅度相位不一致帶來的相位差系統偏差。王強等[10]提出了僅利用相位差變化率單站無源定位方法,只需要兩陣元構成的長基線干涉儀測量來波的相位差變化率,便可實現對三維空間中固定輻射源的無源定位,定位系統設備相對簡單且成本低。

基于相位差變化率的單站無源定位實質上是一個非線性最優估計問題,一般很難直接獲得目標位置的解析解。目前定位算法主要分為迭代法和搜索法兩類,其中迭代法應用最多的是非線性最小二乘法和擴展卡爾曼濾波(EKF)以及性能更好的容積卡爾曼濾波(CKF)、不敏卡爾曼濾波(UKF)算法等。文獻[11]采用非線性最小二乘估計的方法,獲取多個時刻的相位差變化率迭代遞推目標的位置；單月暉等[12]將EKF和修正協方差擴展卡爾曼濾波(MVEKF)方法用于測相位差變化率的無源定位中,驗證了MVEKF方法有良好的收斂性能；郭福成等[13]首次將MVEKF應用于無源定位跟蹤中,無需計算測量方程的修正函數,適用范圍更廣；文獻[14]利用一種基于改進的變尺度最小斜度不敏卡爾曼濾波(SMSS-UKF)估計目標位置,在低信噪比條件下具有較好的適應性,但這一類方法非常依賴于初始化的結果,在狀態估計誤差較大的情況下,尤其是濾波初期,受測量噪聲影響較大,估計過程中協方差易出現病態, 致使濾波結果不穩定。此外,最大似然估計(ML)[15]是一種漸進無偏的估計算法,被廣泛用于定位問題中,可以達到克拉美羅下界(CRLB)。迭代方法也可以用來解決ML問題,常用的方法有牛頓法、共軛梯度法等。這些方法對初始估計同樣很敏感,在高測量噪聲下可能會導致收斂至局部最優或發散。另一個解決方案是將ML問題松弛為凸優化問題,如半定規化松弛(SDP)[16]和二階錐規劃松弛(SOCP)[17]。搜索類算法優勢在于當缺少目標位置先驗信息也可以通過構建代價函數,搜索代價函數極值點獲得目標位置。文獻[10]在沒有目標位置先驗信息的前提下,用ML代替最大后驗估計構造似然函數,將目標位置估計轉化為一個非線性非凸優化問題,然后利用網格搜索法得到全局的最優解,算法穩定性較高。然而,要想實現較高的定位精度,需要足夠大的網格密度,無法避免算法的運算量隨著網格分辨率增加呈指數增長這一固有矛盾,同時難以實現工程化的實時應用。

近年來,蒙特卡洛(MC)技術[18-19]廣泛應用于信號處理領域,其中重要性抽樣是計算多維積分強有力的工具。文獻[20]提出一種基于重要性抽樣的最大似然估計方法,利用重要性抽樣(IS)技術實現了DOA的ML估計,大大降低了ML估計的計算量；文獻[21]采用馬爾科夫鏈蒙特卡洛(MCMC)方法求解似然函數的全局極大值,得到時差-頻差聯合估計,估計精度優于傳統的互模糊函數(CAF)算法,而且克服了迭代算法的初值依賴和收斂問題,同時又避免了網格搜索大量的計算過程；文獻[22]利用一種采用蒙特卡洛重要性抽樣(MCIS)方法近似多維積分的方法解決到達時間差(TDOA)的源定位問題,通過概率形式選取合適的樣本,以樣本均值代替目標位置估計,取得理想的效果。

基于此,本文提出一種基于相位差變化率的高精度,低復雜度的單站無源定位方法。首先給出高斯測量噪聲下的目標輻射源位置的最大似然估計,然后引入Pincus定理導出多維積分下的全局最優解,最后通過MCIS技術可以得到一個近似的全局解。這一過程需要目標輻射源位置粗略的初始估計值。一般的迭代方法本身也需要初值估計,而最小二乘算法能夠根據定位方程直接解算目標位置。所以本文利用加權最小二乘(WLS)算法[23-24]獲得目標位置封閉形式解作為本文定位方法的初值估計。本文將相位差變化率測量方程在目標位置估計初值處一階泰勒級數展開得到線性化方程,構造一個高斯分布的PDF作為重要性函數。以重要性函數的概率分布抽取足夠多的樣本,根據大數定律,以樣本估計的期望值作為全局最優解,從而得到輻射源位置。結果表明,MCIS方法能夠提供最優的性能,不需要大量計算,在一定的噪聲水平下定位誤差逼近CRLB,并且對初始估計誤差的敏感性較低。

1 基于相位差變化率的定位原理

1.1 定位模型

假設單個觀測平臺上垂直機身安裝一組二陣元的一維單基線干涉儀,在一段觀測時間內,干涉儀多次測量來波方向的相位差變化率,定位場景如圖1所示。

空中觀測平臺的位置速度和姿態角信息可以從慣導系統或GPS中獲取。假設固定目標輻射源的位置狀態為xT=[xT,yT,zT]T,在t1～tn

圖1 單基線測相位差變化率的定位模型示意圖Fig.1 Schematic diagram of positioning model for measuring change rate of phase difference with single baseline

(1)

式中:κ=2πdf/c,f為輻射源信號頻率,d為基線長度,c為光速;Lk=rk/rk。相應的tk時刻的相位差變化率為

(2)

(3)

(4)

將K次測量的相位差變化率用向量形式表示

(5)

1.2 僅測相位差變化率的三維可觀測性分析

選機腹下機身軸和機翼軸的交點O作為參考原點,沿著機身軸方向并指向機頭的方向為參考X軸,OY軸與OX軸在同一水平面上且指向左側機翼方向,按照右手關系作OZ軸垂直于OXY平面且指向上方,建立三維直角坐標系O-XYZ,如圖2所示。

假設觀測平臺水平直線勻速運動,對于靜止目標,記k時刻,以OY方向為基準,機載觀測平臺和輻射源之間的方位角為βk,以平面O-XY為基準,機載觀測平臺和輻射源之間的俯仰角為εk;觀測平臺到輻射源徑向距離在水平面的投影距離為rpk。

所以k時刻相位差用方位角和俯仰角表示為

(6)

相位差變化率進一步表示為

(7)

(8)

(9)

圖2 觀測平臺與輻射源相對位置示意圖Fig.2 Schematic diagram of relative position between observing station and radiation source

(10)

(11)

其中：rpk=(xT-xok)2+(yT-yok)2。

將式(8)～式(11)代入式(7)可得

f(xk,yk,zk)

(12)

2 重要性抽樣定位方法

無源定位實質上是一個非線性最優估計問題,在基于相位差變化率的單站無源定位中,相位差變化率與目標位置存在非線性的約束關系,通過多個時刻觀測值獲取目標位置狀態信息。目前的研究中以非線性最小二乘(NLS)、EKF等其他改進的非線性跟蹤濾波算法為主體,但這類算法的性能嚴重依賴初始狀態估計,初始階段定位誤差較大會出現發散的現象,同時受測量噪聲影響大,導致濾波結果不穩定。除此之外,也可以通過構建代價函數,通過網格搜索代價函數最小值獲取目標位置信息,ML就是其中之一。

目標輻射源位置ML可表示為矢量形式

(13)

式中：Q=E(eeT)為測量誤差的協方差。

2.1 ML估計的全局解

ML問題是一種非線性、多維的極值優化問題,需要全局極值的多維搜索,因此其計算量巨大。因此在實際的操作中往往將其降維處理,這里將ML問題轉化為一維優化問題,這里介紹Pincus定理。

(14)

(15)

在實際的測量中,需確定有足夠大的范圍保證x在范圍邊界以內。重要性抽樣的關鍵在于重要性函數的選取,合適的重要性函數可以減少計算量,大大提高算法的計算速度,較差的重要性函數則會影響仿真時間,同樣影響優化結果的質量。這里,定義-λf(x)的歸一化函數為

(16)

(17)

利用式(17)計算x的最佳估計值,實際上是求解一個多維多重積分問題,直接計算非常復雜且困難。由大數定律可知,當樣本數量足夠多時,樣本均值的方差越小,樣本均值逼近真實積分值。因此針對上述問題可以通過p(x)生成多個x的樣本,用樣本均值代替式(17)的復雜積分。

2.2 蒙特卡洛重要性抽樣技術

蒙特卡洛重要性抽樣技術是解決多維積分的強有力工具。蒙特卡洛方法建立簡單的概率統計模型可以求解某個隨機變量的數字特征,通過合適的抽樣方法抽取樣本,求取樣本的統計量,從而得到參數估計。顯而易見,蒙特卡洛方法的關鍵在于抽樣環節,抽樣過程中并不能對服從任意概率分布的隨機變量直接進行抽樣,然而重要性抽樣可以把難以直接抽樣的概率分布,轉換為從易于抽取樣本的重要函數進行抽樣,簡單方便,且可以降低估計方差。其基本思想是通過一個相對簡單的分布函數生成一定數量的隨機數,并將其加權平均來近似計算目標分布函數的數學期望。下面以蒙特卡洛積分為例介紹重要性抽樣過程：

假設一個復雜積分

(18)

式中：p(θ)為θ的概率密度函數。如果q(θ)為積分域上的概率密度函數,那么

(19)

(20)

根據大數定律,式(20)的積分被近似為

(21)

2.3 測相位差變化率定位的重要性函數選擇

對于重要性函數的選擇,這里希望能夠找到一種近似于式(16)的概率密度函數來直接計算積分。考慮到相位差變化率的測量誤差服從高斯分布,同時包含目標的位置信息,除此之外,高斯分布因為只需要均值和協方差,容易構造,樣本的抽取比較簡單,所以本文選擇高斯分布的概率密度函數作為重要性函數。

由式(4)可得

(22)

(23)

(24)

為了簡化公式,令

所以式(23)可以改寫成

(25)

因為是關于xT的線性函數,可以將K次觀測的線性方程用矩陣形式表示

e=HxT-F

(26)

式中：

根據式(13)得到目標位置的最大似然估計近似為

(27)

因此,Pincus定理構造重要性函數為

q(xT)=

(28)

λ的大小關系著抽樣區域的大小,為了構造高斯分布的概率密度函數,這里變換似然函數：

(HxT-F)TQ-1(HxT-F)=[xT-

(HTQ-1H)-1HTQ-1F]T·(HTQ-1H)·[xT-(HTQ-1H)-1HTQ-1F]+FTQ-1F-FTQ-1H(HTQ-1H)-1HTQ-1F=

(29)

其中：

(30)

R=(HTQ-1H)-1

(31)

ζ=FTQ-1F-FTQ-1H(HTQ-1H)-1HTQ-1F

(32)

這里因為ζ項與xT無關,所以q(xT)分子中常數項exp(-λ1ζ)與分母對應項相互抵消。重要性函數q(xT)可以改寫為

q(xT)=

(33)

所以式(33)進一步簡化為

(34)

從式(34)中可以看出高斯分布的協方差與密切相關,λ1的大小決定高斯分布重要性函數的幅度和寬度,當λ1很大時,重要性函數很窄,容易生成大量樣本接近高斯分布均值,當均值明顯偏離輻射源真實位置時,會使最終的定位結果產生偏差。參考文獻[19],這里取λ1=0.5。第4節中討論了參數λ1對MCIS算法定位精度的影響。

3 重要性抽樣定位

3.1 輻射源位置估計

在確定重要性函數之后,可以通過逆變換采樣的方法選擇樣本集,這里不作贅述。繼而很容易獲得目標輻射源位置信息的最大似然估計,可以通過足夠多的樣本均值來近似計算,即

(35)

λf(xTm)]

(36)

顯然似然函數是一個多峰函數,存在局部極值。為了明確區分全局最優和局部極值,使樣本均值更加趨近于真實值,在這里對其指數歸一化處理得到標準化的重要性系數。

(37)

從式(30)中可以看出,參數λ的選擇對于標準化重要性系數ω′(xTm)估計的準確度有著重要影響。根據Pincus定理,當λ→∞時,式(16)的似然函數呈現n維的Dirac-delta函數形狀,使全局最大值和其他極值差異明顯。實際上,對于本文的定位模型,輻射源真實位置附近區域只存在單一的極值點,所以λ的取值不必太大,只需要保證滿足相應的計算要求。經過實驗分析,當λ=0.2時已經能夠滿足計算要求。第4節也給出了參數λ對MCIS算法估計性能的影響。

3.2 算法框架

步驟1通過加權最小二乘(WLS)算法獲取初始目標位置估計,計算系數矩陣H和F。

步驟3設置樣本數量M,利用逆變換采樣定理(Probability Integral Transformation Theorem)生成樣本,生成M個獨立且服從[0,1]均勻分布的隨機向量um,對于每一個um,在q(xT)找到最接近的值,對應的xTm為生成的樣本。

步驟4計算標準化重要性系數ω′(xTm)并根據式(38)估計目標位置

(38)

3.3 定位誤差的CRLB

因為CRLB清楚地刻畫了參數測量誤差下的定位誤差及其分布,所以推導相位差變化率定位誤差的CRLB作為定位方法性能的評價指標。根據CRLB的定義：

(39)

(40)

因此,可得

(41)

4 仿真分析

4.1 實驗1：MCIS算法估計輻射源位置

4.2 實驗2：MCIS參數對定位性能的影響

后面的討論將建立在實驗1模擬場景下進行。同時,為了更好地衡量算法的定位精度,將目標位置估計的均方根誤差(RMSE)作為指標來衡量定位的精度。RMSE定義為

(42)

4.2.1 參數λ和λ1對定位精度的影響

改變參數λ和λ1,進行200次蒙特卡洛實驗,討論MCIS的定位精度。圖4和圖5為不同參數λ和λ1設置下的定位誤差的RMSE。從圖4可以看出,當λ≥0.2時定位結果趨于穩定,說明無需達到理論上的無窮大,很容易滿足本文定位模型的要求。圖5中λ1在[0,5]區間內改變,發現λ1的大小對定位精度幾乎沒有影響。也進一步驗證了實驗1中λ和λ1的選取是合理有效的。

圖4 λ對MCIS定位精度的影響Fig.4 Influence of λ on MCIS positioning accuracy

圖5 λ1對MCIS定位精度的影響Fig.5 Influence of λ1 on MCIS positioning accuracy

4.2.2 抽樣數量對定位精度的影響

圖6為200次蒙特卡洛試驗下,不同樣本數量定位誤差的RMSE,顯然抽樣數量在10～100時始終能夠達到3 km以內的定位誤差,說明樣本數量不是影響定位精度的主要因素。這是因為構造的重要性采樣函數和似然函數相似度很高,生成樣本均服從某一個確定的高斯分布。即使生成不同數量樣本,抽樣值在目標真實位置周圍的分布規律是相同的,所以生成樣本數量對定位結果影響不大。因此,可以選擇較少的樣本達到一個理想的定位結果,進一步減小計算的復雜度,提高定位的收斂速度。

圖6 樣本數量對MCIS定位精度的影響Fig.6 Influence of sample number on MCIS positioning accuracy

4.3 實驗3：不同觀測條件下的定位性能

觀測平臺飛行速度保持300 m/s不變,圖7為不同的相位差變化率測量間隔下的定位相對誤差R隨觀測時長的變化曲線,圖中縱坐標表示相對于觀測平臺和目標歐氏距離的相對誤差。可見,當相位差變化率測量間隔一定,增加觀測時長可以提高定位精度,當觀測時間T在30 s左右時能夠達到R=5%的定位誤差；此外,在相同的觀測時長下,測量間隔越小,也就是數據率越高時,收斂速度越快且定位精度越高。這是因為觀測時間的延長和數據率的提高使得系數矩陣維度增加,會在一定程度上修正IS估計帶來的誤差,但會導致算法的計算量增加的弊端。

圖8驗證了不同運動速率下定位的RMSE,觀測平臺作直線飛行。當平臺速度為100 m/s時收斂速度較慢,觀測至60 s左右時才能達到R=10%的定位精度,觀測平臺速度越大,定位效果越好,尤其是當速度高于400 m/s左右時,定位精度最高,而且觀測至10 s左右時就達到R=10%定位誤差,收斂速度較快。

圖7 不同測量間隔下定位性能Fig.7 Positioning performance at different measurement intervals

圖8 不同平臺速度下的定位性能Fig.8 Positioning performance at different platform speeds

對比分析延長觀測時間和增大觀測平臺的速度對算法的定位性能的影響,兩種情況都會延長觀測平臺的位移。所以,在本文的定位模型中,增大觀測平臺位移能夠有效提高算法的定位性能。這是由于單次測量相位差變化率對應一個定位曲面,多次測量后只有經過目標輻射源的定位曲面才能交于一點,前后觀測平臺位移越大,使定位曲面的差異性越明顯,利于算法在迭代過程中快速排除虛假點,有效提高收斂速度和估計精度,定位結果的可信度更高。但是延長觀測時間與定位收斂速度之間相互沖突,因此可以考慮通過增大飛行速度提升定位性能。

從圖10中可以看出觀測平臺作機動時,MCIS算法的定位性能得到改善,定位精度更高且收斂速度更快。產生這一情況的原因是觀測平臺姿態發生變化的過程中,干涉儀基線矢量方向相對于目標輻射源產生旋轉,這種旋轉產生的“大”相位差變化率明顯改善了定位精度,這足以說明觀測平臺作機動轉彎時可觀測性更強。

圖9 直線運動和機動運動軌跡圖Fig.9 Trajectories of straight and maneuvering motion

圖10 直線運動和機動運動下的定位性能Fig.10 Positioning performance in straight and maneuvering motion

4.4 實驗4：測量噪聲對定位精度的影響

圖11 不同測量噪聲下幾種定位算法的RMSE對比Fig.11 Comparison of RMSE of several positioning algorithms with different measurement noise

4.5 實驗5：算法對初始估計的敏感性分析

圖12 初始估計誤差對算法定位性能的影響Fig.12 Influence of initial estimation error on algorithm positioning performance

4.6 實驗6：算法復雜度分析

本文除了將定位精度和收斂速度作為衡量本文算法的性能標準之外,也對算法的復雜度進行了對比分析。表1對比了本文算法與EKF算法、NLS算法以及網格搜索(GS)的計算復雜度。K表示觀測次數;D表示目標狀態維度;M表示樣本數量;Ngrid表示網格密度。

分析表1中的各算法的計算復雜度,EKF算法和NLS算法復雜度與迭代次數成正比,采用網格搜索(GS)求解最大似然估計的計算復雜度與網格密度密切相關,而且計算量隨估計量維度D的增加呈指數型爆炸增長。這里選擇樣本數量M=20,觀測次數K=100,網格密度Ngrid=25,計算復雜度比,本文提出的MCIS算法計算復雜度上明顯優于其他幾種算法,GS算法計算量最大,耗時最長,可以通過減小網格密度或觀測次數降低其他幾種算法的復雜度,但必然會影響定位精度。

表1 算法復雜度對比Table 1 Complexity assessment of existing algorithms

5 結 論