何山，吳盤龍，李星秀，惲鵬，薄煜明

(1.南京理工大學(xué) 自動化學(xué)院, 江蘇 南京 210094; 2.南京理工大學(xué) 理學(xué)院, 江蘇 南京 210094)

0 引言

隨著網(wǎng)絡(luò)化、信息化技術(shù)的快速發(fā)展，防空反導(dǎo)作為現(xiàn)代戰(zhàn)爭的重要任務(wù)之一，發(fā)揮著越來越重要地位[1]。高炮以其靈活性高、射速快、機(jī)動性強(qiáng)和火力密集等優(yōu)點(diǎn)成為了當(dāng)前低空防御的主力[2]。

火控系統(tǒng)作為武器系統(tǒng)的神經(jīng)中樞，其解算出的射擊諸元精度直接影響武器系統(tǒng)的命中概率[3]。一般來說，火控系統(tǒng)解算的難點(diǎn)在于目標(biāo)航跡和彈丸軌跡的準(zhǔn)確預(yù)測[4]。而隨著末端突防能力的日漸提升，目標(biāo)強(qiáng)機(jī)動性使得在航跡的建模和狀態(tài)估計(jì)上變得尤為困難；另一方面，復(fù)雜的戰(zhàn)場環(huán)境和彈丸的一致性直接影響彈丸零升阻力、誘導(dǎo)阻力、升力以及馬格努斯力等，從而給彈丸的彈道模型建立帶來了極大的挑戰(zhàn)。

現(xiàn)行防空火控系統(tǒng)中的彈道模型多以射表函數(shù)逼近為主，極大縮短了解算時間[5-6]。但在現(xiàn)今信息化作戰(zhàn)的背景下，這一方法存在著精度有限、通用性差、設(shè)計(jì)過程繁瑣等缺點(diǎn)，尤其是只能對各個彈種的具體彈道和各個修正量分別進(jìn)行解算，極大地制約了火控系統(tǒng)的開發(fā)[7]。對此，基于彈道微分方程組的火控外彈道模型在該背景下應(yīng)運(yùn)而生，它具有通用性好、精度高等特點(diǎn)，但解算處理復(fù)雜，耗費(fèi)時間較多，尤其是對動目標(biāo)射擊時，需以逆解法[6]為工具，通過試探法不斷迭代求解出命中點(diǎn)的估計(jì)值，因此計(jì)算量更大。

目前工程應(yīng)用中所選用的彈道微分方程組主要是質(zhì)點(diǎn)彈道方程組和修正的質(zhì)點(diǎn)彈道方程組[8]。質(zhì)點(diǎn)彈道方程組是將彈丸作為質(zhì)點(diǎn)處理，忽略了作用在彈丸上的升力、馬氏力和攻角引起的誘導(dǎo)阻力，使彈道計(jì)算結(jié)果產(chǎn)生了一定的誤差；修正質(zhì)點(diǎn)彈道模型在質(zhì)點(diǎn)彈道模型的基礎(chǔ)上綜合了這些因素，有效地提高了計(jì)算精度。然而對于高轉(zhuǎn)速彈丸，修正質(zhì)點(diǎn)彈道模型在偏流的計(jì)算上存在著一定的誤差。對此有學(xué)者將彈丸在空中的運(yùn)動視為一般剛體的運(yùn)動，考慮彈丸飛行期間產(chǎn)生的進(jìn)動、章動和自轉(zhuǎn)等，極大地提高了對彈丸軌跡的預(yù)測精度[9]，但復(fù)雜的剛體彈道方程解算很難在實(shí)際火控系統(tǒng)的一個工作周期內(nèi)完成。

鑒于此，本文將剛體彈道方程引入到防空火控系統(tǒng)的彈道模型中，極大地提升了火控系統(tǒng)整體的通用性與精確性。同時為了避免剛體彈道方程計(jì)算時給火控帶來的實(shí)時性問題，本文從火控解算的方法入手，將無偏轉(zhuǎn)換得到的噪聲協(xié)方差矩陣做解耦，降低了目標(biāo)狀態(tài)估計(jì)時濾波的復(fù)雜度；并在割弦迭代法中，將虛擬脫靶量引入到射角和射向的迭代- 修正過程中，減少了火控解算中對彈道方程計(jì)算次數(shù)，加快射擊諸元求解的收斂速度。從而在算法執(zhí)行機(jī)制上，有效地克服了基于剛體彈道方程火控解算實(shí)時性不足的缺陷。

1 彈道建模

由于影響彈丸在空氣中運(yùn)動的因素極為錯綜復(fù)雜，本文針對某型防空高炮，根據(jù)各個因素對彈丸飛行影響的大小，在六自由度剛體彈道方程的基礎(chǔ)上[10],忽略科式慣性力和赤道阻尼力矩；同時為了增大龍格- 庫塔法求解時的步長，忽略彈軸高低角和方位角的2階導(dǎo)項(xiàng)，建立如下彈道方程：

(1)

2 脫靶量

如圖1所示，以武器回轉(zhuǎn)中心O為坐標(biāo)原點(diǎn)構(gòu)建地面基準(zhǔn)坐標(biāo)系Oxyz，其中：Ox軸為射面與炮口水平面的交線，指向射擊方向?yàn)檎?，Oy軸在射擊面內(nèi)，垂直O(jiān)x軸，向上為正；彈目偏差觀測坐標(biāo)系TxEyEzE以目標(biāo)中心T為原點(diǎn)，沿著觀目矢量OT為彈目偏差觀測坐標(biāo)系的縱深軸yE，過OT的鉛錘面與迎光面(也叫Q平面)的交線為高低軸zE，向上為正，在迎光面內(nèi)，垂直于zE軸的規(guī)定為方位軸xE，方向滿足右手法則[11-12]；當(dāng)彈丸以速度v發(fā)射時，記T*為目標(biāo)T在Oxz面上的投影，Z*為彈跡穿過迎光面的穿越點(diǎn)，Δβ和Δε分別表示觀目矢量OT與觀炸矢量OZ*之間張角的水平分量和高低分量，βh為觀目矢量在東、北、天地理坐標(biāo)系下的方位角，εh為觀目矢量在東、北、天地理坐標(biāo)系下的高低角。則彈目偏差TZ*為

TZ*=(xE,0,yE,0,zE,0)=
(|OT|tan Δβ,0,|OT|tan Δε).

(2)

圖1 脫靶量示意圖Fig.1 Schematic diagram of miss distance

而在火控系統(tǒng)中，彈目偏差的數(shù)據(jù)處理都是建立在東、北、天地理坐標(biāo)系下的[6]。因此，實(shí)際地面基準(zhǔn)坐標(biāo)系下的脫靶量Δβ′、Δε′分別為

(3)

3 快速防空火控解算方法

3.1 解耦無偏轉(zhuǎn)換卡爾曼濾波

為了準(zhǔn)確描述目標(biāo)狀態(tài)隨著時間變化的過程，將目標(biāo)位置、速度和加速度作為狀態(tài)量，構(gòu)建目標(biāo)的運(yùn)動方程

Xk=Fk-1Xk-1+Wk-1，

(4)

火控系統(tǒng)中對目標(biāo)狀態(tài)的估計(jì)離不開探測設(shè)備對目標(biāo)的跟蹤，在探測系統(tǒng)中，傳感器的量測信息通常是建立在三維球坐標(biāo)中，目標(biāo)與傳感器之間的相對位置如圖2所示。圖2中，rk、ηk、βk分別為徑向距離、高低角、方位角。

圖2 傳感器測量坐標(biāo)系Fig.2 Sensor measurement coordinate system

傳感器的量測值Zk主要包括徑向距離rk、高低角ηk和方位角βk，其測量方程[13-14]為

(5)

式中：Vrk、Vηk和Vβk是相互獨(dú)立、均值為0且恒定方差為σrk、σηk和σβk的高斯白噪聲。

球坐標(biāo)可通過(6)式轉(zhuǎn)化到笛卡爾坐標(biāo)系下的量測：

(6)

因此，根據(jù)(5)式的量測值Zk，對真實(shí)均值和協(xié)方差矩陣求數(shù)學(xué)期望得到無偏量測偏差μk和協(xié)方差Rk[15]：

(7)

(8)

式中：

(9)

3.2 基于剛體彈道模型的快速射擊諸元求解

3.2.1 命中方程建立

假設(shè)觀瞄設(shè)備中心O′與武器回轉(zhuǎn)中心O不在同一位置，記J=OO′為觀炮基線，則實(shí)際防空武器的命中多邊形如圖3所示。

圖3 命中多邊形Fig.3 Hit polygon

相應(yīng)的命中矢量方程為

Dh(t)-D(t)-J-Sh(t)=0，

(10)

式中：Dh(t)=D(t+tf)，tf為彈丸飛行時間。

為了計(jì)算方便，假設(shè)目標(biāo)在彈丸飛行時間內(nèi)做勻速直線運(yùn)動，則

(11)

通過將(10)式中的矢量在東、北、天地理坐標(biāo)系上進(jìn)行投影，得到標(biāo)量命中方程組

(12)

式中：(xoo′,yoo′,zoo′)和(xTh,yTh,zTh)分別為瞄具中心O和命中點(diǎn)Th的坐標(biāo)；瞄準(zhǔn)矢量D(t)=(T,T,T).

由(1)式可知，當(dāng)射角φh和彈丸飛行時間tf已知時，利用4階龍格- 庫塔法可求得彈丸軌跡上任意一點(diǎn)坐標(biāo)，即

(13)

再根據(jù)射向方位角βh，命中點(diǎn)坐標(biāo)可寫為

(14)

由(12)式、(13)式和(14)式即可求得射擊諸元(φh,βh,tf)。

3.2.2 基于虛擬脫靶量的迭代——修正算法

考慮到火控解算的精度與計(jì)算速度等要求，可利用割弦迭代法計(jì)算出命中點(diǎn)坐標(biāo)，在每次割弦過程中，需利用多次修正龍格- 庫塔法求解彈道模型時的初值，才能計(jì)算出符合當(dāng)前預(yù)測未來點(diǎn)的射擊諸元。因此，本文引入虛擬脫靶量概念，即目標(biāo)真實(shí)航跡與計(jì)算機(jī)仿真彈丸軌跡間偏差量，對初始射擊射角和射向進(jìn)行修正，極大地提高了計(jì)算當(dāng)前預(yù)測未來點(diǎn)射擊諸元的收斂速度。具體求解過程如下：

步驟2割弦迭代法初值的選取。假設(shè)在彈丸飛行時間內(nèi)目標(biāo)處于靜止?fàn)顟B(tài)，即目標(biāo)運(yùn)動時間Tf=0 s，則命中點(diǎn)坐標(biāo)為(xTh,k,yTh,k,zTh,k)，有

(15)

這時可將火控解算問題轉(zhuǎn)化為一個外彈道學(xué)中的兩點(diǎn)邊值問題，即從上述邊值條件出發(fā)，求出射擊諸元。

步驟3彈道方程初值給定。根據(jù)所給的目標(biāo)未來點(diǎn)坐標(biāo)預(yù)估初始射角φh，0和射向方位角βh，0:

(16)

步驟4彈道方程求解。利用4階龍格- 庫塔法求解微分方程組(1)式，得到地面基準(zhǔn)坐標(biāo)系下彈道諸元的序列值：

{xi,yi,zi,vx,i,vy,i,vz,i,tf,i},i=0,1,2,…，

(17)

式中：i為龍格- 庫塔法定步長求解時的序列號。

步驟5虛擬脫靶量計(jì)算。

1)將步驟4得到的彈道諸元序列通過(18)式轉(zhuǎn)化到彈目偏差觀測坐標(biāo)系中，

(18)

2)利用拉格朗日插值，可計(jì)算出在Q平面內(nèi)的二維彈目偏差(xE,0,zE,0)和彈丸飛行時間tf,i,0，利用(3)式將(xE,0,zE,0)轉(zhuǎn)換到東、北、天地理坐標(biāo)系下的虛擬脫靶量(Δβ′0,Δε′0)。

步驟6射角與射向迭代修正。判斷(Δβ′0,Δε′0)是否滿足系統(tǒng)所需的精度(一般取0.1 mil)。若滿足，則射擊諸元為φh=φh,0，βh=βh,0，tf=tf,i,0，并轉(zhuǎn)到步驟7；否則對初始射角和射向進(jìn)行修正，得到新的射角和射向方位角：

(19)

將修正后的射角和射向取代原來的變量，并返回步驟4.

至此，一個周期內(nèi)的火控解算完成，相比傳統(tǒng)火控解算的迭代修正算法，步驟5中利用虛擬脫靶量進(jìn)行射角和射向的迭代修正，減少了火控解算中對彈道方程計(jì)算次數(shù)，極大提高了基于剛體彈道方程的防空火控解算方法求解速度。

4 仿真分析

4.1 仿真場景

本文以某型防空高炮為研究對象，構(gòu)建三維仿真場景，對所提火控解算方法進(jìn)行仿真驗(yàn)證。假設(shè)目標(biāo)做勻速運(yùn)動，起始位置為(3 500 m, 3 500 m, 1 000 m)，初始速度為(-300 m/s, -200 m/s, -25 m/s)；測量傳感器采用雷達(dá)，采樣周期為0.02 s，距離誤差為5 m，方位和俯仰角誤差為1 mrad，目標(biāo)在雷達(dá)界面的顯示如圖4所示。

圖4 目標(biāo)在雷達(dá)界面顯示Fig.4 Target on radar interface

彈丸初速v0=1 050 m/s，龍格- 庫塔法的解算步長為0.01 s，在標(biāo)準(zhǔn)氣象條件下，不同射角的彈丸軌跡如圖5所示。

圖5 彈丸軌跡Fig.5 Projectile trajectory

4.2 仿真分析

4.2.1 狀態(tài)濾波

在實(shí)際工程中，通過傳感器測量的目標(biāo)運(yùn)動參數(shù)包含眾多干擾噪聲，這些噪聲不僅影響系統(tǒng)的計(jì)算精度，嚴(yán)重時系統(tǒng)將失穩(wěn)。因此對目標(biāo)進(jìn)行狀態(tài)濾波是準(zhǔn)確預(yù)測目標(biāo)未來點(diǎn)位置的關(guān)鍵技術(shù)。本文通過無偏轉(zhuǎn)換將雷達(dá)對目標(biāo)運(yùn)動狀態(tài)的非線性量測方程圍繞目標(biāo)運(yùn)動狀態(tài)進(jìn)行展開，并將無偏轉(zhuǎn)換得到的噪聲協(xié)方差矩陣做解耦，降低了濾波算法計(jì)算的復(fù)雜度，其位置均方根誤差(RMSE)、速度RMSE和濾波處理時間如圖6、圖7和圖8所示。由圖6～圖8可以看出，解耦無偏轉(zhuǎn)換卡爾曼濾波與無偏轉(zhuǎn)換卡爾曼濾波算法[15]的誤差差值不大于總誤差的1%，而濾波時間減少了15%. 因此，該方法具有較高的精度和良好的性能。

圖6 位置RMSEFig.6 RMSEs of position

圖7 速度RMSEFig.7 RMSEs of velocity

圖8 濾波的處理時間Fig.8 Processing time of filtering

4.2.2 火控解算誤差

根據(jù)上述仿真場景，以2 s穩(wěn)定建航后，分別用傳統(tǒng)的火控解算方法[6]和本文所提方法對不同時間點(diǎn)進(jìn)行火控解算，其解算點(diǎn)的射擊諸元對比如表1所示，火控解算的誤差如圖9所示。

從圖9仿真結(jié)果可以看出，本文所提方法的有效性和可行性。由表1可以看出：當(dāng)目標(biāo)在較遠(yuǎn)處，其解算誤差很大，但隨著目標(biāo)越來越接近，彈丸飛行時間越來越短，誤差呈現(xiàn)減小趨勢，這是由于測量傳感器所得目標(biāo)數(shù)據(jù)存在著一定的量測誤差，導(dǎo)致解算誤差大部分來源于濾波誤差；當(dāng)目標(biāo)距離3 000 m以內(nèi)時，單軸方向上的解算誤差不大于0.5 m，該方法可以很好地滿足實(shí)際防空作戰(zhàn)需求。

表1 射擊諸元對比Tab.1 Comparison of firing data

圖9 火控解算誤差Fig.9 Errors of fire control calculation

4.2.3 解算時間

本文通過對不同距離目標(biāo)的火控解算進(jìn)行仿真實(shí)驗(yàn)來有效地說明所提出方法的實(shí)時性優(yōu)勢。在Matlab2014環(huán)境下(Windows 7, i5-3210M CPU, 2.5 GHz主頻, 8 GB RAM)將傳統(tǒng)的火控解算方法[6]與本文方法進(jìn)行比較，解算時間對比如表2所示。

表2 Windows平臺下的解算時間對比Tab.2 Comparison of calculation times based on Windows s

從表2中可以看出，本文所提方法在計(jì)算量上遠(yuǎn)遠(yuǎn)低于傳統(tǒng)的火控解算方法[6]，這是由于火控諸元的解算主要包括狀態(tài)濾波與命中點(diǎn)解算兩大步驟。對于狀態(tài)濾波，通過將無偏轉(zhuǎn)換卡爾曼濾波的噪聲協(xié)方差做解耦，一定程度上降低了濾波計(jì)算的復(fù)雜度；對于命中點(diǎn)解算，本文提出基于虛擬脫靶量的迭代——修正算法，利用虛擬脫靶量進(jìn)行射角和射向的迭代修正，減少了火控解算中對彈道方程計(jì)算次數(shù)，大幅度提高火控求解速度。

為了有效說明該方法在工程時是可行的，將代碼改編成C語言并移植到配置為i7-4790CPU、SylixOS實(shí)時操作系統(tǒng)的工控機(jī)上運(yùn)行，解算時間對比如表3所示。從表3實(shí)驗(yàn)的數(shù)據(jù)中可以看出，本文所提方法在工程上具有可觀的應(yīng)用前景。

5 結(jié)論

為了提高防空火控解算的精確性與實(shí)時性，本文提出了基于剛體彈道方程的防空火控解算方法，仿真結(jié)果證明了該方法的有效性和可行性。然而，對于實(shí)際的防空系統(tǒng)，目標(biāo)的強(qiáng)機(jī)動特性使得基于傳統(tǒng)運(yùn)動模型的濾波算法出現(xiàn)精度下降甚至發(fā)散現(xiàn)象，因此后續(xù)的工作是考慮在該類情況下，提高目標(biāo)狀態(tài)估計(jì)與軌跡預(yù)測的準(zhǔn)確性，并在一定的基礎(chǔ)上降低算法的復(fù)雜度，使得該算法可用于未來防空武器的火控系統(tǒng)中。

表3 SylixOS下的解算時間對比Tab.3 Comparison of calculation times based on SylixOS ms