郭 躍,劉新學(xué),蔣 鳴
(1.第二炮兵工程大學(xué) 初級(jí)指揮學(xué)院,西安 710025;2.第二炮兵裝備研究院,北京 100089)
彈道導(dǎo)彈在主動(dòng)段不存在機(jī)動(dòng)變軌等突防措施,而且彈道變化單一[1],因此攔截處于主動(dòng)段的彈道導(dǎo)彈的作戰(zhàn)效果十分明顯。對(duì)主動(dòng)段彈道導(dǎo)彈的準(zhǔn)確跟蹤和定位是攔截主動(dòng)段彈道導(dǎo)彈的重要一環(huán)[2-3]。申鎮(zhèn)[4]等基于運(yùn)動(dòng)方程對(duì)單星無(wú)源的主動(dòng)段彈道射向進(jìn)行了估算,張濤、安瑋[5-7]等給出了基于推力加速度模板的主動(dòng)段跟蹤算法,但是需要獲得主動(dòng)段的推力加速度信息。陳映和程臻等[8-9]給出了一種參變量助推彈道導(dǎo)彈時(shí)不變運(yùn)動(dòng)模型,但也是基于情報(bào)數(shù)據(jù)庫(kù)為先驗(yàn)知識(shí),同時(shí)還認(rèn)為導(dǎo)彈主動(dòng)段的運(yùn)動(dòng)模型可以用加加速度模型來(lái)近似計(jì)算,并采用交互多模型(interacting multiple model,IMM)的不敏卡爾曼濾波(IMM-UKF)算法進(jìn)行了仿真。Li X R和Jilkov V P給出了全彈道的跟蹤模型和濾波模型,將主動(dòng)段的跟蹤模型簡(jiǎn)化為變加速度模型或重力轉(zhuǎn)彎模型,但與導(dǎo)彈主動(dòng)段實(shí)際的運(yùn)動(dòng)模型相比還存在較大差異[10-12]。本文從彈道導(dǎo)彈主動(dòng)段的實(shí)際運(yùn)動(dòng)模型入手,建立了一種基于在線修正的主動(dòng)段彈道跟蹤模型,該模型可以對(duì)導(dǎo)彈參數(shù)進(jìn)行在線估計(jì),也可用于離線彈道導(dǎo)彈參數(shù)的采集。
設(shè)地基雷達(dá)坐標(biāo)系定義與NUE(North-Up- East)坐標(biāo)系相同,坐標(biāo)原點(diǎn)經(jīng)、緯度及當(dāng)?shù)馗叱虨?φr0,λr0,hr0),則從地基雷達(dá)坐標(biāo)系到地心直角坐標(biāo)系(earth right angle coordinate system)的轉(zhuǎn)換矩陣為
(1)
如圖1所示,圖中,Fp為導(dǎo)彈推力,θtra為推力方向與軌跡坐標(biāo)系x軸的夾角,FR為總空氣動(dòng)力,θR為總空氣動(dòng)力與軌跡坐標(biāo)系x軸的夾角。AT為射面與北向的夾角,θg為導(dǎo)彈在雷達(dá)坐標(biāo)系中的矢徑與其在OradarXradarYradar平面內(nèi)投影的夾角。在主動(dòng)段的跟蹤過(guò)程中應(yīng)選擇軌跡坐標(biāo)系來(lái)描述,其原點(diǎn)Otraj為導(dǎo)彈質(zhì)心,OtrajXtraj軸與導(dǎo)彈速度方向一致,OtrajYtraj位于射擊平面內(nèi)且垂直于OtrajXtraj,OtrajZtraj符合右手定則。
因此軌跡坐標(biāo)系到雷達(dá)坐標(biāo)系的轉(zhuǎn)換矩陣為
設(shè)t時(shí)刻導(dǎo)彈在雷達(dá)坐標(biāo)系的位置為(xrmyrmzrm)T,則導(dǎo)彈所受引力加速度為[1]
式中:grx,gry,grz分別為重力在雷達(dá)坐標(biāo)系內(nèi)的各軸的分量;Gm為地心引力常數(shù);rm為導(dǎo)彈在地心直角坐標(biāo)系內(nèi)的距離;ω為地球自轉(zhuǎn)角速度;B0為發(fā)射點(diǎn)大地緯度;φs為導(dǎo)彈在(xrmyrmzrm)T所對(duì)應(yīng)的地心緯度;R0,rx,R0,ry,R0,rz分別為雷達(dá)坐標(biāo)系原點(diǎn)的地心矢徑在雷達(dá)坐標(biāo)系各軸的分量;μ為地球扁率系數(shù);a,e1,αe分別為橢球地球體長(zhǎng)半軸,偏心率和地球扁率;(ωrxωryωrz)T為地球自轉(zhuǎn)角速度在雷達(dá)坐標(biāo)系內(nèi)的向量表示。
科氏加速度arc在雷達(dá)坐標(biāo)系的投影為
牽連加速度are在雷達(dá)坐標(biāo)系的投影為
(7)
式中:vrmx,vrmy,vrmz分別為導(dǎo)彈速度在雷達(dá)坐標(biāo)系內(nèi)各軸向的分量;ωrx,ωry,ωrz分別為地球自轉(zhuǎn)角速度在雷達(dá)坐標(biāo)系內(nèi)各軸向的分量。在彈道導(dǎo)彈目標(biāo)的跟蹤過(guò)程中,由科氏力和牽連力引起的位置跟蹤誤差會(huì)隨著時(shí)間不斷累積,對(duì)洲際彈道導(dǎo)彈(ICBM)而言主動(dòng)段飛行時(shí)間在幾百s以上,因此在對(duì)主動(dòng)段飛行的導(dǎo)彈建模時(shí)應(yīng)考慮科氏力和牽連力對(duì)導(dǎo)彈加速度的影響。
導(dǎo)彈主動(dòng)段氣動(dòng)力的計(jì)算由于涉及導(dǎo)彈姿態(tài)、氣動(dòng)力系數(shù)以及飛行程序角等,因此,在計(jì)算導(dǎo)彈氣動(dòng)力時(shí)所建立的六自由度空氣動(dòng)力模型不僅參數(shù)多、模型復(fù)雜,而且估算也很難準(zhǔn)確。
設(shè)i為導(dǎo)彈的級(jí)數(shù),一般i小于3級(jí),導(dǎo)彈的實(shí)時(shí)質(zhì)量可表示為
(9)
推力FP在主動(dòng)段對(duì)導(dǎo)彈作用力最大,因此可以通過(guò)加速度大小的變化準(zhǔn)確推算出導(dǎo)彈發(fā)動(dòng)機(jī)是否工作,這樣對(duì)導(dǎo)彈當(dāng)前質(zhì)量m(tc)的估計(jì)就變?yōu)榕c導(dǎo)彈前一時(shí)刻的質(zhì)量m(tc-1)及常數(shù)C有關(guān)。
主動(dòng)段導(dǎo)彈的推力公式為
FP(t)=FP0,i-Se(p0-pH)=FP0,i-ΔFP
(11)
式中:FP0,i為第i級(jí)的地面額定推力,Se為主動(dòng)段發(fā)動(dòng)機(jī)橫截面積,p0為地面的大氣壓,pH為高度H的大氣壓。多級(jí)火箭各級(jí)推力不同,主動(dòng)段導(dǎo)彈推力的大小可以視為由地面額定推力FP0和由大氣壓強(qiáng)引起的推力變化ΔFP共同決定。
推力加速度矢量表示為
arP=FP(t)/m(t)
(12)
在射面內(nèi)推力FP(t)可以視為
式中:θP為推力FP(t)與速度方向的夾角,ΔFP為推力的變化量。因此推力FP(t)的估算就變?yōu)閷?duì)FP0,θP,ΔFP的估算。
導(dǎo)彈在主動(dòng)段的加速度方程可以表示為
am=arP+ara+arc+are+gr
(13)
式中的加速度gr,are,arc可以根據(jù)上述模型準(zhǔn)確地計(jì)算,推力引起的加速度arP,空氣動(dòng)力引起的加速度ara是預(yù)測(cè)的難點(diǎn)。
1)加加速度模型(Jerk)[2-3,7]。
陳映等通過(guò)研究多種不同類型彈道導(dǎo)彈的加速度隨時(shí)間的變化曲線,認(rèn)為在主動(dòng)段導(dǎo)彈的加速度近似為線性,而將加加速度近似為常數(shù)。主動(dòng)段的狀態(tài)向量可以表示為
主動(dòng)段導(dǎo)彈跟蹤的Jerk模型為
Xr(k+1)=AXr(k)+W(k)
(15)
式中:T為探測(cè)周期,A為轉(zhuǎn)換矩陣,W(k)為噪聲干擾。
2)當(dāng)前統(tǒng)計(jì)模型(CSM)。
導(dǎo)彈在主動(dòng)段所受重力、科氏力、牽連慣性力可以通過(guò)計(jì)算獲得,而推力和空氣動(dòng)力很難有準(zhǔn)確的模型進(jìn)行跟蹤,因此對(duì)推力和空氣動(dòng)力引起的加速度可以視為修正的瑞利馬爾科夫過(guò)程。選擇系統(tǒng)的狀態(tài)變量:
Xr=(xryrzrvrxvryvrzarxaryarz)T
(17)
W(t)=(O1×6ωPax(t)ωPay(t)ωPaz(t))T
(21)
IMM算法是一種動(dòng)態(tài)多模型算法,通過(guò)設(shè)計(jì)和建立多種彈道導(dǎo)彈模型的主動(dòng)段跟蹤模型來(lái)構(gòu)造IMM算法是提高導(dǎo)彈跟蹤精度的有效手段[8,11]。本文根據(jù)彈道導(dǎo)彈主動(dòng)段的運(yùn)動(dòng)特點(diǎn)選擇了CS模型和Jerk模型作為跟蹤模型,同時(shí)結(jié)合不敏卡爾曼濾波(unscented Kalmar filter,UKF)和不敏粒子濾波(unscented particle filter,UPF)算法構(gòu)建了IMM-UKF-UPF算法。算法流程圖如圖2所示。
圖2 IMM-UKF-UPF算法流程圖
CS和Jerk跟蹤模型同前,UKF和UPF算法的實(shí)現(xiàn),限于篇幅,不做贅述。
基于彈道模板的算法一方面可以通過(guò)對(duì)彈道的量測(cè)獲得彈道相關(guān)參數(shù),另一方面在實(shí)際的測(cè)量過(guò)程中可以為彈道預(yù)測(cè)提供參考依據(jù)。算法如圖3所示。圖中,Isp為導(dǎo)彈發(fā)動(dòng)機(jī)的比沖。
圖3 基于彈道模板的算法結(jié)構(gòu)圖
對(duì)彈道采用四階Runge-Kutta法進(jìn)行解算。
設(shè)k時(shí)刻導(dǎo)彈加速度濾波后的矢量為
寫出k和k+1時(shí)刻的向量展開式:
以2級(jí)洲際導(dǎo)彈TD2為模型進(jìn)行算法的仿真驗(yàn)證[12],其在雷達(dá)坐標(biāo)系內(nèi)各參數(shù)的變化情況如圖4~圖7和表1所示。
圖4 主動(dòng)段在全彈道中的示意圖
圖5 雷達(dá)探測(cè)的主動(dòng)段各參數(shù)隨時(shí)間的變化
圖6 解析方法對(duì)主動(dòng)段各力的估算結(jié)果分析
圖7 空氣動(dòng)力系數(shù)的估算情況分析
表1 基于彈道模板的各常值參數(shù)的估算
把估算參數(shù)得出的在線修正的彈道跟蹤模型(ONPT)與交互多模型方法(IMM)分別結(jié)合UPF和UKF算法進(jìn)行濾波,仿真結(jié)果如圖8、圖9所示。
圖8 OLPT與IMM算法跟蹤導(dǎo)彈位置均方差比較圖
圖9 OLPT與IMM算法跟蹤導(dǎo)彈速度均方差比較
從仿真結(jié)果來(lái)看,基于在線修正彈道模板的主動(dòng)段彈道跟蹤方法能夠?qū)椀缹?dǎo)彈的各常值參數(shù)進(jìn)行準(zhǔn)確的估計(jì),同時(shí)對(duì)動(dòng)態(tài)變化的參數(shù)也能進(jìn)行比較好的跟蹤,但是導(dǎo)彈級(jí)間轉(zhuǎn)換時(shí)(如圖6)對(duì)推力和發(fā)動(dòng)機(jī)質(zhì)量流量的跟蹤稍顯滯后。該方法可以用于平時(shí)彈道導(dǎo)彈數(shù)據(jù)的采集和分析并構(gòu)建準(zhǔn)確的彈道模板,也可以對(duì)導(dǎo)彈進(jìn)行在線跟蹤。在對(duì)導(dǎo)彈進(jìn)行跟蹤的過(guò)程中,基于在線修正的彈道模板(ONPT)方法與交互多模型(IMM)方法相比有著更好的跟蹤精度,顯然導(dǎo)彈模板的參數(shù)越準(zhǔn)確跟蹤效果越好。
[1] 張毅,楊輝耀,李俊莉.彈道導(dǎo)彈彈道學(xué)[M].長(zhǎng)沙:國(guó)防科大出版社,1999:111-157.
ZHANG Yi,YANG Hui-yao,LI Jun-li.The ballistic trajectory[M].Changsha:National University of Defense Technology Press,1999:111-157.(in Chinese).
[2] 郭尊華,謝維信.彈道導(dǎo)彈跟蹤技術(shù)進(jìn)展[J].信號(hào)處理,2009,25(8a):578-581.
GUO Zun-hua,XIE Wei-xin.Ballistic missile tracking:a review[J].Signal Processing,2009,25(8a):578-581.(in Chinese)
[3] 周宏仁,敬忠良,王培德.機(jī)動(dòng)目標(biāo)跟蹤[M].北京:國(guó)防工業(yè)出版社,1991.
ZHOU Hong-ren,JING Zhong-liang,WANG Pei-de.Mobile target tracking[M].Beijing:National Defense Industry Press,1991.(in Chinese)
[4] 申鎮(zhèn),強(qiáng)勝,易東云.基于運(yùn)動(dòng)方程的單星無(wú)源主動(dòng)段射向估計(jì)方法[J].彈道學(xué)報(bào),2010,22(1):7-10.
SHEN Zhen,QIANG Sheng,YI Dong-yun.Course head estimation method based on movement equation by single satellites passive detection in boost phase[J].Journal of Ballistics,2010,22(1):7-10.(in Chinese)
[5] 張濤,安瑋,周一宇,等.基于推力加速度模板的主動(dòng)段彈道跟蹤方法[J].宇航學(xué)報(bào),2006,27(3):385-389.
ZHANG Tao,AN Wei,ZHOU Yi-yu,et al.Periodic relative motion condition for satellites formations considering nonlinearity[J].Journal of Astronautics,2006,27(3):385-389.(in Chinese)
[6] 張濤,安瑋,周一宇.基于UKF的主動(dòng)段彈道跟蹤算法[J].彈道學(xué)報(bào),2006,18(2):15-18.
ZHANG Tao,AN Wei,ZHOU Yi-yu.Trajectory tracking method in boost phase by using UKF[J].Journal of Ballistics,2006,18(2):15-18.(in Chinese)
[7] 張濤,安瑋,周一宇.主動(dòng)段彈道定位與跟蹤算法[J].彈道學(xué)報(bào),2005,17(4):11-16.
ZHANG Tao,AN Wei,ZHOU Yi-yu.The trajectory locating and tracking algorithm in boost phase[J].Journal of Ballistics,2005,17(4):11-16.(in Chinese)
[8] 陳映,程臻,文樹梁.彈道導(dǎo)彈助推段同時(shí)跟蹤和類型識(shí)別算法研究[J].信號(hào)處理,2011,27(5):749-754.
CHEN Ying,CHENG Zhen,WEN Shu-liang,Study on method for simultaneously tracking and classifying ballistic missile in boost and post-boost phase[J].Signal Processing,2011,27(5):749-754.(in Chinese)
[9] 陳映,文樹梁,程臻.一種適用于助推段彈道導(dǎo)彈的跟蹤方法研究[J].系統(tǒng)仿真學(xué)報(bào),2012,24(5):1 063-1 067.
CHEN Ying,WEN Shu-liang,CHENG Zhen.Method for tracking ballistic missile on boost phase[J].Journal of System Simulation,2012,24(5):1 063-1 067.(in Chinese).
[10] LI X R,JILKOV V P.A survey of maneuvering target tracking,part Ⅱ:ballistic target models[J].SPIE,2001,4 473:559-581.
[11] LI X R,JILKOV V P.A survey of maneuvering target tracking,part Ⅲ:measurements models[J].SPIE,2001,4 473:423-446.
[12] JOHN A,LUKACS I V.Hit-to-kill guidance for the interception of ballistic in the boost phase[D].Monterey,CA:Naval Postgraduate Shool,2006.