無再入觀測(cè)彈道導(dǎo)彈氣動(dòng)參數(shù)和落點(diǎn)聯(lián)合預(yù)報(bào)

吳 楠,王 鋒,孟凡坤

(信息工程大學(xué) 數(shù)據(jù)與目標(biāo)工程學(xué)院,河南 鄭州 450001)

彈道導(dǎo)彈落點(diǎn)預(yù)報(bào)是導(dǎo)彈預(yù)警系統(tǒng)中的一個(gè)重要環(huán)節(jié),落點(diǎn)預(yù)報(bào)精度是預(yù)警雷達(dá)的重要技術(shù)指標(biāo),提高落點(diǎn)預(yù)報(bào)和誤差分析精度具有很強(qiáng)的工程意義[1-2]。提高落點(diǎn)預(yù)報(bào)精度的關(guān)鍵是建立準(zhǔn)確的預(yù)報(bào)模型[3],如孫瑜等[4]建立了考慮地球J2攝動(dòng)的彈道預(yù)報(bào)模型,潘曉剛等[5]建立了考慮J4攝動(dòng)的彈道預(yù)報(bào)模型,但這2種方法均沒有考慮導(dǎo)彈再入時(shí)所受氣動(dòng)阻力的影響。由于再入氣動(dòng)阻力是影響落點(diǎn)預(yù)報(bào)精度僅次于地球自轉(zhuǎn)的因素[6],較多的研究方法是在再入過程中通過引入觀測(cè)數(shù)據(jù)快速地辨識(shí)氣動(dòng)參數(shù),然后再用于落點(diǎn)預(yù)報(bào),如陳映等[7]利用交互多模型方法估計(jì)彈道系數(shù),史金光等[8]建立了氣動(dòng)復(fù)合系數(shù)模型并借助擴(kuò)展卡爾曼濾波進(jìn)行辨識(shí),趙捍東等[9]借助觀測(cè)數(shù)據(jù)用徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)逼近彈道方程,張金鳳等[10]將再入狀態(tài)估計(jì)和彈道系數(shù)辨識(shí)進(jìn)行聯(lián)合優(yōu)化,均進(jìn)一步提高了落點(diǎn)預(yù)報(bào)精度。然而在導(dǎo)彈的早期預(yù)警階段,由于導(dǎo)彈并未再入,無法獲得再入觀測(cè)數(shù)據(jù),且氣動(dòng)參數(shù)信息也無法預(yù)先準(zhǔn)確獲知,此時(shí)如何考慮再入氣動(dòng)的影響,并建立高精度的落點(diǎn)預(yù)報(bào)和分析模型是研究的難點(diǎn)。

本文在無再入觀測(cè)數(shù)據(jù)條件下,基于彈道導(dǎo)彈打擊目標(biāo)時(shí)對(duì)終端速度變化范圍的限制,反推出氣動(dòng)阻力參數(shù)的變化范圍,然后將氣動(dòng)阻力參數(shù)作為隨機(jī)過程建模,并增廣為目標(biāo)的狀態(tài),建立考慮氣動(dòng)阻力、地球自轉(zhuǎn)和J2攝動(dòng)的高精度落點(diǎn)預(yù)報(bào)方程,采用數(shù)值積分和無跡變換的方法對(duì)落點(diǎn)經(jīng)緯度和落地時(shí)刻進(jìn)行預(yù)報(bào)和誤差分析。

1 考慮氣動(dòng)阻力的導(dǎo)彈落點(diǎn)計(jì)算方程

(1)

(2)

(3)

(4)

2 氣動(dòng)阻力參數(shù)變化范圍預(yù)估

(5)

(6)

(7)

(8)

式中:δβmax(k)為相對(duì)于βmax(k)的小量,可取為δβmax(k)=10-8βmax(k)。

3 基于無跡變換的落點(diǎn)均值預(yù)報(bào)和誤差傳播計(jì)算

氣動(dòng)阻力參數(shù)是導(dǎo)彈的固有屬性,具有確定性,但防御方無法獲取該參數(shù)準(zhǔn)確的變化規(guī)律,因此在彈道和落點(diǎn)預(yù)報(bào)時(shí),只能將氣動(dòng)力參數(shù)作為隨機(jī)量處理。假設(shè)用Wiener過程描述氣動(dòng)力參數(shù)[11],相應(yīng)的氣動(dòng)阻力參數(shù)模型可表示為

(9)

式中:w(t)為白噪聲過程。將式(9)離散化,得:

β(k+1)=β(k)+wk

(10)

式中:wk為離散化驅(qū)動(dòng)白噪聲序列,方差取決于隨機(jī)序列{βk}的分布模型。由于沒有任何氣動(dòng)阻力參數(shù)的先驗(yàn)信息,且為了表征氣動(dòng)力參數(shù)的劇烈變化情況,假設(shè)氣動(dòng)阻力參數(shù)在變化范圍內(nèi)服從均勻分布,則氣動(dòng)阻力參數(shù)的均值為

(11)

方差為

(12)

(13)

當(dāng)系統(tǒng)中含有隨機(jī)量時(shí),對(duì)其統(tǒng)計(jì)特性的分析方法有Monte Carlo抽樣統(tǒng)計(jì)法和協(xié)方差分析描述函數(shù)法等,其中前者具有推導(dǎo)計(jì)算簡(jiǎn)便、結(jié)果直觀可靠等優(yōu)點(diǎn),但缺點(diǎn)是耗時(shí)過長(zhǎng)且結(jié)果不是精確解,后者則具有省時(shí)、高精度且便于理論分析等優(yōu)點(diǎn),但缺點(diǎn)是統(tǒng)計(jì)線性化過程過于繁瑣,且要求系統(tǒng)模型全域一階可微,當(dāng)系統(tǒng)模型較為復(fù)雜且不滿足可微性條件時(shí),協(xié)方差分析描述函數(shù)法將不再適用。無跡變換(unscented transform,UT)是用于計(jì)算經(jīng)過非線性變換的隨機(jī)變量統(tǒng)計(jì)特性的新方法[12],UT不需要對(duì)非線性方程線性化,而是對(duì)狀態(tài)矢量的概率密度函數(shù)進(jìn)行近似,因此將UT應(yīng)用于非線性隨機(jī)系統(tǒng)的狀態(tài)預(yù)報(bào),與協(xié)方差分析描述函數(shù)法相比具有2個(gè)優(yōu)點(diǎn):①不需要求解Jacobian矩陣,不要求系統(tǒng)模型具有一階可微性,算法更易實(shí)現(xiàn),應(yīng)用范圍進(jìn)一步擴(kuò)大;②預(yù)報(bào)獲得的后驗(yàn)均值和方差都能夠精確到二階,具有更高的預(yù)報(bào)精度。本節(jié)采用UT方法,基于再入段增廣狀態(tài)方程式(13),對(duì)考慮再入氣動(dòng)阻力的落點(diǎn)進(jìn)行預(yù)報(bào)和誤差傳播分析。定義落點(diǎn)信息矢量z=(tfLfBf)T,利用式(13)和式(2)～式(4)可以獲得落點(diǎn)信息矢量和初始增廣狀態(tài)矢量間的函數(shù)關(guān)系H(·):

(14)

(15)

式中:O為零矩陣。利用UT進(jìn)行均值和誤差傳播計(jì)算,首先基于式(15)構(gòu)造2n+1個(gè)Sigma點(diǎn)和相應(yīng)權(quán)重,對(duì)于本文情況,y0的維數(shù)n=7。

(16)

z(i)=H(χ(i))i=0,1,…,2n

(17)

最后利用統(tǒng)計(jì)的方法計(jì)算落點(diǎn)信息矢量的均值和協(xié)方差矩陣:

(18)

4 落點(diǎn)誤差橢圓計(jì)算

(19)

然后將計(jì)算得出的誤差橢圓旋轉(zhuǎn)平移至大地經(jīng)緯度平面。

(20)

式中:φ∈[0,2π),l為標(biāo)準(zhǔn)差的倍數(shù),通常取l=3(落入該誤差橢圓的概率為97.3%)。

5 仿真分析

5.1 仿真條件設(shè)置

基于某雷達(dá)針對(duì)某型彈道導(dǎo)彈的探測(cè)數(shù)據(jù),進(jìn)行濾波和彈道估計(jì),以濾波穩(wěn)定時(shí)刻t0=540 s時(shí)的彈道估計(jì)值作為預(yù)報(bào)初值,此時(shí)時(shí)間為2011年3月7日8時(shí)53分19秒(協(xié)調(diào)世界時(shí)),狀態(tài)矢量(地固系位置和速度矢量,單位分別為m和m/s)估計(jì)均值為

對(duì)應(yīng)的估計(jì)誤差協(xié)方差矩陣Px0為

表1 導(dǎo)彈的真實(shí)落點(diǎn)信息

5.2 未知?dú)鈩?dòng)阻力參數(shù)預(yù)估

圖1 氣動(dòng)阻力參數(shù)預(yù)測(cè)均值與真實(shí)值

從結(jié)果可以看出,再入過程中真實(shí)氣動(dòng)阻力參數(shù)是時(shí)變的,變化范圍為6.2×10-5～6.8×10-5m2/kg,但變化幅度不是很大,且由于再入過程持續(xù)時(shí)間較短,所以將其作為常值進(jìn)行近似處理是較為合理的,而本文基于先驗(yàn)信息預(yù)測(cè)獲得的氣動(dòng)阻力參數(shù)均值與真實(shí)值相比,處于同一數(shù)量級(jí),且相對(duì)誤差約為10%,說明本文算法通過先驗(yàn)信息預(yù)測(cè)的氣動(dòng)阻力參數(shù)是較為準(zhǔn)確而有效的。

5.3 落點(diǎn)預(yù)報(bào)誤差分析

表2 3種方法落點(diǎn)預(yù)報(bào)結(jié)果比較

從計(jì)算結(jié)果可以看出:①解析法雖然耗時(shí)最短,約為其他2種算法的1/10,但是無論是落點(diǎn)還是時(shí)間預(yù)報(bào)精度都較差;②積分法落點(diǎn)預(yù)報(bào)精度較高,但由于其時(shí)間預(yù)報(bào)誤差較大,造成時(shí)空不匹配,大大增加了落點(diǎn)預(yù)報(bào)誤差;③本文算法耗時(shí)與積分法相當(dāng),由于考慮了氣動(dòng)阻力的影響而且對(duì)阻力參數(shù)進(jìn)行了有效的預(yù)測(cè),落地時(shí)間精度顯著提高,所以時(shí)空匹配較好,具有較高的落點(diǎn)和時(shí)間聯(lián)合預(yù)報(bào)精度;④本文方法預(yù)報(bào)時(shí)間的標(biāo)準(zhǔn)差要大于另2種方法,這說明當(dāng)考慮了氣動(dòng)阻力參數(shù)后,再入段速度變化較為劇烈,導(dǎo)致導(dǎo)彈的預(yù)報(bào)落地時(shí)間具有較大的散布,較另2種方法與實(shí)際情況更為符合;⑤3種算法計(jì)算的誤差橢圓大小和形狀類似,橢圓的長(zhǎng)軸(表征縱向標(biāo)準(zhǔn)差)近似沿射向方向,縱向標(biāo)準(zhǔn)差(長(zhǎng)軸)顯著大于橫向標(biāo)準(zhǔn)差(短軸),前者約為后者的5倍,這是因?yàn)閷?dǎo)彈飛行時(shí)速度縱向分量顯著大于橫向分量。與前2種方法相比,本文算法計(jì)算的誤差橢圓縱向標(biāo)準(zhǔn)差(長(zhǎng)軸)有所增加,這是因?yàn)榭紤]氣動(dòng)阻力的影響后增加了導(dǎo)彈速度縱向分量變化的不確定性,從而增大了落點(diǎn)縱向散布。

以Monte Carlo打靶仿真2 000次獲得的落點(diǎn)散布為基準(zhǔn),分別與3種算法計(jì)算的誤差橢圓進(jìn)行對(duì)比,以驗(yàn)證算法的準(zhǔn)確性和有效性,比較結(jié)果如圖2所示。

圖2 打靶(2 000次)落點(diǎn)散布與計(jì)算誤差橢圓的比較

從結(jié)果可以看出:①解析幾何法計(jì)算的誤差橢圓精度最差,縱向和橫向都存在偏差,縱向偏差是由于沒有考慮再入氣動(dòng)阻力導(dǎo)致,因此整體向前偏移(實(shí)際落點(diǎn)由于受到阻力作用射程會(huì)更近一些),橫向偏差則主要是沒有考慮地球引力J2攝動(dòng)影響導(dǎo)致;②數(shù)值積分法的誤差橢圓精度優(yōu)于解析幾何法,由于考慮了地球引力J2攝動(dòng)影響,橫向偏差基本消除,但是由于沒有考慮再入氣動(dòng)阻力的影響,仍然存在縱向偏差;③本文提出的考慮再入氣動(dòng)阻力模型法預(yù)報(bào)精度最高,與實(shí)際的打靶結(jié)果最為符合,這是因?yàn)樵摲椒瓤紤]了地球引力J2攝動(dòng)影響,又考慮了再入氣動(dòng)阻力影響,有效消除了縱向偏差和橫向偏差。

同樣,以Monte Carlo仿真打靶2 000次獲得的落地時(shí)間散布為基準(zhǔn),分別與3種算法計(jì)算的落地時(shí)間均值和±3σ范圍進(jìn)行比較,如圖3所示,圖中,N為打靶次數(shù)。

圖3 打靶(2 000次)落地時(shí)間與預(yù)測(cè)均值方差的比較

結(jié)合表2和圖3可以看出,解析法和積分法的結(jié)果幾乎相同,由于其預(yù)報(bào)時(shí)間均值顯著小于真值,且標(biāo)準(zhǔn)差結(jié)果偏小(與實(shí)際不符),導(dǎo)致打靶結(jié)果幾乎都超出其預(yù)報(bào)均值±3σ范圍之外,這說明誤差分析結(jié)果與實(shí)際統(tǒng)計(jì)特性不匹配,這是因?yàn)?種方法均沒有考慮氣動(dòng)阻力的減速作用,導(dǎo)致預(yù)報(bào)的導(dǎo)彈落地時(shí)間均有所提前,從而帶來更大的落點(diǎn)預(yù)報(bào)誤差;而本文算法由于考慮了再入氣動(dòng)阻力減速的作用,預(yù)報(bào)的落地時(shí)間與打靶結(jié)果符合較好,預(yù)報(bào)均值與打靶均值較為接近(相差1 s),且打靶結(jié)果均落在其預(yù)報(bào)的均值±3σ范圍以內(nèi),說明本文算法的誤差分析結(jié)果與實(shí)際統(tǒng)計(jì)特性更為匹配,可靠度和準(zhǔn)確度更高。

6 結(jié)論

在彈道導(dǎo)彈早期預(yù)警階段,氣動(dòng)阻力參數(shù)對(duì)落點(diǎn)預(yù)報(bào)精度具有顯著影響,但通常無法獲知且無再入觀測(cè)數(shù)據(jù)對(duì)其進(jìn)行辨識(shí)。本文基于導(dǎo)彈末端速度變化范圍的限制,有效地預(yù)估了氣動(dòng)阻力參數(shù)的變化范圍,然后建立了考慮氣動(dòng)阻力、地球自轉(zhuǎn)和J2攝動(dòng)項(xiàng)的落點(diǎn)預(yù)報(bào)和誤差分析模型,并計(jì)算了落點(diǎn)誤差散布橢圓。仿真結(jié)果表明,與傳統(tǒng)的解析法和積分法相比,本文方法由于有效地預(yù)估出氣動(dòng)阻力參數(shù)范圍,落點(diǎn)預(yù)報(bào)精度尤其是預(yù)報(bào)時(shí)間精度顯著提高,計(jì)算的誤差散布橢圓與打靶結(jié)果最為吻合。對(duì)于導(dǎo)彈早期預(yù)警階段的彈道和落點(diǎn)預(yù)報(bào)工作具有較好的應(yīng)用價(jià)值。