李 悅,2,柳 征3,伍渝江
(1.蘭州大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,蘭州 730000;2.解放軍61541部隊(duì),北京 100094;3.國防科技大學(xué) 電子科學(xué)與工程學(xué)院,長沙 410073)
三星座時(shí)差定位系統(tǒng)[1]通過三顆星測量同一脈沖的TOA(Time of Arrival),形成兩路時(shí)差曲面,加上地球面的約束,實(shí)現(xiàn)對(duì)地面輻射源的定位。但在實(shí)際中,被偵察輻射源往往在衛(wèi)星過頂時(shí)段關(guān)機(jī),規(guī)避偵察,達(dá)到反偵察的目的。這種情況下,三星座定位系統(tǒng)只能在過頂前、后的一小段時(shí)間捕獲有關(guān)目標(biāo)的少量數(shù)據(jù),這種數(shù)據(jù)被稱為弱觀測數(shù)據(jù),并將此類目標(biāo)源所處的區(qū)域稱為弱觀測區(qū)?,F(xiàn)有三星座定位方法多采用牛頓迭代或解析法[2],利用兩路時(shí)差數(shù)據(jù)和地球模型進(jìn)行定位解算,當(dāng)過頂時(shí)段無法獲得有效數(shù)據(jù)時(shí),現(xiàn)有方法不能充分利用弱觀測數(shù)據(jù),對(duì)弱觀測區(qū)目標(biāo)定位精度較差,并且受目標(biāo)高程誤差影響較大,當(dāng)目標(biāo)沒有先驗(yàn)高程信息支持時(shí),定位精度不理想[2]。
針對(duì)上述問題,本文重點(diǎn)研究弱觀測數(shù)據(jù)定位綜合處理問題,提出一種基于非線性最小二乘(NLS)的時(shí)差定位綜合算法。
如圖1所示,三星座僅在過頂前和過頂后的實(shí)線區(qū)域?qū)δ繕?biāo)輻射源進(jìn)行了有效觀測,而過頂時(shí)的虛線區(qū)域輻射源關(guān)機(jī)規(guī)避。根據(jù)文獻(xiàn)[2]的定位誤差分析可知,采用牛頓迭代法定位在過頂時(shí)段定位效果較好,過頂前、后定位效果相對(duì)較差。因?yàn)槲墨I(xiàn)[2]的定位方程采用了地球模型作為約束,所以高程假設(shè)在很大程度上影響了定位精度。時(shí)差定位綜合算法將過頂前的數(shù)據(jù)同過頂后的數(shù)據(jù)配對(duì)綜合處理,在一次定位運(yùn)算中通過同時(shí)運(yùn)用過頂前、后的兩組觀測數(shù)據(jù),提高定位精度,同時(shí)定位方程不依賴地球模型作為約束,對(duì)近地目標(biāo)定位精度穩(wěn)定,不受高程影響。
圖1 弱觀測數(shù)據(jù)定位示意圖
假設(shè)目標(biāo)真實(shí)位置為X(x,y,z),三星座在過頂前和過頂后的觀測位置分別為
(1)
(2)
對(duì)目標(biāo)進(jìn)行觀測可得到4個(gè)有效時(shí)差值[3-5],據(jù)此可進(jìn)行三維定位[1]。假設(shè)三星座在過頂前位置得到的兩組時(shí)差測量值為t1、t2,在過頂后位置觀測得到的兩組時(shí)差測量值為t3、t4。在沒有測量誤差時(shí),ti與真實(shí)時(shí)差相等,當(dāng)引入誤差后:
cti=fi(X)+ni,i=1,2,3,4
(3)
式中,c為電磁波傳播速度,ni為誤差項(xiàng),為方便將cti記為ri:
ri=fi(X)+ni,i=1,2,3,4
(4)
方程組可以重寫為
(5)
式中 ‖·‖2表示2-范數(shù)。可以將方程(4)寫成一個(gè)向量等式:
r=f(X)+n
(6)
測量誤差n為多元隨機(jī)變量,其協(xié)方差矩陣記為4×4的正定矩陣N:
N=E[(n-E(n))(n-E(n))T]
(7)
因?yàn)镹為對(duì)稱的正定矩陣,所以N有逆。在三星座定位系統(tǒng)中,我們認(rèn)為真實(shí)位置X為未知的非隨機(jī)變量,測量誤差n為零均值高斯噪聲。由文獻(xiàn)[6-8],r關(guān)于X的條件密度函數(shù)為
則式(8)的極大似然估計(jì)為使P值最大的X,令
Q(X)=[r-f(X)]TN-1[r-f(X)]
(9)
可知式(8)的極大似然估計(jì)等價(jià)為Q(X)最小。
事實(shí)上,這種最小化Q(X)的估計(jì)方法對(duì)非高斯噪聲也是合理的,因?yàn)镹-1在方程(9)中僅被視為權(quán)系數(shù)矩陣。
距離方程f(X)為非線性方程組,將f(X)線性化,在迭代點(diǎn)X0附近泰勒展開,保留前兩項(xiàng)有:
f(X)=f(X0)+f′(X0)(X-X0)
(10)
(11)
將式(10)代入式(9)得:
迭代格式為
本節(jié)給出時(shí)差定位綜合算法同文獻(xiàn)[2]中牛頓迭代法的誤差理論分析及GDOP分布。
式(13)經(jīng)過簡單的變換可得:
由式(15)可以看出,估計(jì)誤差受到線性化和噪聲的影響,假設(shè)線性化誤差為零,則對(duì)真值的估計(jì)是無偏的,且誤差為
(16)
令
(17)
則式(16)可以寫成:
(18)
定位誤差協(xié)方差矩陣為
(20)
定位誤差的GDOP分布[9]為
(21)
由文獻(xiàn)[2]可知,牛頓法迭代格式為
Xk+1=Xk-[g′(Xk)]-1g(Xk)
(22)
式中,g(X)為定位方程組,處理方法同3.1,得:
(23)
誤差協(xié)方差矩陣為
定位誤差的GDOP分布為
(25)
本節(jié)給出兩組仿真結(jié)果:4.1節(jié)將時(shí)差定位綜合算法與文獻(xiàn)[2]提出的牛頓迭代法的定位誤差GDOP圖進(jìn)行對(duì)比分析,來驗(yàn)證時(shí)差定位綜合算法對(duì)弱觀測數(shù)據(jù)的處理能力;4.2節(jié)中給出了兩個(gè)相差3 km高程的真實(shí)定位點(diǎn),通過比對(duì)分析,給出高程對(duì)兩種算法定位精度的影響。比對(duì)指標(biāo)主要考慮定位誤差。鑒于三星時(shí)差定位的應(yīng)用背景是對(duì)全球可能輻射源的定位,因此應(yīng)當(dāng)采用全球平均誤差較小的地心大地坐標(biāo)下的地球橢球。目前,以GPS為代表的應(yīng)用中用的都是WGS-84坐標(biāo)系[10],因此,本文仿真中所用地球橢球都采用WGS-84坐標(biāo)系的地球橢球,方程為
(26)
式中a=6 378 137 m,e2=0.006 694 379 990 13。衛(wèi)星的高度約為1 000 km,星間距約為110 km,成三角形分布,過頂前后觀測時(shí)間間隔為700 s。
過頂前衛(wèi)星真實(shí)位置為A1(-2.2552×106,5.7649×106,4.3409×106)m,B1(-2.1702×106,5.7321×106,4.4235×106)m,C1(-2.2578×106,5.706×106,4.4286×106)m。
過頂后衛(wèi)星真實(shí)位置為A2(-4.6821×106,5.9947×106,2.5667×105)m,B2(-4.6411×106,6.0212×106,3.5363×105)m,C2(-4.7180×106,5.9603×106,3.8311×105)m。
在式(14)、(22)的迭代過程中,初始迭代點(diǎn)的選擇會(huì)影響解的精度[11],我們不進(jìn)行討論。在仿真過程中使用衛(wèi)星A觀測時(shí)刻星下點(diǎn)作為初始迭代點(diǎn),在實(shí)際應(yīng)用中可以采用文獻(xiàn)[2]的解析解作為時(shí)差定位綜合算法的初始迭代點(diǎn)。
本小節(jié)給出3種定位方法在不同誤差條件下的定位誤差GDOP分布。其中定位方法1為時(shí)差定位綜合算法;定位方法2為牛頓法,在圖3中我們將牛頓法過頂前后的兩次定位誤差GDOP分布圖合在一起給出;定位方法3是將牛頓法在過頂前、后的兩次定位結(jié)果的算術(shù)平均作為定位結(jié)果。3種方法的誤差GDOP圖如圖2~4所示。
圖2 方法1定位誤差GDOP分布圖
圖3 方法2定位誤差GDOP分布圖
圖4 方法3定位誤差GDOP分布圖
圖2~4給出了3種定位方法對(duì)零高程目標(biāo)在不同的測距誤差和星址測量誤差條件下定位產(chǎn)生的誤差GDOP分布,圖中黑色實(shí)心三角表示三星座在觀測位置的地面投影,三角間的連接虛線為衛(wèi)星航跡的地面投影。弱觀測區(qū)域?yàn)閳D中衛(wèi)星航跡的中垂線附近。當(dāng)α=0 m、δ=24 m,方法1對(duì)弱觀測數(shù)據(jù)的定位誤差主要分布在2 000 m以內(nèi),少量分布在2 000~2 500 m;方法2定位誤差主要分布在2 500~9 500 m;方法3產(chǎn)生的誤差為2 500~3 000 m。當(dāng)α=150 m、δ=24 m時(shí),方法1對(duì)弱觀測數(shù)據(jù)的定位誤差主要分布在7 000~9 000 m,方法2定位誤差主要分布在9 000~45 000 m,方法3產(chǎn)生的誤差為8 000~12 000 m??梢钥闯?,在兩種誤差條件下,方法1對(duì)弱觀測區(qū)目標(biāo)的定位精度要高于方法2和方法3。
為了更清楚地表現(xiàn)3種方法的各自定位優(yōu)勢,我們給出相對(duì)誤差GDOP圖,即將兩種不同方法的誤差值做比值,則GDOP圖的1等高線將定位好壞明顯區(qū)分。如圖5所示,我們僅就α=150 m、δ=24 m一組誤差條件做相對(duì)GDOP對(duì)比。
圖5 相對(duì)定位誤差GDOP分布圖
圖5(a)是方法3與方法2的相對(duì)GDOP圖,小于等高線1的區(qū)域表示方法3定位精度高于方法2,可以看出,方法3對(duì)弱觀測區(qū)目標(biāo)定位要好于方法2。圖5(b)是方法1與方法2的相對(duì)GDOP圖,可以看出,方法1對(duì)弱觀測區(qū)目標(biāo)定位要好于方法2。圖5(c)是方法1與方法3的相對(duì)GDOP圖,表明方法1對(duì)弱觀測區(qū)定位好于方法3。通過對(duì)相對(duì)誤差GDOP圖的分析,顯然,方法1對(duì)弱觀測區(qū)目標(biāo)的定位精度是最高的。
運(yùn)用4.1節(jié)中的3種定位方法對(duì)兩個(gè)有高程差的弱觀測區(qū)目標(biāo)在不同誤差條件下定位,為方便理解,目標(biāo)改用經(jīng)緯高坐標(biāo),WGS-84坐標(biāo)系的坐標(biāo)轉(zhuǎn)換參照文獻(xiàn)[2]。目標(biāo)A(128.8,25.492,0)、B(128.8,25.492,3000),這兩點(diǎn)僅在高程有差異。定位誤差結(jié)果如表1所示。其中方法2有兩組結(jié)果,分別為過頂前、后兩次定位的誤差。
表1 定位誤差對(duì)比Table 1 Comparison of position error m
通過對(duì)表1數(shù)據(jù)縱向比較,顯然方法1對(duì)弱觀測數(shù)據(jù)在兩種誤差條件下定位效果都要好于方法2和方法3。通過橫向數(shù)據(jù)比較,發(fā)現(xiàn)方法1定位精度并不受高程影響,在0高程和3 km高程上,定位誤差基本保持穩(wěn)定。方法2和方法3的定位誤差都隨高程增加而增加,這是由于定位方程依賴地球模型方程,高程對(duì)方法2和方法3定位誤差的影響文獻(xiàn)[2]中有詳細(xì)的討論。當(dāng)觀測目標(biāo)可以確定高程時(shí),例如對(duì)海面艦船目標(biāo)的觀測,認(rèn)為目標(biāo)零高程,可以將方法1的結(jié)果在高程面做投影得到結(jié)果,這樣可以提高定位精度。
基于WGS-84模型的時(shí)差定位方法在電子戰(zhàn)、導(dǎo)航等領(lǐng)域有較為廣泛的應(yīng)用。實(shí)際應(yīng)用中,因?yàn)樾l(wèi)星軌道相對(duì)固定,所以目標(biāo)往往對(duì)衛(wèi)星偵察進(jìn)行規(guī)避,導(dǎo)致衛(wèi)星在過頂時(shí)段不能獲得有效數(shù)據(jù)。因此,研究弱觀測區(qū)的時(shí)差定位方法應(yīng)該是具有實(shí)際意義的。
本文以處理三星座弱觀測區(qū)數(shù)據(jù)為背景,提出了一種基于NLS的時(shí)差定位綜合方法,分析了該方法對(duì)目標(biāo)定位的原理,并對(duì)其定位誤差進(jìn)行了理論推導(dǎo)與數(shù)值分析。通過與其它方法的比較表明,在弱觀測區(qū)該方法能獲得更優(yōu)的定位結(jié)果,同時(shí)該方法不受目標(biāo)高程影響,對(duì)近地目標(biāo)能夠提供穩(wěn)定的定位結(jié)果。
參考文獻(xiàn):
[1] 孫仲康,周一宇,何黎星. 單多基地有源無源定位技術(shù)[M].北京:國防工業(yè)出版社,1996.
SUN Zhong-kang , ZHOU Yi-yu , HE Li-xing. Sigle/Mulit-Bulti-Base Active and Passive Location Technology[M].Beijing: The National Defense Industry Press,1996. (in Chinese)
[2] 鐘丹星. 低軌三星星座測時(shí)差定位方法若干問題研究[D].長沙:國防科學(xué)技術(shù)大學(xué),2002.
ZHONG Dan-xing. Research on some technology problems of Low-Earth-Orbit 3-satellites location system using TDOA measurements[D]. Changsha: National University of Defense Technology, 2002. ( in Chinese)
[3] Schmidt R O. A new approach to geometry of range difference location[J]. IEEE Transactions on Aerospace Electronic Systems,1972, 8(4): 821-835.
[4] Hahn W R,Tretter S A.Optimum processing for delay-vectorestimation in passive signal arrays[J]. IEEE Transactions on Information Theory, 1973, 19(9): 608-614.
[5] Spirito M A.On the accuracy of cellular mobile station location estimation[J]. IEEE Transactions on Vehicle Technology,2001, 50(5): 674-685.
[6] Mensin C,Plass S.Positioning algorithms for cellular networks using TDOA[C]//Proceedings of Internationnl Conferen ce Acoustics Speech Signal Processing (ICASSP).Toulouse,France:[s.n.],2006: 513-516.
[7] So H C,Chan Y T,Chan F K W.Closed-form formulae foroptimum time-difference-of-arrival based localization[J]. IEEE Transactions on Signal Process, 2008, 56(8): 2614-2620.
[8] Torrieri D J.Statistical theory of passive location systems[J]. IEEE Transactions on Aerospace Electronic System,1984, 20(2): 183-197.
[9] Levanon N. Lowest GDO Pin 2-D scenarios[J]. IEE Proceedings of Radar Sonar Navigation, 2000, 147(3): 149-155.
[10] 郗曉寧,王威,高玉東.近地航天器軌道基礎(chǔ)[M].長沙:國防科技大學(xué)出版社,2003.
XI Xiao-ning, WANG Wei,GAO Yu-dong. Fundamentals of Near-earth Spacecraft Orbit [M]. Changsha: National University of Defense Technology Press, 2003 .(in Chinese)
[11] 袁亞湘,孫文瑜.最優(yōu)化理論與方法[M].北京:科學(xué)出版社,1997.
YUAN Ya-xiang,SUN Wen-yu. Optimization Theory and Methods[M]. Beijing: Science Press,1997.(in Chinese)