一種基于無跡卡爾曼濾波和炮彈運動軌跡模型的落/發(fā)點外推新算法

顏波濤，于 瓊，李同亮

(1.中國電子科技集團公司第三十八研究所，安徽 合肥 230031；2.中國人民解放軍陸軍工程大學(xué)，江蘇 南京 210001)

0 引 言

現(xiàn)代戰(zhàn)場的炮兵作戰(zhàn)，依然是一種強有力的作戰(zhàn)方式。炮位偵校雷達一方面具備偵察功能，通過對敵軍射出的飛行彈丸上升段進行搜索、跟蹤、彈道外推來獲取敵軍火炮陣地的位置，為我方武器對其進行打擊提供指導(dǎo)；另一方面，其具備校射功能，通過對我方射出彈丸的落點位置進行快速預(yù)測、比較偏差，從而校正火炮發(fā)射，很大程度上提高了我方火炮作戰(zhàn)性能[1-2]。

從炮位偵校雷達的工作過程來看，決定其效果的主要有3個階段：(1)目標(biāo)搜索確認(rèn)階段；(2)目標(biāo)跟蹤測量階段；(3)目標(biāo)數(shù)據(jù)處理階段，即彈丸落/發(fā)點外推。本文重點研究對測量數(shù)據(jù)的處理，即彈道外推算法。

現(xiàn)今，主要有2類彈道外推算法：一是基于先驗樣本的機器學(xué)習(xí)算法，其原理是利用彈道先驗數(shù)據(jù)對彈道模型進行統(tǒng)計建模，再對當(dāng)前彈道進行辨識。但是，炮位偵校雷達的實時性要求較高，且機器學(xué)習(xí)算法基于的大量準(zhǔn)確的先驗樣本一般難以提供，因此，一般認(rèn)為炮位偵校雷達的數(shù)據(jù)處理不宜采用機器學(xué)習(xí)算法；二是卡爾曼濾波估計算法，該類算法以彈道模型為基礎(chǔ)，先利用卡爾曼濾波方法從包含噪聲的量測數(shù)據(jù)中對系統(tǒng)參數(shù)進行辨識，再進行落/發(fā)點預(yù)測。其中較為常用的濾波算法有擴展卡爾曼濾波(EKF)、無跡卡爾曼濾波(UKF)等。應(yīng)雷達實際需求的發(fā)展，各種非線性濾波方法逐漸成為炮位偵校雷達數(shù)據(jù)處理的研究重點[3-5]。

根據(jù)作者文獻調(diào)研的結(jié)果，本文首先選擇了較為先進的無跡卡爾曼濾波(UKF)作為彈道外推新算法的濾波方式，減輕了量測噪聲的影響、辨識系統(tǒng)參數(shù)，再基于一種先進的炮彈運動軌跡模型，利用四階龍格-庫塔方程進行落/發(fā)點外推。期間，針對偵察模式加入了反向濾波和“七態(tài)”等策略，進一步提升算法效果。仿真測試結(jié)果顯示，相比于其它2種國內(nèi)雷達常用的經(jīng)典彈道外推算法，該算法的定位精度大大提高。

1 算法及測試

1.1 算法介紹

前文已述，所提出的算法可以分成2個部分，利用無跡卡爾曼濾波方法(UKF)對雷達量測數(shù)據(jù)進行濾波，再利用濾波結(jié)果，基于一種先進的炮彈運動軌跡模型和四階龍格-庫塔方程進行外推，獲得彈丸落/發(fā)點坐標(biāo)。針對偵察模式，附加了反向濾波和“七態(tài)”等策略。下文將對各個部分做詳細介紹。

1.1.1 濾波模型

1.1.1.1 炮彈運動彈道模型

(1) 校射所用彈道方程

校射時，雷達跟蹤彈道降弧段。由于彈丸的參數(shù)均已知，因而所采用的彈道模型擬考慮利用彈丸自身的氣動力系數(shù)，同時要涉及非標(biāo)準(zhǔn)條件(包括非標(biāo)準(zhǔn)彈道條件、非標(biāo)準(zhǔn)氣象條件)。但是，顧及到計算時間，彈道方程又不能太復(fù)雜。對彈道降弧段，也沒有必要采用6自由度剛體彈道方程，不過對偏流則應(yīng)適當(dāng)計及[6-8]。

以x表示射距離，y為高度，z為側(cè)偏，非標(biāo)準(zhǔn)條件的彈丸運動方程為：

(1)

(2)

式中：vr是相對速度，其表達式為：

(3)

Cs是聲速，其表達式為：

(4)

式中：G=6.328×10-3(K/m)；τ0=273.15+t0，t0是地面氣溫(℃)。

ρ是空氣密度，當(dāng)y<10 km時，可取下列公式：

ρ=ρ0exp(-β′y)

(5)

(6)

地面密度ρ0與地面氣壓p0的關(guān)系由下式?jīng)Q定：

(7)

式中：M為空氣1 kmol的質(zhì)量，為28.964 4 kg/kmol；R為普適氣體恒量，為8.314 32×103 J/(kmol·K)。

(8)

(2) 偵察所用彈道方程

由于雷達跟蹤仰角不會小于2°～3°，因此即使敵方發(fā)射戰(zhàn)術(shù)火箭，其彈道也將超過主動段而屬彈道被動段。與校射彈道模型比較而言，其主要不同點在于：(a)偵察中的彈箭參數(shù)是未知量；(b)偵察彈道升弧段的前部分接近直線，偏流很小，可略去不計。由此仍以式(1)為基礎(chǔ)，并做適當(dāng)變更。首先令：

(9)

(10)

式中：C為彈道系數(shù)；Cxon(Ma)為標(biāo)準(zhǔn)彈的阻力系數(shù)，通常稱之為阻力定律。

在彈道計算中，彈道系數(shù)對給定射角為常量，那么，N與N1所依賴的狀態(tài)變量之間的函數(shù)關(guān)系是相同的。

其次，設(shè)置式(1)中的Kz=0。

最后，增設(shè)狀態(tài)變量C，對彈道系數(shù)進行估計，得到偵察彈道的狀態(tài)方程。

1.1.1.2 無跡卡爾曼濾波

無跡卡爾曼濾波(UKF)的產(chǎn)生和發(fā)展是基于不敏變換。其原理是，首先利用1組確定的樣本點對狀態(tài)向量分布進行反映，然后利用任意非線性函數(shù)進行轉(zhuǎn)換后，這些樣本點依然能夠良好反映函數(shù)分布，通過整個過程對函數(shù)真實值精確估計，精度幾乎可以達到2階以上。UKF相對于一階擴展卡爾曼濾波(EKF)具有更高的精度，而且避開雅可比矩陣的計算，不過相比EKF計算量仍有所增加[7-9]。

在外推算法的數(shù)據(jù)處理中，校射模式下，由于彈道系數(shù)已知，估計的狀態(tài)向量維數(shù)為6，即六態(tài)，包括位置和速度；偵察模式下，增加了對彈道系數(shù)C的估計，即七態(tài)。

UKF的處理過程：

(1) 選取不敏變換采樣點，以標(biāo)準(zhǔn)對稱采樣為例。在已知非線性函數(shù)的均值Xk-1|k-1和協(xié)方差Pk-1|k-1時，選擇2L+1個采樣點(L為狀態(tài)向量的維度)進行下列計算:

χ0=Xk-1|k-1

(11)

(12)

(13)

(14)

(15)

(16)

(17)

(2) 計算狀態(tài)預(yù)測

(18)

χk|k-1=f(Xk|k-1,k-1)

(19)

(3) 計算預(yù)測協(xié)方差矩陣

[χi,k|k-1-Xk|k-1]T+Qk-1

(20)

(4) 計算量測預(yù)測

(21)

ξk|k-1=h(χk|k-1,k)

(22)

在上述步驟完成后，即可完成濾波估計。UKF的過程通常是逐漸收斂的，也即是說，一般濾波末點的精確度最高，因此，彈道外推算法通常選擇濾波終點作為外推起點。

1.1.2 彈道解算

炮位偵校雷達的彈道外推算法包括2個步驟：(1)對跟蹤目標(biāo)的量測數(shù)據(jù)進行非線性濾波，辨識彈道參數(shù)；(2)沿彈道起點方向進行解算，當(dāng)解算至既定高度時結(jié)束計算，獲取外推結(jié)果[9-10]。

目前較為廣泛的解算方法是通過既定彈道模型建立彈道方程組，再采用數(shù)值計算方法解算彈道，最后推算落/發(fā)點。

解彈道微分方程組的方法很多，目前工程上應(yīng)用較為廣泛的是4階龍格-庫塔法[11]。

龍格-庫塔解法是一種以函數(shù)y(x)的泰勒級數(shù)分解為基礎(chǔ)的解算方法，形式如下，首先建立微分方程組并獲取初值:

(23)

yi(t0)=yi0(i=1,2,…,m)

(24)

假設(shè)方程組在某時間點n處各變量的值為(t,y1n,y2n,…,ymn)，則下一時間點n+1處，可通過4階龍格-庫塔算法計算各變量的值：

(25)

(26)

式中：h為迭代時間間隔。

采用此算法求解彈道落/發(fā)點的迭代過程如下：

步驟1：令n=0,t0=0,濾波末點狀態(tài)向量作為yi,0的初值代入；

步驟2：利用tn和yi,n，計算過程值Ki1，Ki2，Ki3，Ki4；

步驟3：計算

yi,n+1=yi,n+(Ki1+2Ki2+2Ki3+Ki4)/6

(27)

步驟4：判斷循環(huán)是否可結(jié)束，若否，令n=n+1，返回步驟2；若是，輸出結(jié)果，結(jié)束。

1.1.3 偵察模式下的附加策略

1.1.3.1 “七態(tài)”向量

在所有的彈道外推算法中，一個影響彈道外推的重要參數(shù)就是彈道系數(shù)。校射模式下，彈道系數(shù)為已知，濾波算法通常構(gòu)造六態(tài)狀態(tài)向量，包括彈丸的位置坐標(biāo)和速度，再以此建立狀態(tài)方程，這可以認(rèn)為是合理的。偵察模式下，敵方炮彈的彈道系數(shù)顯然未知，于是出現(xiàn)了2種做法:一是算法認(rèn)為彈道系數(shù)是常數(shù)，以相鄰兩點的速度、加速度濾波值代入彈道方程進行彈道系數(shù)粗略估計，實驗證明誤差較大；二是另一種彈道系數(shù)估計方法，即通過馬赫數(shù)、阻力系數(shù)、彈丸速度等參數(shù)求取彈道系數(shù)，不過這種方法對先驗知識要求極高，實際中難以使用[12-14]。

有文獻提出，可將該參數(shù)納入濾波過程，不斷迭代，最終獲取較為精確的參數(shù)估計結(jié)果。選取彈丸位置(x,y,z)及其在三維坐標(biāo)軸上的速度分量(VX,VY,VZ)和彈道系數(shù)C形成7個參數(shù)的狀態(tài)向量，簡稱“七態(tài)”：

X=(x1,x2,x3,x4,x5,x6,x7)=

(x,y,z,vx,vy,vz,C)

(28)

此時，彈道估計狀態(tài)方程為：

(29)

1.1.3.2 反向濾波

反向濾波的原理在于減小外推距離。很多彈道外推算法，在偵察模式下進行濾波時，依然選擇濾波時沿時間軸的正方向，即起始波束Ss至截止波束Se，外推時從濾波終點Se開始，按照估計的狀態(tài)向量進行外推[15-17]。但這種方法可以通過反向外推減小外推距離，從而提高外推精度。

為了提高偵察模式下的發(fā)點外推精度，有文獻提出，可以在濾波時選擇相反的方向。即，將截止波束Se作為濾波起點，起始波束Ss作為濾波終點。該濾波算法對時間“方向”取反，將跟蹤波束中的最后一點Ss作為外推時的起點，各坐標(biāo)軸的速度分量為(-VX,-VY,-VZ),于是狀態(tài)向量為X=[x,y,z,-Vx,-Vy,-Vz,C]。

實際證明，該方法的外推距離更短，誤差更小，可以加以利用[18-19]。

1.2 其它2種對比算法

本文選擇了其它2種較為常見的彈道外推算法作為對比，介紹如下。

算法1：國內(nèi)某雷達實際使用算法。彈道運動模型與本文一致，但是選擇的濾波方法為擴展卡爾曼濾波(EKF)；偵察模式下利用了反向濾波的思想，但是對彈道系數(shù)的估計是基于對現(xiàn)有彈道量測點擬合的思想計算得到的。

算法2：國內(nèi)某雷達實際使用算法。選擇的濾波方法與本文一致，即UKF，但是利用的是不一樣的彈道運動模型;偵察模式下也利用了反向濾波的思想，但是對彈道系數(shù)的處理是選取1個固定的經(jīng)驗值。

1.3 測試結(jié)果及分析

外推算法對比試驗選取同一炮彈的雷達量測軌跡，選取同一起點，分別截取的量測軌道點數(shù)為50、100、150、200，每次外推次數(shù)為100。校射和偵察誤差實驗結(jié)果分別如圖1和圖2所示。

圖1 3種算法校射誤差

圖2 3種算法偵察誤差

顯然，本文算法的外推誤差要明顯小于另外2種算法，且原始軌跡采樣點數(shù)越多，即外推起始點越接近外推目標(biāo)點，外推誤差越小。

2 結(jié)束語

炮彈彈道落發(fā)點外推算法的好壞是決定炮位偵察校射雷達的關(guān)鍵因素之一。本文在對文獻梳理和現(xiàn)有算法仿真研究的基礎(chǔ)上，提出了一種基于無跡卡爾曼濾波和炮彈運動軌跡模型的落發(fā)點外推新算法。在所有算法中，該算法首次嘗試將這兩者結(jié)合起來，并加入了其它一些先進的外推策略，如“七態(tài)”和反向濾波。經(jīng)過大量仿真和雷達實測數(shù)據(jù)的測試，結(jié)果顯示，相比于國內(nèi)其它2種常用的彈道外推算法，該算法的外推效果得到了大大提高。