熊 杰,陳 俊,寧 靜,曹師齊,楊 海,萬洪容
(中國西南電子技術(shù)研究所 成都 610036)
輻射源無源定位是雷達[1]、聲納[2]、無線通信網(wǎng)絡[3]以及電子對抗[4]等領(lǐng)域中的一個基本問題。針對某一外源輻射源,多平臺協(xié)同無源定位通??梢圆捎没诘竭_角(AOA)的定位體制[5-6]、基于到達時間差(TDOA)的定位體制[7],以及基于AOATDOA 聯(lián)合的定位體制[8-9]。當源輻射源運動或者協(xié)同平臺運動時,到達頻率差(FDOA)通常也被加入無源定位算法以提高定位精度[10-12]。另一方面,由于上述幾類觀測量與輻射源位置成非線性關(guān)系,因此輻射源無源定位問題不是一個簡單的問題。近年來涌現(xiàn)了多種算法試圖解決輻射源無源定位問題,如基于泰勒展開的迭代搜索算法[13]、閉合解定位算法[6-12]等。雖然閉合解定位算法的定位精度容易受到傳感器位置誤差的影響,但是由于其具有不需要輻射源位置的先驗知識且回避了迭代搜索算法常遇到的迭代發(fā)散問題的優(yōu)點,因而備受關(guān)注。本文同樣將焦點聚于移動輻射源無源定位閉合解算法。
現(xiàn)有的很多無源定位閉合解算法能直接用于電子對抗領(lǐng)域中針對敵方移動輻射源無源定位的問題。文獻[8]所提的基于AOA-TDOA 的聯(lián)合定位閉合解算法忽視了傳感器位置誤差所帶來的影響;而文獻[9]所提的基于AOA-TDOA 的聯(lián)合定位閉合解算法雖然考慮了傳感器位置誤差的影響,但所提算法只針對雙平臺協(xié)同定位這一特例;文獻[7]所提的基于TDOA 的定位閉合解算法以及文獻[10-12]所提基于TDOA-FDOA 的聯(lián)合定位閉合解算法,雖然都考慮了協(xié)同平臺的位置誤差,但均要求參與定位任務的我方協(xié)同平臺的個數(shù)不小于5 個。這一要求限制了上述幾種算法在實際作戰(zhàn)場景中的應用范圍。針對上述問題,在建立移動輻射源位置與各種觀測量以及各個傳感器位置的近似等式關(guān)系的基礎(chǔ)上,本文提出一種新的AOA-TDOAFDOA 聯(lián)合定位閉合解算法。所提算法在保持閉合解框架不變的前提下,不僅考慮了傳感器位置誤差所帶來的影響,同時還能將參與定位的我方平臺個數(shù)減小至3 個。進一步對所提算法進行理論分析,在觀測噪聲及傳感器位置誤差均為高斯噪聲并且噪聲強度適中時,所提算法的定位誤差均方根能達到克拉美羅(CRLB)界。3 平臺協(xié)同定位仿真實驗結(jié)果以及對照仿真實驗結(jié)果驗證了上述理論結(jié)論。
假定某時刻有M 個傳感器對一移動輻射源開展聯(lián)合協(xié)同無源定位,定位場景如圖1 所示。
記移動輻射源在該時刻的三維位置向量為uo=[xo,yo,zo]T,三維速度向量為記該時刻第 i(i =1,2,···,M)個傳感器的所在真實位置為,真實速度為一般地,定位算法并不能獲取上述傳感器的真實坐標與真實速度,而只能獲取傳感器所在平臺的導航設備提供的位置向量與速度向量。記該時刻導航設備給出的第i個傳感器的所在位置為 si= [xi,yi,zi]T,速度為
為簡便行文,將各傳感器的位置向量與速度向量收集在一起,并記為, 其中若記傳感器位置誤差向量為 ?β ,那么, 其中,?s=一般地, ?β通常認為是零均值的高斯向量,其協(xié)方差矩陣為 Qβ。
設1 號傳感器為參考傳感器,那么該時刻則有M?1 個TDOA 觀測量作為無源定位算法的輸入。 ri1為第 i(i =2,3,···,M)個傳感器與第1 個傳感器之間的TDOA 觀測量。一般地,其中 ni1為 測量噪聲,為第i 個傳感器與第1 個傳感器之間的TDOA 真值,并且有,其中代 表該時刻移動輻射源與第k(k=1,2,···,M)個傳感器之間的歐氏距離。記在該時刻M?1 個TDOA 觀測向量為r = [r21,r31,···,rM1]T,則r = ro+nt,其中為TDOA 真實值向量,nt=[n21,n31,···,nM1]T為對應的觀測噪聲向量。假設噪聲向量 nt為零均值、協(xié)方差矩陣為 Qt的高斯分布向量,并且與傳感器位置誤差向量 ?β、方位角噪聲向量 nθ以 及高低角噪聲向量 n?統(tǒng)計獨立。
同樣地,設1 號傳感器為參考傳感器,那么該時刻則有M?1 個FDOA 觀測量作為無源定位算法的輸入。為 第 i(i =2,3,···,M)個傳感器與第1 個傳感器之間的 FDOA 觀測量。一般地,其中為 測量噪聲,為第i 個傳感器與第1 個傳感器之間的FDOA 真值,并且有記在該時刻M?1 個FDOA 觀測向量為則為FDOA 真實值向量,為對應的觀測噪聲向量。假設噪聲向量 nν為 零均值、協(xié)方差矩陣為 Qν的高斯分布向量,并且與傳感器位置誤差向量 ?β、方位角噪聲向量 nθ、 高低角噪聲向量 n?以及TDOA噪聲向量 nt統(tǒng)計獨立。
需要說明的是,許多文獻(如文獻[12])在定義TDOA 表達式時將光速c 也包含在內(nèi),在定義FDOA 表達式時同時將光速c 和輻射源信號載頻fo也包含在內(nèi)。而這兩類定義在已知輻射源信號載頻 fo時沒有本質(zhì)區(qū)別。因此,本文不特意區(qū)分上述兩類定義,其含義可由上下文得到。
綜上,為簡化算法推導,將上述3 類觀測量匯集 起 來。記 α= [θT,?T,rT,r˙T]T為 該 時 刻 的AOA、TDOA 及FDOA 測量值向量,為對應的觀測噪聲向量。一般地,觀測噪聲向量?α被 建模為零均值、協(xié)方差矩陣 Qα為高斯向量,其中 Qα的定義為:
那么,移動輻射源AOA-TDOA-FDOA 聯(lián)合無源定位算法所涉問題為:已知某一時刻的AOA、TDOA 及FDOA 測量向量 α與傳感器位置向量 β,以及二者對應的噪聲/誤差向量的協(xié)方差矩陣 Qα與Qβ,盡可能準確地估計該時刻的移動輻射源位置uo。
移動輻射源AOA-TDOA-FDOA 聯(lián)合定位算法采用兩階段估計流程。第一階段,首先建立起移動輻射源位置參數(shù) uo、與 觀測向量 α及傳感器位置向量 β的偽線性關(guān)系,進而得到 uo與粗略估計值;第二階段,修正第一階段的粗估計值,從而精確估計出某時刻移動輻射源的位置。
式中, ηθi為測量噪聲 nθi以及傳感器位置誤差 ?si帶來的近似誤差。記 gθi=[?sinθi,cosθi,0]T,那么式(4)可以等價表達為:
當測量噪聲 nθi很弱時, cos (nθi)≈1及 sin(nθi)≈nθi,及。將上述兩近似等式以及帶入式(4),有。再利用圖1 中的一個幾何關(guān)系式并定義,那么有:
式(5)與式(6)分別建立了向量 uo與 方位角 θi之間的近似等式關(guān)系,以及近似誤差 ηθi與測量噪聲nθi及傳感器位置誤差 ?si之間的函數(shù)關(guān)系。
式中, g?i=[?sin?icosθi,?sin?isinθi,cos?i]T; η?i為測量噪聲 n?i以 及傳感器位置誤差 ?si帶來的近似誤差。當測量噪聲 n?i很 弱時,及近似成立,并將上述兩近似關(guān)系及代入式(7),忽略二階誤差項,有:
如文獻[10-12]及文獻[7]所示,TDOA 觀測量 ri1, i(i =2,3,···,M)可以通過如下方式進行線性化處理。首先,對等式的兩邊同時取平方運算,有:
式中, εt,i代表近似誤差,其定義為:
對于FDOA 觀測量 r˙i1, i(i=2,3,···,M)的線性化處理方式為同時對式(10)及式(11)的右端對時間求導數(shù)[10-12],即:
1)第一階段:粗估計。
記 η=[ηθ1,ηθ2,···,ηθM,η?1,η?2,···,η?M]T,εt=,以及,并定義,那么聯(lián)合式(5)~式(8)、式(10)~式(12),有:
式中,
式(14)中的矩陣 Oi×j代表 i×j維的全零矩陣,矩陣hθ?,Gθ?,ht,hν,Gt,Gν分別定義為:
式中, ⊙表示向量Schur 乘法; 1m×1與 0m×1分別表征m 維取值全1 的列向量與全0 的列向量;Gθ=[gθ1,gθ2,···,gθM]T,G?=[g?1,g?2,···,g?M]T;S=[s1,s2,···,sM]T。
解決上述矛盾的一種折中的處理方式是采取兩步加權(quán)最小二乘計算,即:
2)第二階段:精估計。
對于第一階段的參數(shù)估計值 φ1,其前三維分量φ1(1:3)是 移動輻射源在某時刻的位置向量 uo的估計值,其第五至第七維分量 φ1(5:7)是該輻射源在該時刻的瞬時速度向量的估計值。但根據(jù)最小二乘理論,由式(20)估計得出 φ1的 前提條件是與(uo,u˙o)無 關(guān),但這又與第一階段估計參數(shù)的定義矛盾,因而按照這種方式估計出的移動輻射源的位置向量與速度向量不精確。因此,在第二階段,本文使用第一階段的參數(shù)估計值 φ1(4)與 φ1(8)來構(gòu)建另一個最小二乘模型,以提高 uo與估計精度,從而實現(xiàn)對該時刻移動輻射源的無源定位。
記向量 φ1對 應的估計誤差向量為 ?φ1,由于φ1(1:3)=uo+?φ1(1:3),則有:
類似地,有:
以及:
基于上面的討論,第二階段參數(shù)向量可選擇為:
聯(lián)合式(21)~式(24),有:
式中,
另一方面,針對等式誤差 ε2,有:
式中,矩陣 B2,D2的定義如下所示:W2的選取形式[10],即
一旦得到 φ2, 那么位置向量 uo及 速度向量 u˙o的最終估計值可由以下方式計算得到,即:
式中, Π =diag[φ1(1:3)?s1];表示對應元素相除。矩陣 Π的存在避免了對 φ2元素開平方時的正負符號模糊問題。
式(32)最終給出了某時刻移動輻射源AOATDOA-FDOA 的聯(lián)合定位閉合解。下面給出定位解u與對應的協(xié)方差矩陣,若記算法給出的定位解向量,那么對式(25)求一階全微分有:
在第二階段的估計過程中, s1與被視作常量,因而 ?s1與未出現(xiàn)在式(33)里。式(28)表明當噪聲強度不大時, ε2近似滿足零均值,因此式(31)所示的 ?φ2同樣近似為零均值。所以,本文所提的AOA-TDOA-FDOA 聯(lián)合定位給出的閉合解同樣近似滿足零均值特性。那么,由式(33)可得定位解向量 Θ的理論協(xié)方差矩陣為:
需要指出的是,文獻[10-11]分析了基于TDOA與FDOA 觀測量的聯(lián)合定位閉合解算法,其第一階段估計參數(shù)向量與本文選取的向量一致。是8 維向量,因此文獻[10-11]所提算法至少需要5 個協(xié)同平臺參與定位任務,才能按照式(19)及式(20)估計。而在實際的電子對抗作戰(zhàn)過程中難以滿足這一條件,本文在繼承文獻[10-11]所提算法框架的前提下,將AOA 觀測量引入定位算法,這就使得實際作戰(zhàn)所需的協(xié)同平臺個數(shù)降至3 個,既節(jié)約了定位成本,也提高了無源定位閉合形式算法的戰(zhàn)場生存能力。
上節(jié)詳細討論了某時刻移動輻射源AOA-TDOAFDOA 的聯(lián)合定位閉合解,本節(jié)將對上述算法進行理論分析,并給出其估計誤差的理論誤差下界,即CRLB 界。
AOA、TDOA 以及FDOA 測量值向量 α與傳感器位置和速度向量 β均服從高斯分布,并且兩者獨立,因 此 數(shù) 據(jù) 向 量 μ =[αT,βT]T與 參 數(shù) 向 量的概率密度函數(shù)為:
式中,
由向量求導定義可知,式(37)中的偏導?αo/?Θ與? αo/?β 實質(zhì)是向量 αo對 向量u、 向量、 向量 s及向量求偏導,即有:
以及:
下面逐一給出式(38)、(39)中每個偏導矩陣的具體表達式。根據(jù)定義有OM×3以 及=OM×3M,=OM×3M。記λi,κi,(i=1,2,···,M),具有如下表達形式:
那 么,?θo/?u=[λ1,λ2,···,λM]TAu1,??o/?u=[κ1,κ2,···,κM]TAu2,?θo/?s=diag(Au1)As1,??o/?s=diag(Au2)As2。
記 E3M×3=[I3×3,···,I3×3]T,并定義兩個(M?1)×3M 的矩陣 C1與 矩陣 C2,其第 i(i =1,2,···,M ?1)行的定義為:
綜合上述各個偏導矩陣,那么式(38)與式(39)所表征的偏導矩陣的具體表達形式為:
將式(43)代入式(37),即可得6M+6 維方陣CRLB(?)的具體形式。
由矩陣求逆引理可得移動輻射源位置及速度參數(shù) Θ的CRLB 界,即:
式中, X?1表示當測角與時差/頻差傳感器不存在位置誤差及速度誤差時的參數(shù)Θ 的CRLB 界,因此式(44)中的第二項表征了由測角與時差/頻差傳感器位置及速度誤差引起的參數(shù) Θ的CRLB 界的增量部分。式(44)的跡表征了任一種多平臺AOA-TDOA-FDOA聯(lián)合定位算法的估計誤差均方根的下界。
文獻[10]已證明,當TDOA 與FDOA 測量誤差與傳感器位置與速度誤差都服從高斯分布并且噪聲不大,并且傳感器距離輻射源較遠時,定位解向量 Θ的理論協(xié)方差矩陣 cov (Θ)能達到其理論誤差下界。基于類似的證明過程,當無源定位算法加入AOA 觀測量后,本文所提算法所給出定位解向量Θ的理論協(xié)方差矩陣 cov (Θ)也能達到對應的CRLB界。由于篇幅限制,在此省略這一結(jié)論的推導過程,下面將給出仿真實驗,驗證本文提出算法的均方根誤差在噪聲不大時同樣能達到CRLB 界。
例。移 動 輻 射 源 初 始 位 置 為[57 891, ?12 536,125]T(單位為m,下同),初始速度為 [12 ,30,0]T(單位為m/s,下同);協(xié)同平臺主站初始位置為[5 000,1 500,3 500]T,初始速度為 [12 0,?20,0]T;協(xié)同平臺輔站1 的初始位置為 [21 500,15 000,5 800]T,初始速度為 [13 0,0,0]T;協(xié)同平臺輔站2 的初始位置為[20 500,?45 000,5 800]T, 初 始 速 度 為 [1 10,50,0]T。輻射源載波頻率為1.1 25×109Hz。為評價無源定位算法的定位精度,本文采用均方根誤差為評價標準,其定義為:
為驗證本文提出的移動輻射源AOA-TDOAFDOA 聯(lián)合定位閉合解算法的性能,本節(jié)首先給出一組3 平臺協(xié)同定位一個移動輻射源位置的仿真實
式中, MS E為某時刻無源定位所估計的輻射源位置的均方根誤差; uo為該時刻移動輻射源真實位置;ui為第i 次獨立仿真實驗所估計的輻射源位置;N 為獨立仿真實驗次數(shù)。需要指出的是,在N 次獨立仿真實驗過程中,每次仿真中輻射源移動軌跡以及每個平臺的真實運動軌跡都一致。本節(jié)后續(xù)部分將首先以此為基礎(chǔ)給出幾組仿真實例,驗證本文所提算法在各種噪聲強度下的定位能力;然后再給出兩組對比實驗,說明本文所提算法在定位精度與定位時間等方面較其他已有定位算法所具備的優(yōu)勢。
為驗證本文所提算法對平臺位置誤差的適應能力,在這組仿真中將3 架空中平臺的圖1 中x、y、z 向上的位置誤差標準差變化范圍同時設置為10?2~105m,所有平臺圖1 中x、y、z 向上的速度誤差標準差固定為0.2 m/s,并且將所有AOA 測量噪聲標準差均固定為3°,所有TDOA 測量噪聲標準差均固定為100 ns,以及所有FDOA 測量噪聲標準差均固定為100 Hz。5 000 次獨立仿真形成的均方根定位誤差曲線如圖2 所示。為直觀說明本文所提算法對平臺位置誤差的適應能力,每種設置下的定位誤差理論下界(即式(44)所示的 CR LB(Θ)矩陣的前3 個對角元素和的平方根)也同時繪制在圖2 中。
圖2表明當所有協(xié)同平臺在圖1 中x、y、z這3 個方向的位置誤差標準差均不超過5×103m時,本文所提算法在3 平臺協(xié)同定位場合下的均方根定位誤差能達到對應的CRLB 界。
為驗證本文所提算法對測角噪聲強度的適應能力,在這組仿真中,將每個平臺的方位角測量噪聲標準差變化范圍同時設置為10?2°~102°,而所有平臺在圖1 中x、y、z 向的位置誤差標準差都固定為5 m、速度誤差標準差都固定為0.2 m/s,并且將所有高低角測量噪聲標準差均固定為3°,所有TDOA測量噪聲標準差均固定為100 ns,所有FDOA 測量噪聲標準差均固定為100 Hz。5 000 次獨立仿真形成的均方根定位誤差曲線如圖3 所示。同時為直觀說明本文所提算法對測角噪聲強度的適應能力,每種設置下的CRLB 界也同時繪制在圖3 中。
圖3表明當所有協(xié)同平臺所測得的方位角的測量誤差標準差不大于20°時,本文所提算法在3 平臺協(xié)同定位場合下的均方根定位誤差能達到對應的CRLB 界。
為驗證本文所提算法對TDOA 噪聲強度的適應能力,在這組仿真中,將每個TDOA 測量噪聲標準差變化范圍同時設置為10?2~106ns,所有平臺在x、y、z 向的位置誤差標準差都固定為5 m,速度誤差標準差都固定為0.2 m/s,并且將所有方位/高低角測量噪聲標準差均固定為3°,所有FDOA測量噪聲標準差均固定為100 Hz。5 000 次獨立仿真形成的均方根定位誤差曲線如圖4 所示。同時為直觀說明本文所提算法對TDOA 測量噪聲強度的適應能力,每種設置下的CRLB 界也同時繪制在圖4 中。
圖4表明當所有TDOA 觀測量的測量噪聲標準差不大于105ns 時,本文所提算法在3 平臺協(xié)同定位場合下的均方根定位誤差能達到對應的CRLB 界。
為驗證本文所提算法對FDOA 噪聲強度的適應能力,在這組仿真中,將每個FDOA 測量噪聲標準差變化范圍同時設置為10?2~104Hz,所有平臺在x、y、z 向的位置誤差標準差都固定為5 m,速度誤差標準差都固定為0.2 m/s,并且將所有方位/高低角測量噪聲標準差均固定為3°,所有TDOA測量噪聲標準差均固定為100 ns。5 000 次獨立仿真形成的均方根定位誤差曲線如圖5 所示。同時為直觀說明本文所提算法對FDOA 測量噪聲強度的適應能力,每種設置下的CRLB 界也同時繪制在圖5 中。
圖5表明當所有FDOA 觀測量的測量噪聲標準差不大于200 Hz 時,本文所提算法在3 平臺協(xié)同定位場合下的均方根定位誤差能達到對應的CRLB 界。
針對移動輻射源TDOA-FDOA 聯(lián)合無源定位閉合解算法需要至少5 個協(xié)同平臺參與一次定位任務的問題,本文將AOA 觀測量引入,提出一種新的基于AOA-TDOA-FDOA 的聯(lián)合定位閉合解算法。本文所提算法在兼顧平臺位置不確定性以及保持代數(shù)閉合解框架不變的前提下,能將參與定位的平臺個數(shù)有效地減少至3 個。同時,理論分析與仿真實驗結(jié)果都表明所提算法在面對工程上常見的噪聲強度時,定位性能可以達到理論誤差下界。這表明,本文所提算法既能拓寬算法的戰(zhàn)場應用范圍,也能提高算法的戰(zhàn)場生存能力。