顏波濤,于 瓊,李同亮
(1.中國電子科技集團公司第三十八研究所,安徽 合肥 230031;2.中國人民解放軍陸軍工程大學(xué),江蘇 南京 210001)
現(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)典彈道外推算法,該算法的定位精度大大提高。
前文已述,所提出的算法可以分成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]。
本文選擇了其它2種較為常見的彈道外推算法作為對比,介紹如下。
算法1:國內(nèi)某雷達實際使用算法。彈道運動模型與本文一致,但是選擇的濾波方法為擴展卡爾曼濾波(EKF);偵察模式下利用了反向濾波的思想,但是對彈道系數(shù)的估計是基于對現(xiàn)有彈道量測點擬合的思想計算得到的。
算法2:國內(nèi)某雷達實際使用算法。選擇的濾波方法與本文一致,即UKF,但是利用的是不一樣的彈道運動模型;偵察模式下也利用了反向濾波的思想,但是對彈道系數(shù)的處理是選取1個固定的經(jīng)驗值。
外推算法對比試驗選取同一炮彈的雷達量測軌跡,選取同一起點,分別截取的量測軌道點數(shù)為50、100、150、200,每次外推次數(shù)為100。校射和偵察誤差實驗結(jié)果分別如圖1和圖2所示。
圖1 3種算法校射誤差
圖2 3種算法偵察誤差
顯然,本文算法的外推誤差要明顯小于另外2種算法,且原始軌跡采樣點數(shù)越多,即外推起始點越接近外推目標(biāo)點,外推誤差越小。
炮彈彈道落發(fā)點外推算法的好壞是決定炮位偵察校射雷達的關(guān)鍵因素之一。本文在對文獻梳理和現(xiàn)有算法仿真研究的基礎(chǔ)上,提出了一種基于無跡卡爾曼濾波和炮彈運動軌跡模型的落發(fā)點外推新算法。在所有算法中,該算法首次嘗試將這兩者結(jié)合起來,并加入了其它一些先進的外推策略,如“七態(tài)”和反向濾波。經(jīng)過大量仿真和雷達實測數(shù)據(jù)的測試,結(jié)果顯示,相比于國內(nèi)其它2種常用的彈道外推算法,該算法的外推效果得到了大大提高。