李小元,王中原,常思江
(1.海軍研究院,北京 100161;2.南京理工大學(xué) 能源與動(dòng)力工程學(xué)院,江蘇 南京 210094)
炮彈在發(fā)射飛行過程中,存在火炮初始擾動(dòng)、初速或然誤差、炮彈結(jié)構(gòu)參數(shù)誤差、氣象條件誤差等,造成射彈散布。
彈道修正技術(shù)是指對(duì)發(fā)射出去的炮彈的一段飛行彈道參數(shù)進(jìn)行實(shí)時(shí)測量,同預(yù)定彈道參數(shù)進(jìn)行比較,進(jìn)行彈道偏差計(jì)算,形成控制指令,通過彈上控制機(jī)構(gòu)作用調(diào)節(jié)彈道向預(yù)定彈道逼近,從而大幅度減小炮彈散布的一種技術(shù)。如果僅對(duì)縱向彈道進(jìn)行調(diào)節(jié),稱為一維彈道修正彈。當(dāng)彈上控制機(jī)構(gòu)采用阻力環(huán)裝置來調(diào)節(jié)阻力,則為打遠(yuǎn)修近、阻力環(huán)體制的一維彈道修正彈,它是目前廣泛采用、結(jié)構(gòu)簡潔的一種一維彈道修正技術(shù)[1-3]。
一維彈道修正技術(shù)在實(shí)際中要取得好的應(yīng)用效果,一個(gè)重要問題是如何根據(jù)實(shí)測的一段飛行彈道參數(shù)快速準(zhǔn)確地辨識(shí)出該彈的彈道特征參數(shù),進(jìn)而準(zhǔn)確地預(yù)報(bào)后續(xù)彈道、落點(diǎn),并準(zhǔn)確地知道這發(fā)彈同預(yù)定彈道的偏差,從而實(shí)施彈道修正。目前,采用卡爾曼濾波、特征參數(shù)辨識(shí)方法是一種較為有效的方法。文獻(xiàn)[4-5]分別采用卡爾曼濾波、速度-時(shí)間序列匹配等方法對(duì)該問題進(jìn)行了研究。文獻(xiàn)[6]通過建立彈道濾波方程,采用卡爾曼濾波法外推彈道,對(duì)彈道修正彈的落點(diǎn)進(jìn)行了預(yù)測。為了進(jìn)一步提高彈箭密集度,文獻(xiàn)[7]比較了線性卡爾曼模型和擴(kuò)展卡爾曼模型在GPS量測數(shù)據(jù)下的差異,并分析了2種模型在不同條件下的優(yōu)缺點(diǎn)。但是上述文獻(xiàn)在數(shù)據(jù)測量精度、基礎(chǔ)阻力系數(shù)變化等因素對(duì)預(yù)報(bào)彈道落點(diǎn)的影響方面研究較少。
在卡爾曼濾波、特征參數(shù)辨識(shí)方法中,在對(duì)應(yīng)整個(gè)可提供的實(shí)時(shí)測量彈道段上,分為預(yù)報(bào)和校正2個(gè)環(huán)節(jié)交替實(shí)施,核心目的是快速、準(zhǔn)確辨識(shí)出廣義氣動(dòng)力系數(shù)。在預(yù)報(bào)環(huán)節(jié),彈的阻力系數(shù)狀況對(duì)預(yù)報(bào)結(jié)果極為敏感,影響很大;在校正環(huán)節(jié),雷達(dá)實(shí)時(shí)測量的彈道參數(shù)的精度對(duì)校正結(jié)果影響極大。它們對(duì)最終辨識(shí)出的廣義氣動(dòng)力系數(shù)影響很大,直接影響在線彈道預(yù)報(bào)精度。本文對(duì)阻力環(huán)修正體制下的這些問題進(jìn)行分析、討論,給出相關(guān)建議。
針對(duì)本文研究的一維彈道修正彈,根據(jù)其工作原理,可建立一維修正彈道方程:
(1)
控制方程為
(2)
式中:Fx2,Fy2,Fz2為作用在彈丸上的合力F在彈道坐標(biāo)系x2y2z2中的投影,在計(jì)算Fx2,Fy2時(shí)需要使用阻力系數(shù)符合系數(shù)kcx和升力系數(shù)符合系數(shù)kcy,上述2個(gè)系數(shù)為待辨識(shí)的參數(shù);Mξ,Mη,Mζ為作用在彈丸質(zhì)心處的合力矩M在彈軸坐標(biāo)系ξηζ中的投影;MD為阻力環(huán)附加的極阻尼力矩;FD為阻力環(huán)張開后增加的阻力;ΔCx為阻力系數(shù)增量;v為速度;θa為彈道傾角;ψ2為彈道偏角;ωξ,ωη,ωζ為彈丸角速度在彈軸坐標(biāo)系ξηζ中的分量;φa為彈軸高低擺動(dòng)角;φ2為彈軸側(cè)向擺動(dòng)角;γ為滾轉(zhuǎn)角;x,y,z為彈丸在地面坐標(biāo)內(nèi)的位置分量;JA,JC分別為極轉(zhuǎn)動(dòng)慣量和赤道轉(zhuǎn)動(dòng)慣量;m為彈丸質(zhì)量;ρ為大氣密度;S為彈丸特征面積。各坐標(biāo)系的定義參見文獻(xiàn)[8]。
如前所述,如何根據(jù)實(shí)測的一段飛行彈道參數(shù)快速、準(zhǔn)確地辨識(shí)出該彈的彈道特征參數(shù)是一個(gè)核心問題。本文選取阻力系數(shù)符合系數(shù)kcx和升力系數(shù)符合系數(shù)kcy作為彈道特征參數(shù),分別用于縱向射程與橫向距離的彈道、落點(diǎn)預(yù)報(bào)。為了能夠利用一段實(shí)測彈道數(shù)據(jù)盡可能準(zhǔn)確地估計(jì)出這2個(gè)彈道特征參數(shù),本文采用卡爾曼濾波方法,將kcx和kcy作為狀態(tài)變量進(jìn)行估計(jì)。卡爾曼濾波方法的優(yōu)點(diǎn)是考慮了物理系統(tǒng)本身的動(dòng)態(tài)模型和測量系統(tǒng)的量測誤差,利用從某時(shí)刻k-1以前直到現(xiàn)在的測量數(shù)據(jù)Y(j個(gè))逐次遞推,可估計(jì)出k時(shí)刻最優(yōu)的系統(tǒng)狀態(tài)X??柭鼮V波基本方程簡介如下。
預(yù)測方程:
(3)
(4)
(5)
濾波方程:
(6)
卡爾曼濾波增益:
(7)
協(xié)方差矩陣:
Pk=(I-KkHk)Pk/k-1
(8)
式中:Δt為卡爾曼濾波的采樣間隔;Φk為基本矩陣,H為量測矩陣,Qk為過程噪聲方差矩陣,Rk為量測噪聲方差矩陣,其他各符號(hào)定義參見文獻(xiàn)[9]。
下面以上述卡爾曼濾波方程為基礎(chǔ),結(jié)合具體的彈道模型,構(gòu)建一維修正彈的彈道濾波模型。
由于在利用卡爾曼濾波辨識(shí)彈道特征參數(shù)的過程中,不涉及阻力環(huán)機(jī)構(gòu)作用。因此,在構(gòu)建一維修正彈彈道濾波模型時(shí),令方程(1)中FD=0,并將參數(shù)kcx和kcy分別置于阻力系數(shù)和升力系數(shù)之前作為乘子,并補(bǔ)充如下微分方程:
(9)
則狀態(tài)變量取為
X=(vθaψ2xyzωξωηωζφaφ2γkcxkcy)T
(10)
將方程(1)改寫成狀態(tài)空間形式為
(11)
式中:fi(i=1,2,…,14)為方程(1)和方程(9)對(duì)應(yīng)的右端函數(shù)。
對(duì)方程(11)進(jìn)行線性化,可得:
(12)
式中:F為動(dòng)力矩陣,其元素為fi(i=1,2,…,14)對(duì)每個(gè)狀態(tài)變量的偏導(dǎo)數(shù),由于具體形式冗長,為節(jié)約篇幅,這里略去。
在求出了動(dòng)力矩陣F的具體表達(dá)式后,可求得彈道濾波系統(tǒng)的基本矩陣Φk為
Φk≈I14×14+F14×14·Δt
(13)
式中:I14×14為單位矩陣。
假設(shè)可直接或間接測得的彈道參數(shù)為
Z=(vθaψ2xyz)T
(14)
則彈道濾波系統(tǒng)的量測方程可以表達(dá)為
Z=h(X)+n=H·X+n
(15)
式中:n為量測噪聲。
由于狀態(tài)變量X為14×1維矩陣,Z為6×1維矩陣,故量測矩陣H為6×14維矩陣,n為6×1維矩陣。
(16)
量測噪聲n的方差矩陣為
(17)
式中:σ為測量值的標(biāo)準(zhǔn)差。
另外,與狀態(tài)方程相對(duì)應(yīng)的還有過程噪聲方差矩陣Qk,由于很難描述彈道模型與實(shí)際彈丸飛行狀態(tài)的差異,故一般取Qk=O14×14。
數(shù)值迭代求解如下形式的Ricatti方程:
(18)
式中:Pk/k-1為k時(shí)刻數(shù)據(jù)更新前的狀態(tài)估計(jì)值協(xié)方差矩陣;Pk為數(shù)據(jù)更新后的狀態(tài)估計(jì)值協(xié)方差矩陣,相應(yīng)地,Pk-1為上一時(shí)刻的狀態(tài)估計(jì)值協(xié)方差矩陣;Kk為卡爾曼增益矩陣。
為了使卡爾曼濾波系統(tǒng)正常工作,狀態(tài)估計(jì)值協(xié)方差矩陣的初始值P0根據(jù)量測噪聲的標(biāo)準(zhǔn)偏差來選取,有:
P0=diag(P0,11,P0,22,…,P0,1414)
(19)
卡爾曼濾波的基本方程為
(20)
根據(jù)卡爾曼濾波方程(20)可以實(shí)時(shí)估計(jì)出彈道方程組對(duì)應(yīng)的狀態(tài)變量,將測量數(shù)據(jù)段最后一點(diǎn)的狀態(tài)變量最優(yōu)估計(jì)值作為彈道預(yù)報(bào)的起點(diǎn),對(duì)彈道方程組進(jìn)行數(shù)值積分直至落點(diǎn),可預(yù)報(bào)出阻力環(huán)機(jī)構(gòu)全彈道不作用條件下的該發(fā)彈的無控彈道及落點(diǎn)。將事先獲取的阻力環(huán)阻力系數(shù)ΔCx代入彈道方程組中(即FD≠0),并預(yù)設(shè)一個(gè)阻力環(huán)啟動(dòng)作用時(shí)刻,可預(yù)報(bào)出該發(fā)彈在某一時(shí)刻張開阻力環(huán)后的修正彈道及落點(diǎn)。如果將無控彈道預(yù)測落點(diǎn)與目標(biāo)彈的射程相減,取得該發(fā)彈的射程修正量,利用該修正量進(jìn)行迭代計(jì)算,可很容易求出該射程修正量對(duì)應(yīng)的阻力環(huán)機(jī)構(gòu)啟動(dòng)時(shí)刻。
為了驗(yàn)證上述模型與方法的正確性和有效性,本節(jié)以某大口徑一維彈道修正彈的實(shí)際數(shù)據(jù)為例,開展在線彈道預(yù)報(bào)分析。在該一維彈道修正彈的某次試驗(yàn)中,連續(xù)射擊6發(fā)無控彈,采用兩部彈道跟蹤雷達(dá)同時(shí)對(duì)每發(fā)彈丸實(shí)時(shí)跟蹤,分別測得對(duì)應(yīng)的飛行彈道數(shù)據(jù)。這兩部雷達(dá)的跟蹤精度有差異,其中一部為高精度彈道跟蹤雷達(dá)(記為HA雷達(dá)),而另一部的跟蹤精度相對(duì)較低(記為LA雷達(dá))。采用不同精度雷達(dá)同時(shí)跟蹤的目的是為了研究測量數(shù)據(jù)精度對(duì)彈道預(yù)報(bào)的影響。需要說明的是,由于一維彈道修正彈的彈道預(yù)報(bào)主要關(guān)注縱向距離x,故下面分析均對(duì)應(yīng)預(yù)報(bào)縱向距離x。
具體分析包括幾個(gè)方面:研究不同阻力系數(shù)對(duì)彈道預(yù)報(bào)精度的影響;分析測量飛行彈道數(shù)據(jù)對(duì)彈道預(yù)報(bào)的影響。
利用HA雷達(dá)實(shí)測5 s數(shù)據(jù),包括彈丸速度和位置,進(jìn)行在線彈道預(yù)報(bào),所得結(jié)果如圖1~圖3和表1所示。
其中,圖1~圖3分別為其中1發(fā)彈丸的位置、速度及阻力系數(shù)符合系數(shù)kcx濾波曲線,表1為6發(fā)彈丸預(yù)報(bào)射程與實(shí)測射程之差。從圖1~圖3可以看出,采用卡爾曼濾波方法得到了彈丸位置、速度及阻力系數(shù)符合系數(shù)的最優(yōu)估計(jì)值。由表1所示的在線射程預(yù)報(bào)結(jié)果看,6發(fā)預(yù)報(bào)射程與實(shí)測射程之差出現(xiàn)了正負(fù)號(hào),即對(duì)于這一組射彈,有的射程預(yù)報(bào)值大于射程實(shí)測值,有的則小于實(shí)測值。這表明由于各種隨機(jī)誤差因素的存在,在線彈道預(yù)報(bào)誤差本身也存在隨機(jī)性,而一維彈道修正彈的主要功能是提高射程密集度,因此應(yīng)盡可能減小彈道預(yù)報(bào)誤差的跳差范圍,即盡量控制極差。從數(shù)值上看,試驗(yàn)彈的平均射程為45 700 m,而6發(fā)彈的射程預(yù)報(bào)極差為73 m(約為平均射程的0.16%),射程預(yù)報(bào)誤差平均值為21 m(不到平均射程的0.05%),能夠滿足實(shí)際工程應(yīng)用的要求。
圖1 彈丸位置實(shí)測值與最優(yōu)估計(jì)值
圖2 彈丸速度實(shí)測值與最優(yōu)估計(jì)值
圖3 參數(shù)kcx的最優(yōu)估計(jì)結(jié)果
表1 6發(fā)彈丸預(yù)報(bào)射程與實(shí)測射程之差
在彈丸阻力系數(shù)標(biāo)準(zhǔn)值基礎(chǔ)上變化±2.5%和±5%,將變化后的阻力系數(shù)代入前面建立的在線彈道預(yù)報(bào)模型中,所得預(yù)報(bào)射程與實(shí)測射程之差ΔXc列于表2中。
表2 不同基礎(chǔ)阻力系數(shù)對(duì)應(yīng)的預(yù)報(bào)射程與實(shí)測射程之差
表2中結(jié)果說明,彈丸阻力系數(shù)標(biāo)準(zhǔn)值變化對(duì)預(yù)報(bào)結(jié)果有一定影響。按照相對(duì)于標(biāo)準(zhǔn)值的變化百分比,預(yù)報(bào)射程誤差呈線性變化。從一組數(shù)據(jù)來看,預(yù)報(bào)射程誤差相對(duì)于標(biāo)準(zhǔn)值條件下的預(yù)報(bào)誤差產(chǎn)生了一個(gè)系統(tǒng)偏差,6發(fā)彈丸預(yù)報(bào)射程誤差的極差變化很小。從這個(gè)意義上講,彈丸阻力系數(shù)標(biāo)準(zhǔn)值發(fā)生變化時(shí),主要影響在線彈道預(yù)報(bào)結(jié)果的準(zhǔn)確度,而不影響反映彈丸發(fā)與發(fā)之間差異的散布,即對(duì)一維修正落點(diǎn)的密集度無影響。從武器使用角度,應(yīng)兼顧準(zhǔn)確度和密集度,故有必要采用盡可能準(zhǔn)確的阻力系數(shù)。由于本文通過大量試驗(yàn)取得了較為精確的彈丸阻力系數(shù),因此預(yù)報(bào)準(zhǔn)確度也相對(duì)較高。
試驗(yàn)中LA雷達(dá)斜距離標(biāo)準(zhǔn)差比HA雷達(dá)斜距離標(biāo)準(zhǔn)差大2 m,LA雷達(dá)高低角標(biāo)準(zhǔn)差比HA雷達(dá)高低角標(biāo)準(zhǔn)差大0.65 mrad,LA雷達(dá)方位角標(biāo)準(zhǔn)差比HA雷達(dá)方位角標(biāo)準(zhǔn)差大0.5 mrad。表3所示為采用HA和LA兩部雷達(dá)跟蹤數(shù)據(jù)預(yù)報(bào)出的射程值,其他條件參數(shù)均取為相同。Xm為實(shí)測射程,Xc,HA為HA雷達(dá)預(yù)報(bào)射程,Xc,LA為LA雷達(dá)預(yù)報(bào)射程,ΔXc,HA=Xc,HA-Xm,ΔXc,LA=Xc,LA-Xm。
表3 利用兩部雷達(dá)跟蹤數(shù)據(jù)預(yù)報(bào)的射程
從平均值上看,基于LA雷達(dá)跟蹤數(shù)據(jù)的預(yù)報(bào)誤差要小于基于HA雷達(dá)跟蹤數(shù)據(jù)的預(yù)報(bào)誤差,特別是對(duì)于2號(hào)彈、3號(hào)彈以及6號(hào)彈。但對(duì)于一維彈道修正彈,其原理主要是減小彈丸的縱向散布,故對(duì)于在線彈道預(yù)報(bào)而言,更重要的是盡可能減小一組射彈射程預(yù)報(bào)誤差的極差。這一極差主要受跟蹤數(shù)據(jù)精度的影響,由表3可知,HA雷達(dá)對(duì)應(yīng)的射程預(yù)報(bào)誤差極差為73 m,而LA雷達(dá)對(duì)應(yīng)的射程預(yù)報(bào)誤差極差為116 m。因此,應(yīng)盡可能采用跟蹤精度高的雷達(dá)。
測量數(shù)據(jù)的長度對(duì)預(yù)報(bào)結(jié)果也有影響。前面計(jì)算中,都是采用了5 s跟蹤數(shù)據(jù)。下面針對(duì)HA雷達(dá)數(shù)據(jù),以1號(hào)彈為例,研究不同數(shù)據(jù)長度對(duì)預(yù)報(bào)結(jié)果的影響,如圖4所示。
圖4 不同數(shù)據(jù)長度對(duì)應(yīng)的射程預(yù)報(bào)誤差
由圖4可以看出跟蹤數(shù)據(jù)長度對(duì)預(yù)報(bào)結(jié)果的影響。當(dāng)數(shù)據(jù)時(shí)長大于5 s,預(yù)報(bào)誤差略有波動(dòng),但基本上維持在25 m左右,數(shù)據(jù)時(shí)長越長對(duì)預(yù)報(bào)誤差的影響較小,而數(shù)據(jù)越多,在線計(jì)算量也越大;當(dāng)數(shù)據(jù)時(shí)長小于5 s,預(yù)報(bào)誤差隨著數(shù)據(jù)時(shí)長的減小而迅速增大,3 s時(shí)長對(duì)應(yīng)的預(yù)報(bào)誤差已接近70 m。由于圖4所示是1號(hào)彈對(duì)應(yīng)的結(jié)果,如果是對(duì)于6號(hào)彈,3 s時(shí)長對(duì)應(yīng)的預(yù)報(bào)誤差接近150 m,對(duì)修正彈最終密集度指標(biāo)的實(shí)現(xiàn)極為不利。上述結(jié)果表明,當(dāng)綜合考慮在線彈道預(yù)報(bào)的精度和計(jì)算量(實(shí)時(shí)性),選取5 s數(shù)據(jù)是比較合適的。
前面所建基于卡爾曼濾波的在線彈道預(yù)報(bào)模型中,參數(shù)kcx也是作為一個(gè)狀態(tài)變量。實(shí)際研究中發(fā)現(xiàn),該狀態(tài)變量對(duì)應(yīng)的協(xié)方差矩陣元素初始值Pk,0對(duì)參數(shù)kcx的最優(yōu)估計(jì)結(jié)果有影響,并且較為敏感,進(jìn)而影響射程預(yù)報(bào)的結(jié)果。本節(jié)仍采用HA雷達(dá)跟蹤數(shù)據(jù),研究不同Pk,0值對(duì)預(yù)報(bào)結(jié)果的影響,結(jié)果如表4所示。表中,ΔXc為預(yù)報(bào)射程與實(shí)測射程之差。
表4 不同Pk,0值對(duì)應(yīng)的預(yù)報(bào)射程與實(shí)測射程之差
圖5為不同Pk,0值對(duì)應(yīng)的射程預(yù)報(bào)誤差極差(ΔXc,m)曲線圖。
圖5 不同Pk,0值對(duì)應(yīng)的射程預(yù)報(bào)誤差極差
從表4和圖5中不難看出,盡管Pk,0的量值不大,但對(duì)預(yù)報(bào)結(jié)果卻具有較大影響。前面研究表明,對(duì)于一維彈道修正彈,一組射程預(yù)報(bào)誤差的極差大小可在一定程度上反映在線彈道預(yù)報(bào)的性能。當(dāng)Pk,0>0.01,隨著Pk,0值的增大,射程預(yù)報(bào)誤差的平均值和極差都增大;當(dāng)Pk,0<0.01,隨著Pk,0值的減小,射程預(yù)報(bào)誤差的平均值和極差也增大。圖5所示極差曲線中,顯然Pk,0=0.01是一個(gè)最優(yōu)值。需要指出的是,對(duì)于不同彈丸、不同測量設(shè)備或不同模型,該最優(yōu)值也不盡相同,可采用試驗(yàn)手段繪制類似圖5的極差曲線,從而確定Pk,0的最優(yōu)值。工程實(shí)踐表明,這一方法是比較有效的。
本文針對(duì)阻力環(huán)體制的一維彈道修正彈,建立了飛行控制彈道模型,并基于卡爾曼濾波理論建立了彈道濾波模型,給出了彈道預(yù)報(bào)的方法,并通過算例對(duì)本文研究結(jié)果的正確性和有效性進(jìn)行了分析和驗(yàn)證。通過對(duì)實(shí)際工程應(yīng)用結(jié)果的分析, 提出以下建議:為了提高無控彈道和修正彈道的預(yù)報(bào)精度,應(yīng)盡可能準(zhǔn)確地獲取全彈的阻力系數(shù),系數(shù)值應(yīng)覆蓋全彈道馬赫數(shù)。在確定某一實(shí)時(shí)測量彈道參數(shù)的基礎(chǔ)上,進(jìn)行若干發(fā)炮射試驗(yàn),通過反復(fù)調(diào)整預(yù)報(bào)-校正分段、測量參數(shù)方差同實(shí)測彈道的對(duì)比,確定出較佳匹配(協(xié)方差、分段等)關(guān)系,來提高“廣義阻力系數(shù)”的辨識(shí)效果,提高后續(xù)彈道預(yù)報(bào)精度。