葉 飚,曾占魁、,馮 剛,曹喜濱
(1.上海宇航系統(tǒng)工程研究所,上海 201109;2.哈爾濱工業(yè)大學(xué) 航天學(xué)院,黑龍江 哈爾濱 15001)
天基預(yù)警衛(wèi)星由高軌部分的(DSP)和低軌部分的STSS組成。DSP衛(wèi)星主要通過(guò)短波紅外探測(cè)器掃描探測(cè)彈道導(dǎo)彈助推段的高熱尾焰,先確定威脅導(dǎo)彈的初始彈道信息,之后將相關(guān)目標(biāo)信息傳遞至STSS,引導(dǎo)STSS衛(wèi)星控制長(zhǎng)波紅外探測(cè)器由平時(shí)的掃描待命狀態(tài)進(jìn)入捕獲跟蹤狀態(tài),實(shí)現(xiàn)對(duì)自由飛行段彈頭的凝視跟蹤與彈道重構(gòu)。
根據(jù)STSS星座部署的構(gòu)型分析,若STSS星座采用24星的構(gòu)型配置,臨邊探測(cè)模式下可保證軌道高度230~1 500km范圍內(nèi)的全球單重覆蓋,并保證赤緯大于30°的二重以上覆蓋(兩顆星同時(shí)可見(jiàn)),兩極上空的預(yù)警覆蓋達(dá)到極限,最多可有6~7顆的可見(jiàn)星。根據(jù)我國(guó)實(shí)際,多數(shù)對(duì)我國(guó)構(gòu)成威脅的洲際導(dǎo)彈從北極上空經(jīng)過(guò),彈道的大部分均在超過(guò)2顆以上預(yù)警衛(wèi)星的覆蓋范圍內(nèi)。通常情況下,對(duì)多于2顆預(yù)警衛(wèi)星的角度測(cè)量信息,可用最小二乘估計(jì)技術(shù)進(jìn)行處理,獲得更高精度的彈道重構(gòu)信息。但單顆預(yù)警衛(wèi)星的跟蹤目標(biāo)數(shù)有限,過(guò)多占用不必要的監(jiān)視資源,會(huì)削弱整個(gè)系統(tǒng)對(duì)多個(gè)威脅目標(biāo)的預(yù)警能力。因此,在保證彈道重構(gòu)精度條件下,盡可能降低所占用的預(yù)警衛(wèi)星數(shù)。因此,當(dāng)觀測(cè)星數(shù)多于2顆時(shí),需引入選星算法,選擇兩顆觀測(cè)幾何條件最好的預(yù)警衛(wèi)星執(zhí)行跟蹤測(cè)量操作。為此,本文對(duì)基于低軌預(yù)警衛(wèi)星測(cè)量數(shù)據(jù)的彈道重構(gòu)與選星算法進(jìn)行了研究。
J2000地心赤道慣性坐標(biāo)系O-XYZ:坐標(biāo)原點(diǎn)O為地心,基本平面為J2000地球平赤道面,X軸在基本平面內(nèi)指向J2000平春分點(diǎn),Z軸垂直于基本平面,由地球中心指向地球平赤道自轉(zhuǎn)北極方向,X、Y、Z軸構(gòu)成右手系。
傳感器測(cè)量坐標(biāo)系O′-UEN:以傳感器中心位置O′為坐標(biāo)原點(diǎn),地心O與O′的連線為U軸,N軸位于地球自轉(zhuǎn)軸與U軸組成的平面內(nèi),指向正北方向?yàn)檎珽軸依右手法則確定。
O′-UEN,O-XYZ系的關(guān)系如圖1所示。

圖1 坐標(biāo)系示意Fig.1 Coordination system
雙星交匯定位幾何關(guān)系如圖2所示。根據(jù)文中的參數(shù)定義,設(shè)預(yù)警衛(wèi)星S1,S2在J2000坐標(biāo)系中的位置分別為r1,r2。對(duì)彈道導(dǎo)彈(ICBM)探測(cè),每顆預(yù)警衛(wèi)星有兩個(gè)測(cè)量參數(shù):O′-UEN系中星載探測(cè)傳感器測(cè)得視線矢量s1,s2的方位角α和高低角δ。
在O′-UEN系中,由星載探測(cè)傳感器可測(cè)得視線矢量s1,s2的單位矢量分別為

式中:δi,αi分別為si在O′-UEN系中的高低角及方位角(i=1,2)。
假設(shè)ICBM在J2000坐標(biāo)系中的位置為rm。用O′-UEN系至J2000慣性系的變換矩陣,可將導(dǎo)彈在J2000坐標(biāo)系中的坐標(biāo)用預(yù)警衛(wèi)星測(cè)得的坐標(biāo)表示。

圖2 雙星交匯定位幾何關(guān)系Fig.2 Geometry of solving ICBM position using two visible stars observation
設(shè)d1,d2為u1,u2轉(zhuǎn)換至J2000坐標(biāo)系的分量,即

由正弦定理可得ICBM與預(yù)警衛(wèi)星S1間的距離

式中:s為兩顆預(yù)警衛(wèi)星之間的相對(duì)距離;θ為兩觀測(cè)視線矢量交匯夾角,且θ=π-β1-β2。此處:

由此,確定ICBM在J2000坐標(biāo)系中的位置為

在ICBM自由飛行段,飛行空氣動(dòng)力和其他攝動(dòng)力可忽略,考慮地球非球型引力帶諧系數(shù)J2的動(dòng)力學(xué)方程精度即可。取狀態(tài)變量為ICBM在J2000坐標(biāo)系中的位置/速度

則系統(tǒng)狀態(tài)方程可寫(xiě)為

式中:W為系統(tǒng)狀態(tài)噪聲向量,其協(xié)方差矩陣為Qk;g(r)為地球引力加速度矢量,且

觀測(cè)量取雙星幾何定位算法確定的ICBM位置,根據(jù)式(8),有

式中:k為觀測(cè)序列編號(hào);v為測(cè)量噪聲,是零均值的高斯隨機(jī)過(guò)程。
線性化動(dòng)力學(xué)方程和觀測(cè)方程。式(10)對(duì)狀態(tài)向量求偏導(dǎo)

式中:I3×3為單位陣。
離散化處理可得

式中:T為濾波采樣周期;n為泰勒展開(kāi)級(jí)數(shù)。
觀測(cè)方程線性化后有

式中:Hk=[I3×303×3]。
EKF濾波方程為

用低軌預(yù)警衛(wèi)星觀測(cè)數(shù)據(jù)幾何定位算法,先給出3個(gè)不同時(shí)刻t1,t2,t3的ICBM 位置矢量r1,r2,r3,用Gibbs方法確定速度矢量2,有

其中:τ1,τ2,τ3分別為時(shí)刻t1,t2,t3的時(shí)間間隔,且τ1=t1-t2;τ3=t3-t2;τ13=t1-t3[1]。
則,由式(21)~(23)可算出時(shí)刻t2的速度矢量后,綜合時(shí)刻t2的位置矢量r2,就確定出了ICBM在時(shí)刻t2的軌道狀態(tài)。
因低軌預(yù)警衛(wèi)星的常規(guī)姿態(tài)為對(duì)地定向,故只有當(dāng)目標(biāo)與預(yù)警衛(wèi)星、地球間的幾何關(guān)系滿足紅外臨邊探測(cè)的原則時(shí),星上的長(zhǎng)波紅外探測(cè)設(shè)備才可能探測(cè)到目標(biāo)。但完全按可見(jiàn)原則對(duì)24顆低軌預(yù)警衛(wèi)星逐一篩選的效率低,因此有必要建立快速判斷可見(jiàn)預(yù)警衛(wèi)星的算法。
設(shè)當(dāng)前時(shí)刻可見(jiàn)的預(yù)警衛(wèi)星數(shù)目為n,其中的一顆預(yù)警衛(wèi)星Sj與彈道導(dǎo)彈幾何觀測(cè)方程式(12)移項(xiàng)后可得


在式(24)中,rm的三個(gè)分量與ρj為4個(gè)待定量,其余為已知量或觀測(cè)量。對(duì)n顆可見(jiàn)預(yù)警星,可得n組式(24)的方程,則3×n個(gè)方程求解3+n個(gè)未知數(shù),寫(xiě)成觀測(cè)方程的形式有



當(dāng)同時(shí)可見(jiàn)2顆以上預(yù)警衛(wèi)星時(shí),用最小二乘解算多星的觀測(cè)數(shù)據(jù)可得較高的目標(biāo)定位精度。但因一顆預(yù)警衛(wèi)星無(wú)法同時(shí)對(duì)兩個(gè)及以上的目標(biāo)進(jìn)行凝視跟蹤觀測(cè),故對(duì)同一目標(biāo)占用過(guò)多的預(yù)警衛(wèi)星將降低對(duì)其他目標(biāo)的發(fā)現(xiàn)與跟蹤能力。
在雙星定位模式下,不考慮預(yù)警衛(wèi)星的位置誤差,根據(jù)式(24),可寫(xiě)出以彈道導(dǎo)彈位置矢量rm為狀態(tài)變量,觀測(cè)視線的高低角δ和方位角α為觀測(cè)量的觀測(cè)方程。在雙星幾何解算時(shí),觀測(cè)方程為非線性方程,經(jīng)線性化后可得

式中:ΔX為待求狀態(tài)偏差量,且ΔX=Δrm;ΔZ為觀測(cè)偏差量,且

推導(dǎo)可得

因在選星前不可能獲得觀測(cè)信息δ,α,以及彈道導(dǎo)彈rm。因此,可用先驗(yàn)估計(jì)的彈道位置矢量替代式(27)~(31)中的rm,并由計(jì)算出估計(jì)的替代δ,α。
當(dāng)可見(jiàn)預(yù)警星大于2顆時(shí),可用最小二乘法求解式(27),得

本文假設(shè)觀測(cè)量δ,α的誤差ΔZ為零均值的聯(lián)合高斯分布。當(dāng)幾何布局固定時(shí),由ΔZ與ΔX的函數(shù)關(guān)系可得ΔX亦為零均值的高斯分布,計(jì)算ΔX(ΔX)T的期望值可得ΔX的協(xié)方差。則有

式(33)利用了矩陣(HTH)-1為對(duì)稱矩陣的特性。一般情況下,假定觀測(cè)角誤差Δδ,Δα是分布相同且相互獨(dú)立的,其方差等于(σang)2,則

當(dāng)可見(jiàn)兩顆預(yù)警星時(shí),n=2。將式(34)代入式(33),可得

ΔX的3個(gè)分量分別對(duì)應(yīng)的是估計(jì)彈道導(dǎo)彈在ECI坐標(biāo)系中的三方向的位置誤差,其協(xié)方差陣可展開(kāi)為

定義彈道導(dǎo)彈位置估計(jì)的精度因子

(HTH)-1分量寫(xiě)作

由于定位解算的結(jié)果將作為下一步卡爾曼濾波觀測(cè)量,定位誤差方差矩陣可由式(35)給出,提取方差陣covΔX的跡就是幾何定位解算的誤差估計(jì)均值,并作為卡爾曼濾波的觀測(cè)噪聲方差陣引入濾波器,有

設(shè)模擬的ICBM自由飛行段彈道歷元秒0時(shí)刻的軌道根數(shù)為:半長(zhǎng)軸5 462 756.648m,偏心率0.45,傾角124.108°,升交點(diǎn)赤經(jīng)9.96°,近地點(diǎn)幅角265.995°,平近點(diǎn)角98.431°;低軌預(yù)警衛(wèi)星按24顆星三軌道面部署,軌道高度1 600km,傾角102°,預(yù)警衛(wèi)星的軌道位置三方向誤差由均值0、均方根50m的高斯白噪聲;δ,α的測(cè)量誤差由均值0、均方根0.01°的高斯白噪聲;卡爾曼濾波的濾波采樣周期10s,預(yù)警衛(wèi)星多星最小二乘幾何定位的位置矢量作為濾波觀測(cè)量和濾波初值的位置分量,由前30s三次定位解算位置分量的Gibbs差分得到初值的速度分量。仿真結(jié)果如圖3~8所示。

圖3 ICBM自由飛行段可見(jiàn)STSS預(yù)警衛(wèi)星數(shù)Fig.3 Number of visible STSS satellites on ICBM passive phase

圖4 未選星時(shí)雙星幾何解算定位誤差Fig.4 Position solving error using two stars observation without selection

圖5 噪聲方差最小原則選星后雙星定位誤差Fig.5 Position solving error using two stars observation with selection of least covariance noise

圖6 噪聲方差最小原則選星后雙星定軌EKF濾波的位置誤差Fig.6 Position filtering error using two stars observation with selection of least covariance noise

圖7 噪聲方差最小原則選星后雙星定軌EKF濾波的速度誤差Fig.7 Velocity Filtering error using two stars observation with selection of least covariance noise

圖8 觀測(cè)噪聲方差最小原則選星后雙星定軌EKF濾波位置與速度誤差均方根Fig.8 Root mean square of the position and velocity dealing with EKF on the least covariance noise principle selection using two stars observations
本文對(duì)STSS低軌道預(yù)警衛(wèi)星雙星定軌的選星算法進(jìn)行了研究。結(jié)果發(fā)現(xiàn):按24顆星部署的STSS星座可滿足對(duì)洲際彈道導(dǎo)彈彈道重構(gòu)的最小數(shù)量可見(jiàn)星要求;用EKF濾波解算的ICBM彈道重構(gòu)的位置精度優(yōu)于100m(3σ),速度精度優(yōu)于0.15m/s(3σ);在可見(jiàn)預(yù)警衛(wèi)星多于2顆的條件下,用最小二乘解算可提高初步的彈道定位的精度,但將增加系統(tǒng)占用,削弱對(duì)整個(gè)預(yù)警系統(tǒng)對(duì)多威脅目標(biāo)的跟蹤預(yù)警能力;基于觀測(cè)噪聲估計(jì)方差最小的原則選星后,可在保證雙星定軌精度的條件下,最小限度占用監(jiān)視預(yù)警衛(wèi)星資源。
[1] 王 洋,曲長(zhǎng)文,蔣 波.低軌星座對(duì)自由段空間目標(biāo)跟蹤算法研究[J].電子測(cè)量技術(shù),2009,32(5):157-160.
[2] 秦永元,張洪鉞,汪叔華.卡爾曼濾波與組合導(dǎo)航原理[M].西安:西北工業(yè)大學(xué),1998.
[3] BATTIN R H.An introduction to the mathematics methods of astrodynamics[M].Washinton D C.AIAA,1999.