喻晨龍, 譚賢四, 曲智國, 謝 非, 李陸軍
(1. 空軍預(yù)警學(xué)院研究生管理大隊,湖北 武漢 430019; 2. 空軍預(yù)警學(xué)院防空預(yù)警裝備系, 湖北 武漢430019;3. 長江勘測規(guī)劃設(shè)計研究院, 湖北 武漢 430019)
臨近空間升力式再入飛行器的飛行高度高、飛行速度快、突防能力強、打擊距離遠,具有升力式氣動布局,其由火箭助推至高空后,依靠地球重力再入大氣層,可長距離、無動力滑翔飛行,典型的有HTV-2和AHW兩種[1-3]。這類再入式飛行器在氣動外形、運行環(huán)境、飛行受力、控制方式和機動模式等方面與彈道導(dǎo)彈存在明顯差異,給預(yù)警探測系統(tǒng)帶來了極大的挑戰(zhàn)[4-8]。針對這一殺手锏武器,目前可行的探測手段是以反導(dǎo)預(yù)警系統(tǒng)為依托,以米波雷達、凝視雷達等新型裝備為補充進行協(xié)同探測,從而盡早盡快發(fā)現(xiàn)目標(biāo)并維持較高的跟蹤精度,為后端的攔截系統(tǒng)提供連續(xù)穩(wěn)定的信息支撐。
通常,遠程預(yù)警相控陣?yán)走_靠近邊境線前沿部署,探測距離可達5 000 km,然而受地球曲率的限制,在探測臨近空間目標(biāo)時探測距離最遠只有1 300 km。臨近空間升力式再入飛行器在縱向上存在跳躍機動,在橫向上存在側(cè)向偏移,其中:只有當(dāng)目標(biāo)處于雷達探測范圍內(nèi)時,才能保持連續(xù)穩(wěn)定跟蹤;當(dāng)目標(biāo)飛出探測范圍,即向下機動至地平線以下時會使情報中斷。若能利用目標(biāo)在雷達水平面波束層中飛行時獲得的有限觀測信息,對目標(biāo)的運動包絡(luò)和可能的跳躍空域進行預(yù)測,則可輔助雷達在遠程預(yù)警相控陣?yán)走_丟失目標(biāo)后補充照射,從而保證情報的連續(xù)。
升力式再入飛行器的彈道包括平衡滑翔和跳躍滑翔[9]。研究者們針對目標(biāo)的運動參數(shù)、攻角和傾側(cè)角的變化規(guī)律進行了深入分析,如:李廣華等[10]認(rèn)為,在平衡滑翔狀態(tài)下目標(biāo)的速度與高度、速度與速度傾角之間存在解析關(guān)系,且隨著傾側(cè)角的變化而變化;文獻[11-12]的作者通過分析氣動力加速度對目標(biāo)質(zhì)心運動和繞質(zhì)心運動的影響,說明氣動力改變是目標(biāo)跳躍機動的主要原因;李廣華[13]認(rèn)為目標(biāo)的運動狀態(tài)參數(shù)反映著目標(biāo)的氣動力,只有當(dāng)參數(shù)滿足特定數(shù)學(xué)關(guān)系時,才能維持平衡滑翔軌跡,當(dāng)偏離平衡點時就會產(chǎn)生跳躍軌跡;王路等[14]通過目標(biāo)的三自由度彈道方程計算了目標(biāo)的攻角和最大升阻比。然而,上述研究未利用運動參數(shù)、攻角和傾側(cè)角的變化規(guī)律進行彈道預(yù)測。鑒于此,筆者在分析目標(biāo)的運動特性的基礎(chǔ)上,估算目標(biāo)的運動參數(shù)、攻角和傾側(cè)角,合理假設(shè)攻角和傾側(cè)角的變化規(guī)律,最后通過數(shù)值仿真方法預(yù)測目標(biāo)飛離雷達探測范圍后可能的跳躍空域。
為了描述目標(biāo)的運動過程和雷達的探測過程,定義地心(Earth Canter,EC)坐標(biāo)系、飛行器東北天(Vehicle East-North-Up,VENU)坐標(biāo)系和雷達站東北天(Radar East-North-Up, RENU)坐標(biāo)系3個相關(guān)的參考坐標(biāo)系[14]。這3個坐標(biāo)系各自包括1個直角坐標(biāo)系和1個球坐標(biāo)系。
假設(shè)地球為均勻圓球,不考慮自轉(zhuǎn)影響,則飛行器在無動力滑翔過程中滿足EC球坐標(biāo)系和VENU直角坐標(biāo)系的三自由度運動方程[8,14],即
(1)
(2)
式中:r(t)、χ(t)、φ(t)分別為飛行器在EC球坐標(biāo)系下的地心距、經(jīng)度和緯度;v(t)、γ(t)、ψ(t)分別為飛行器在VENU直角坐標(biāo)系下的速度、航跡傾角和航向角;m為飛行器的質(zhì)量;g為重力加速度;σ(t)為傾側(cè)角;D(t)、L(t)分別為飛行過程中所受阻力和升力,
(3)
(4)
其中CD(t)為阻力系數(shù),CL(t)為升力系數(shù),ρ(t)為大氣密度,S為飛行器氣動參考面積。
文獻[15-18]作者認(rèn)為:當(dāng)飛行器的速度>8 Ma時,升力系數(shù)可簡化為攻角α(t)的線性函數(shù),阻力系數(shù)可簡化為α(t)的二次函數(shù)。對氣動參數(shù)插值擬合后可得
CD(t)=0.044+0.08α(t)+2.2α2(t),
(5)
CL(t)=-0.151 7+3.151α(t)。
(6)
大氣密度滿足指數(shù)模型,即
ρ(t)=ρ0exp[-(r(t)-Re)/hs],
(7)
式中:ρ0=1.752 kg/m3,為參考大氣密度;hs=6.7 km,為參考高度;Re=6 371 km,為地球半徑。
飛行器穿越水平波束屏的過程中,雷達每間隔一段時間Ts(即搜索間隔時間)照射一次目標(biāo),得到了N組觀測值{z1,z2,…,zN},這些觀測值是在RENU球坐標(biāo)系下獲取的,記為
zi=[?r(ti)θr(ti)φr(ti)]T,
(8)
式中:?r(ti)、θr(ti)、φr(ti) 分別為第i次照射時的目標(biāo)徑向距離、方位角、俯仰角,其中ti為第i次照射時對應(yīng)的時刻。將觀測值在RENU直角坐標(biāo)系中表示,有
(9)
通過式(9)位置差分可得到速度信息Vr(ti)=[vrx(ti)vry(ti)vrz(ti)]T。
防御方預(yù)警系統(tǒng)可直接獲取的信息僅包括目標(biāo)在穿越水平波束屏過程中的觀測信息,即RENU直角坐標(biāo)系下的位置和速度。而目標(biāo)的面質(zhì)比、控制量,以及在EC球坐標(biāo)系下的地理位置和VENU直角坐標(biāo)系下的速度是未知的。
由式(1)-(7)可知:目標(biāo)的彈道軌跡取決于運動狀態(tài)初值和攻角/傾側(cè)角的變化規(guī)律。為了實現(xiàn)軌跡預(yù)測,首先把觀測信息由RENU直角坐標(biāo)系轉(zhuǎn)換到EC球坐標(biāo)系和VENU直角坐標(biāo)系,得到目標(biāo)飛離水平波束屏的運動狀態(tài),然后結(jié)合式(1)-(7)估算目標(biāo)的攻角/傾側(cè)角和面質(zhì)比。
1) 將位置觀測信息從RENU直角坐標(biāo)系轉(zhuǎn)換到EC球坐標(biāo)系。首先從RENU直角坐標(biāo)系轉(zhuǎn)換到EC直角坐標(biāo)系,然后從EC直角坐標(biāo)系轉(zhuǎn)換到EC球坐標(biāo)系。
假設(shè)雷達站OR所在地的經(jīng)度為J,緯度為W,海拔高度為H。要從RENU直角坐標(biāo)系變換到EC直角坐標(biāo)系,根據(jù)歐拉角關(guān)系,首先把RENU直角坐標(biāo)系原點先沿地心矢移動Re+H到地心O,然后繞x軸順時針旋轉(zhuǎn)π/2-W,最后繞z軸順時針旋轉(zhuǎn)π/2+J。即
(10)
式中:Xr(ti)=[xr(ti)yr(ti)zr(ti)]T,為RNEU直角坐標(biāo)系下的位置;Xe(ti)=[xe(ti)ye(ti)ze(ti)]T,為EC直角坐標(biāo)系下的位置;
(11)
為從RENU直角坐標(biāo)系到EC直角坐標(biāo)系的坐標(biāo)變換矩陣,其中Tx、Tz分別為繞x、z軸的基元旋轉(zhuǎn)矩陣。
將位置觀測信息從EC直角坐標(biāo)系轉(zhuǎn)換到EC球坐標(biāo)系,即
(12)
假設(shè)在第k次照射后飛行器飛離水平波束屏,雷達測量得到目標(biāo)的質(zhì)心位置為[?r(tk)θr(tk)φr(tk)]T(其中tk為第k次照射時對應(yīng)的時刻),這是一個RENU球坐標(biāo),通過坐標(biāo)變換得到此時質(zhì)心在EC球坐標(biāo)系下位置為[re(tk)χe(tk)φe(tk)]T。
2) 將速度信息從RENU直角坐標(biāo)系轉(zhuǎn)換到VENU直角坐標(biāo)系。先考慮位置坐標(biāo)的變換,將位置坐標(biāo)從RENU直角坐標(biāo)系轉(zhuǎn)換到EC直角坐標(biāo)系,再從EC直角坐標(biāo)系轉(zhuǎn)換到VENU直角坐標(biāo)系。由于從EC直角坐標(biāo)系轉(zhuǎn)換到VENU直角坐標(biāo)系,與從RENU直角坐標(biāo)系轉(zhuǎn)換到EC直角坐標(biāo)系的變換關(guān)系剛好相反,因此可類比式(10)中的方法轉(zhuǎn)換。特別需注意的是,此時坐標(biāo)系的原點變?yōu)榱孙w行器的質(zhì)心。
要從EC直角坐標(biāo)系轉(zhuǎn)換到VENU直角坐標(biāo)系,根據(jù)歐拉角關(guān)系,首先把EC直角坐標(biāo)系先繞z軸逆時針旋轉(zhuǎn)π/2+χe(tk),然后繞x軸逆時針旋轉(zhuǎn)π/2-φe(tk),最后坐標(biāo)原點由地心O沿地心矢移動re(tk)到飛行器的質(zhì)心。即
(13)
式中:Xe(ti)=[xe(ti)ye(ti)ze(ti)]T和Xv(ti)=[xv(ti)yv(ti)zv(ti)]T分別表示第i次照射時目標(biāo)在EC直角坐標(biāo)系和VENU直角坐標(biāo)系下的位置;
(14)
為從EC直角坐標(biāo)系到VENU直角坐標(biāo)系的坐標(biāo)變換矩陣。
那么,從RENU直角坐標(biāo)系到VENU直角坐標(biāo)系的位置坐標(biāo)變換關(guān)系為
(15)
因此,從RENU直角坐標(biāo)系下的速度矢量[vrx(ti)vry(ti)vrz(ti)]T轉(zhuǎn)換到VENU直角坐標(biāo)系下的速度矢量[vvx(ti)vvy(ti)vvz(ti)]T,其變換關(guān)系為
(16)
由于VNEU直角坐標(biāo)系下的速度、航跡傾角與航向角有如下關(guān)系:
(17)
因此,根據(jù)RENU直角坐標(biāo)系下的一組量測信息差分得到飛行器飛離水平波束屏的速度[vrx(tk)vry(tk)vrz(tk)]T,然后根據(jù)坐標(biāo)轉(zhuǎn)換關(guān)系把各軸上的分速度變換到VENU直角坐標(biāo)系下,最后得到飛行器的速度vv(tk)、航跡傾角γv(tk)和航向角ψv(tk)。
式(1)、(2)描述了目標(biāo)的地理位置和速度的迭代關(guān)系,式(1)-(7)隱含了目標(biāo)控制量的解算方法。因此,將式(2)改寫為[20]
(18)
(19)
(20)
式中:ρ(ti)=ρ0exp[-(re(ti)-Re)/hs]。
(21)
(22)
(23)
式中:K(ti)為目標(biāo)的升阻比。
將式(21)代入式(5)、(6),可以解算出目標(biāo)第i次照射時的攻角α(ti)、升力系數(shù)CL(ti)和阻力系數(shù)CD(ti)。
當(dāng)初值確定后,飛行器的彈道由攻角和傾側(cè)角決定,攻角和傾側(cè)角并不是隨機變化的。其中:1) 攻角直接決定飛行過程中的氣動特性,再入過程中若調(diào)整頻繁且幅值較大,將嚴(yán)重影響飛行器的熱防護系統(tǒng),因此攻角的變化要盡量簡單,一般選擇分段線性函數(shù),大小通常控制為0°~ 20°,此外攻角與升阻比緊密相關(guān),最常用的制導(dǎo)律是以最大升阻比飛行;2) 傾側(cè)角是主控制量,其確定升力在縱向和側(cè)向上的分量,并決定飛行器橫向機動程度,改變傾側(cè)角是為了規(guī)避禁飛區(qū)和突防,傾側(cè)角最常用的制導(dǎo)律是方位角誤差走廊,其幅值不變、符號翻轉(zhuǎn),幅值通??刂圃凇?0°以內(nèi)。
因此,升力式再入飛行器在遠端運動時,為了維持最小能耗和最大升阻比飛行,并且實現(xiàn)有效突防,對攻角和傾側(cè)角的變化規(guī)律可做如下假設(shè):
(24)
(25)
式中:α(tk)為飛行器飛出水平玻束屏?xí)r估計的攻角;αmaxK(tk)為最大升阻比對應(yīng)的攻角;σ(tk)為飛行器飛出水平玻束屏?xí)r估計的傾側(cè)角;t為飛行時間;T為攻角改變的時間;T1為傾側(cè)角改變的時間。
將飛行器飛離水平波束屏?xí)r的運動參數(shù)、攻角和傾側(cè)角的預(yù)測值代入式(1)-(7)進行數(shù)值仿真,隨著迭代的進行,當(dāng)式(2)中第1次出現(xiàn)γ(t)=0時停止。
假設(shè)某型遠程預(yù)警相控陣?yán)走_部署地點在東經(jīng)110°、北緯10°、海拔高度0 km,波束寬度1.2°,沿著正東設(shè)置單層水平波束屏,量測噪聲服從高斯分布,測量誤差??=100 m,?θ=?φ=0.001 5°,搜索間隔時間Ts=1 s。飛行器穿越水平波束屏的過程中雷達在RENU球坐標(biāo)系內(nèi)得到了14組觀測數(shù)據(jù),如表1所示。
表1 RENU球坐標(biāo)系內(nèi)觀測數(shù)據(jù)
按照2.1節(jié)的估算方法對觀測數(shù)據(jù)進行坐標(biāo)轉(zhuǎn)換,計算得到目標(biāo)飛離水平波束屏?xí)r的地理位置和速度分別為:re(tk)=6 425 km,χe(tk)=117.472 8°,φe(tk)=9.889 5°;vv(tk)=19.850 9 Ma,γv(tk)=-1.868 6°,ψv(tk)=92.029 4°。進一步,根據(jù)2.2節(jié)的估算方法,計算得到目標(biāo)的面質(zhì)比為S/m=0.002 5,在飛離水平波束屏?xí)rα(tk)≈5°,σ(tk)≈40°。
飛行器飛離水平波束屏后,其攻角和傾側(cè)角按照式(24)、(25)的規(guī)律變化,即攻角和傾側(cè)角先是維持不變,再分別在時刻T和T1發(fā)生改變,攻角由飛離時的攻角α(tk) 變化為最大升阻比攻角αmaxK(tk) ,傾側(cè)角幅值不變、符號翻轉(zhuǎn)。這時,隨著迭代的進行,軌跡發(fā)生較大變化。由于攻角和傾側(cè)角改變的時刻是無法預(yù)知的,當(dāng)式(2)中第1次出現(xiàn)γ(t)=0時迭代停止,此時的軌跡曲線就圍成了一個空域范圍Ω,在EC球坐標(biāo)系下的目標(biāo)跳躍空域范圍如圖1所示。
臨近空間武器的興起將打破常規(guī)作戰(zhàn)規(guī)則,在探測這類目標(biāo)時,遠程預(yù)警雷達因受到地球地平線遮擋而丟失目標(biāo),本文基于目標(biāo)彈道特性,利用有限的觀測數(shù)據(jù)對目標(biāo)可能的跳躍空域進行預(yù)判,能提前指示預(yù)警系統(tǒng)調(diào)度雷達資源探測,具有較強的軍事應(yīng)用價值。
參考文獻:
[1] PAUL L M,VINCENT L R,LUAT T N,et al.NASA hypersonic fight demonstration overview,status,and future plans[J].Acta astronautica,2004,55(7):619-630.
[2] 李益翔.美國高超聲速飛行器發(fā)展歷程研究[D].哈爾濱:哈爾濱工業(yè)大學(xué),2016:18-32.
[3] STEVEN H W,JEFFREY S,DALE S,et al.The DAPPA/AF falcon program:the hypersonic technology vehicle #2(HTV-2) flight
demonstration phase[C]∥15th AIAA International Space Planes and Hypersonic Systems and Technologies Conference.28 April-1 2008,Dayton,Ohio:2008-2539.
[4] RANDALL T V,LAERENCE D H,CHARLES R,M,et al.X-43A Hypersonic vehicle technology development[J].Acta Astronautica,2006(59):181-191.
[5] STEVEN W,MING T,SUE M,et al.Falcon HTV-3X-A reusable hypersonic test bed[C]∥15th AIAA International Space Planes and Hypersonic Systems and Technologies Conference.28 April-1 May,2008,Dayton,Ohio:2008-2544.
[6] 秦雷,李君龍,周荻.國外臨近空間高速滑翔彈頭運動特征研究[J].戰(zhàn)術(shù)導(dǎo)彈技術(shù),2015(1):16-20.
[7] 李君龍,秦雷,謝曉瑛.臨近空間非彈道式機動模式與跟蹤濾波問題[J].現(xiàn)代防御技術(shù),2017,45(1):99-106.
[8] 李慧峰.高超聲速飛行器制導(dǎo)與控制技術(shù)[M].北京:中國宇航出版社,2012:10.
[9] 陳小慶,候中喜,劉建霞.高超聲速滑翔飛行器彈道特性分析[J].導(dǎo)彈與航天運載技術(shù),2011(2):5-9.
[10] 李廣華,張洪波,湯建國.高超聲速滑翔飛行器典型彈道特性分析[J].宇航學(xué)報,2015,36(4):397-403.
[11] OLIVIER D M,ROBERT H B.Tracking and identification of a maneuvering reentry vehicle[C]∥AIAA Guidance,Navigation,and Control Conference and Exhibit.11-14 August 2003,Austin,Texas,USA,2003:2003-5446.
[12] 張凱,熊家軍,韓春耀,等.一種基于氣動力模型的高超聲速滑翔目標(biāo)跟蹤算法[J].宇航學(xué)報,2017,38(2):123-130.
[13] 李廣華.高超聲速滑翔飛行器運動特性分析及彈道跟蹤預(yù)報方法研究[D].長沙:國防科學(xué)技術(shù)大學(xué),2016:17-49.
[14] 王路,邢清華,毛藝帆.基于升阻比變化規(guī)律的再入高超聲速滑翔飛行器軌跡預(yù)測算法[J].系統(tǒng)工程與電子技術(shù),2015,37(10):2335-2339.
[15] PHILIPS T H.A Common Aero Vehicle(CAV) model,description,and employment guide[R].USA:Schafer Corporation,2003.
[16] 趙漢元.飛行器再入動力學(xué)與制導(dǎo)[M].長沙:國防科技大學(xué)出版社,1997:48-97.
[17] 黃長強,國海峰,丁達理.高超聲速滑翔飛行器軌跡優(yōu)化與制導(dǎo)綜述[J].宇航學(xué)報,2014,35(4):369-379.
[18] 王璐璐,秦玉亮,王宏強,等.飛行器升阻比估計誤差分析方法研究[J].雷達科學(xué)與技術(shù),2012,10(2):175-179.