代 霞，方 睿，靳俊峰，劉軍偉

(1.中國(guó)電子科技集團(tuán)公司第三十八研究所，安徽 合肥 230088；2.孔徑陣列與空間探測(cè)安徽省重點(diǎn)實(shí)驗(yàn)室，安徽 合肥 230088)

0 引 言

隨著現(xiàn)代科學(xué)技術(shù)的不斷進(jìn)步，新型戰(zhàn)爭(zhēng)武器層出不窮。其中，彈道導(dǎo)彈被譽(yù)為“殺手锏”武器，自誕生之初就受到了各國(guó)的青睞[1-2]。應(yīng)運(yùn)而生地，導(dǎo)彈防御成為各國(guó)軍事國(guó)防領(lǐng)域亟待攻破的堡壘之一[3-4]。

彈道導(dǎo)彈落點(diǎn)外推精度是衡量導(dǎo)彈防御能力的一個(gè)重要指標(biāo)?？焖佟⒕_地預(yù)報(bào)導(dǎo)彈落點(diǎn)是反導(dǎo)雷達(dá)研究的重點(diǎn)課題之一，其研究核心一方面在于雷達(dá)硬件測(cè)量精度的提高，另一方面在于軟件系統(tǒng)落點(diǎn)預(yù)測(cè)算法的優(yōu)化。

當(dāng)前比較主流的落點(diǎn)推算方法有2種：一種是數(shù)值計(jì)算法，建立彈道導(dǎo)彈動(dòng)力學(xué)方程，采用Runge-Kutta數(shù)值積分方法，從初始點(diǎn)開始，一步一步外推到落點(diǎn)[5]；另一種是軌道解析法，利用一個(gè)或多個(gè)初始點(diǎn)的位置、速度信息，采用幾何方法，計(jì)算出導(dǎo)彈軌道參數(shù)信息，再通過軌道參數(shù)計(jì)算導(dǎo)彈落點(diǎn)[6]。

在選取足夠小的積分步長(zhǎng)的情況下，數(shù)值計(jì)算法可以達(dá)到很高的精度，但計(jì)算時(shí)間較長(zhǎng)，耗費(fèi)資源較多；軌道解析法計(jì)算速度很快，但由于忽略了大氣層內(nèi)空氣阻力對(duì)導(dǎo)彈的作用，越靠近地球表面，空氣阻力越大，因此該方法外推精度較數(shù)值計(jì)算法要低。當(dāng)前計(jì)算機(jī)技術(shù)迅猛發(fā)展，計(jì)算資源已經(jīng)不足以成為雷達(dá)軟件系統(tǒng)發(fā)展的制約因素，因此數(shù)值計(jì)算法成為彈道導(dǎo)彈落點(diǎn)外推的優(yōu)選方法。

在反導(dǎo)雷達(dá)工作過程中，落點(diǎn)外推是伴隨著點(diǎn)跡更新實(shí)時(shí)進(jìn)行的。雷達(dá)點(diǎn)跡經(jīng)凝聚、起始、關(guān)聯(lián)后形成航跡，對(duì)航跡進(jìn)行濾波，最后選取合適的濾波點(diǎn)作為初始點(diǎn)進(jìn)行彈道外推。外推的精度在很大程度上取決于初始點(diǎn)信息的準(zhǔn)確性。

初始點(diǎn)信息的誤差主要來自兩方面：一方面是雷達(dá)測(cè)量誤差，這是由硬件決定的，不可避免[7]；另一方面來源于濾波算法的誤差，常用的濾波算法包括最小二乘法濾波和卡爾曼濾波等等。

理論和經(jīng)驗(yàn)都表明，隨著觀測(cè)時(shí)間的推移，航跡關(guān)聯(lián)點(diǎn)跡越來越多，彈道外推精度越來越高，因此實(shí)時(shí)外推時(shí)大多選取航跡最新關(guān)聯(lián)點(diǎn)作為初始點(diǎn)。最小二乘法濾波對(duì)于中點(diǎn)的平滑效果最好，在需要進(jìn)行高精度實(shí)時(shí)外推時(shí)不適用；而卡爾曼濾波器隨著關(guān)聯(lián)點(diǎn)跡增多越來越穩(wěn)定[8]。

針對(duì)這兩方面的誤差，本文提出了一種基于平滑優(yōu)化思想的彈道導(dǎo)彈落點(diǎn)外推方法，將不敏卡爾曼濾波(UKF)[9]與不敏卡爾曼平滑方法(UKS)相結(jié)合，先對(duì)關(guān)聯(lián)點(diǎn)跡按時(shí)間先后順序進(jìn)行UKF濾波，然后對(duì)濾波值進(jìn)行UKS平滑，之后將平滑值按時(shí)間從后往前反向排序，對(duì)其進(jìn)行二次UKS平滑，取時(shí)間最新的點(diǎn)為初始點(diǎn)進(jìn)行落點(diǎn)外推，最后計(jì)算最新幾次外推落點(diǎn)的平均值，作為最終的結(jié)果。

利用STK軟件進(jìn)行彈道仿真試驗(yàn)，選取最新關(guān)聯(lián)點(diǎn)作為初始點(diǎn)，分別使用濾波信息和本文提出的方法進(jìn)行外推，多次試驗(yàn)結(jié)果均表明，本文提出的方法能大幅度提高落點(diǎn)外推的精度。

1 基于平滑優(yōu)化的彈道導(dǎo)彈落點(diǎn)外推方法

UKF是在卡爾曼濾波(KF)的基礎(chǔ)上發(fā)展起來的一種非線性濾波方法，它利用有限的參數(shù)近似反映隨機(jī)量的統(tǒng)計(jì)特性，然后用加權(quán)統(tǒng)計(jì)線性回歸的方法來估計(jì)隨機(jī)量的均值和協(xié)方差。

UKF按時(shí)間順序從前向后濾波，利用時(shí)間靠前的點(diǎn)跡信息對(duì)時(shí)間靠后的點(diǎn)跡信息進(jìn)行優(yōu)化更新，減小時(shí)間靠后的點(diǎn)跡誤差。從理論上來說，也可以利用時(shí)間靠后的點(diǎn)跡信息對(duì)時(shí)間靠前的點(diǎn)跡進(jìn)行優(yōu)化，從而減小時(shí)間靠前的點(diǎn)跡誤差，UKS正是這種思想的實(shí)現(xiàn)。

1.1 UKF

UKF濾波器的基礎(chǔ)是不敏變換(UT)。不敏變換以雷達(dá)量測(cè)為基礎(chǔ)，選取2n+1個(gè)δ采樣點(diǎn)(其中n為狀態(tài)向量的維數(shù))，并賦予不同的權(quán)值W，通過狀態(tài)方程計(jì)算量測(cè)向量的估計(jì)均值和協(xié)方差。

UKF濾波的矩陣表示如下：

(1) 預(yù)測(cè)

(1)

(2)

(3)

(4)

(2) 更新

(5)

(6)

(7)

(8)

(9)

由此，計(jì)算得濾波增量為：

(10)

狀態(tài)向量和狀態(tài)協(xié)方差更新分別為：

(11)

(12)

1.2 UKS

UKS以UKF計(jì)算得到的狀態(tài)向量和狀態(tài)協(xié)方差更新為基礎(chǔ)，采用和UKF類似的UT方法，結(jié)合狀態(tài)方程和量測(cè)方程進(jìn)行計(jì)算，其第一步和UKF類似：

(13)

(14)

(15)

(16)

(17)

由此計(jì)算得平滑增量為：

(18)

狀態(tài)向量和狀態(tài)協(xié)方差平滑值分別為：

(19)

(20)

用STK軟件仿真生成一條理論彈道，時(shí)長(zhǎng)為100 s，數(shù)據(jù)率為1 Hz，在距離、方位、仰角3個(gè)方向上各加上高斯噪聲，形成包含100個(gè)點(diǎn)的航跡，以航跡點(diǎn)為基礎(chǔ)，按時(shí)間先后順序進(jìn)行UKF濾波，之后以濾波值為基礎(chǔ)，進(jìn)行一次UKS平滑。

圖1是航跡濾波值和平滑值相對(duì)理論值的誤差。從圖中可以看到，無論是從位置誤差還是速度誤差來看，UKF濾波之后，航跡誤差相對(duì)量測(cè)大幅減小，但在量測(cè)點(diǎn)誤差出現(xiàn)峰值的地方依然存在較大的誤差，經(jīng)過一次UKS之后，整條航跡變得更平滑，更接近理論航跡。

圖1 量測(cè)值、濾波值和平滑值誤差對(duì)比圖

1.3 基于平滑優(yōu)化的彈道導(dǎo)彈落點(diǎn)外推方法

UKS是從時(shí)間最新關(guān)聯(lián)點(diǎn)向前進(jìn)行的，以最新關(guān)聯(lián)點(diǎn)進(jìn)行初始化，因此沒有對(duì)最新關(guān)聯(lián)點(diǎn)做平滑，其狀態(tài)向量和協(xié)方差與UKF計(jì)算結(jié)果相同，如圖1中插入的小圖所示。為了減小最新關(guān)聯(lián)點(diǎn)誤差，提高外推精度，本文在一次UKS的基礎(chǔ)上進(jìn)行了二次平滑。

第2次UKS方向和第1次UKS方向相反，從時(shí)間靠前的關(guān)聯(lián)點(diǎn)開始，從前向后進(jìn)行平滑，將第1次UKS優(yōu)化后的信息再次傳遞到最新關(guān)聯(lián)點(diǎn)，減小最新關(guān)聯(lián)點(diǎn)誤差。以更新后的最新關(guān)聯(lián)點(diǎn)作為初始值，進(jìn)行落點(diǎn)外推。圖1畫出了2次UKS誤差的對(duì)比圖，可以看到，經(jīng)過第2次UKS之后，最新關(guān)聯(lián)點(diǎn)位置和速度誤差均有明顯減小。

在雷達(dá)工作過程中，航跡長(zhǎng)度是不確定的，如果外推時(shí)航跡包含n個(gè)點(diǎn)，則2次平滑一共需要處理2n個(gè)點(diǎn)，當(dāng)n很大時(shí)，處理的時(shí)間會(huì)比較長(zhǎng)，不符合外推的實(shí)時(shí)性要求；另一方面，外推只需要高精度的最新關(guān)聯(lián)點(diǎn)的信息，作為第1次UKS的首點(diǎn)和第2次UKS的末點(diǎn)，航跡長(zhǎng)度不影響平滑過程的的進(jìn)行。因此，外推時(shí)設(shè)置合適的窗口寬度，只對(duì)窗口內(nèi)的點(diǎn)進(jìn)行平滑處理，如圖2所示。

圖2 基于平滑優(yōu)化的彈道導(dǎo)彈落點(diǎn)外推方法示意圖

利用STK軟件進(jìn)行彈道仿真試驗(yàn)，選取不同的航跡長(zhǎng)度，以最新關(guān)聯(lián)點(diǎn)作為初始點(diǎn)，采用數(shù)值計(jì)算法進(jìn)行外推。

圖3是外推落點(diǎn)誤差隨平滑窗口寬度的變化圖。從圖3可以看到，采用本文方法進(jìn)行落點(diǎn)外推的誤差小于直接采用濾波結(jié)果進(jìn)行外推(平滑窗口寬度為0)的誤差；隨著平滑窗口寬度的增大，落點(diǎn)誤差逐漸減??；當(dāng)平滑窗口寬度小于5個(gè)點(diǎn)時(shí)，落點(diǎn)誤差和采用濾波結(jié)果進(jìn)行外推的誤差差別很小，這是因?yàn)閁KS和UKF類似，幾步平滑之后才能穩(wěn)定；當(dāng)窗口寬度增大到25時(shí)，落點(diǎn)誤差已經(jīng)趨于穩(wěn)定；再進(jìn)一步增大窗口寬度，落點(diǎn)誤差變化很小，基本保持穩(wěn)定。因此，在本文的外推示例中，選取25個(gè)點(diǎn)作為2次UKS的窗口寬度。

圖3 外推落點(diǎn)誤差隨平滑窗口寬度變化圖

在雷達(dá)實(shí)際工作過程中，每隔一定時(shí)間以最新關(guān)聯(lián)點(diǎn)作為初始點(diǎn)進(jìn)行一次落點(diǎn)外推。由于量測(cè)誤差的存在，最新關(guān)聯(lián)點(diǎn)可能存在較大的偶然誤差。為了降低偶然誤差對(duì)落點(diǎn)外推精度的影響，計(jì)算最新若干(記為p)次外推落點(diǎn)的幾何平均位置，作為最終的落點(diǎn)。

綜上所述，基于平滑優(yōu)化的彈道導(dǎo)彈落點(diǎn)外推方法分為以下幾個(gè)步驟：

(1) 將航跡點(diǎn)按時(shí)間先后順序進(jìn)行UKF，假設(shè)每隔Δt時(shí)間進(jìn)行一次落點(diǎn)外推；

(2) 外推時(shí)間到，當(dāng)前航跡長(zhǎng)度記為n，平滑窗口寬度為w。如果n≥w，則取第(n-w+1)個(gè)點(diǎn)到第n個(gè)點(diǎn)之間(包含2個(gè)端點(diǎn))的航跡進(jìn)行UKS；如果n

(3) 將上一步中經(jīng)第1次UKS優(yōu)化之后的航跡按時(shí)間反向順序，進(jìn)行第2次UKS；

(4) 以第2次UKS優(yōu)化后的最新關(guān)聯(lián)點(diǎn)作為初始點(diǎn)，采用數(shù)值計(jì)算法進(jìn)行落點(diǎn)外推；

(5) 假設(shè)當(dāng)前時(shí)刻及之前外推總次數(shù)為m，如果m≥p，則計(jì)算第m-p+1次外推落點(diǎn)到當(dāng)前時(shí)刻外推落點(diǎn)的幾何平均值；如果m

2 仿真結(jié)果

本文采用STK軟件進(jìn)行仿真試驗(yàn)[10]，模擬了一條射程為1 000 km的標(biāo)準(zhǔn)彈道，發(fā)射時(shí)間為2018年1月1日4點(diǎn)30分0秒，彈道生成函數(shù)為Ballistic。導(dǎo)彈發(fā)射點(diǎn)經(jīng)度為127.433°，緯度為39.166 7°，高度為0 m；導(dǎo)彈落點(diǎn)經(jīng)度為136.452°，緯度為33.843°，高度為0 m。導(dǎo)彈飛行時(shí)間為523.7 s，最大飛行高度為284.7 km。雷達(dá)站中心經(jīng)度為136.452°，緯度為33.843°，高度為0 m，數(shù)據(jù)率為1 Hz。雷達(dá)測(cè)量誤差滿足高斯分布，測(cè)量精度：距離標(biāo)準(zhǔn)差為100 m，方位角標(biāo)準(zhǔn)差為0.3°，俯仰角標(biāo)準(zhǔn)差為0.3°。

試驗(yàn)中，對(duì)理論彈道在斜距、方位、仰角3個(gè)緯度分別單獨(dú)進(jìn)行加噪，生成了50條模擬彈道，計(jì)算過程中取50條彈道計(jì)算結(jié)果的平均值。模擬過程中從第30 s開始，到第500 s結(jié)束，每隔5 s進(jìn)行一次外推，UKS窗口為25個(gè)點(diǎn)，對(duì)落點(diǎn)求平均的窗口寬度為12次外推，分別采用UKF結(jié)果和本文中的方法進(jìn)行落點(diǎn)外推，圖4為外推落點(diǎn)誤差隨時(shí)間變化圖。

圖4 彈道外推落點(diǎn)誤差隨外推時(shí)間變化圖

從圖4中可以看到，總體看來，無論使用哪種外推方法，隨著時(shí)間的推進(jìn)，外推落點(diǎn)誤差均逐漸減小，和理論分析一致。除了第35 s進(jìn)行的外推，每一次外推中，采用2次平滑之后的最新關(guān)聯(lián)點(diǎn)作為初始值進(jìn)行外推的落點(diǎn)誤差均小于直接采用濾波結(jié)果進(jìn)行外推的落點(diǎn)誤差。

第405 s進(jìn)行的外推落點(diǎn)誤差減小最多，達(dá)到23.67%，落點(diǎn)誤差平均減小17.16%；第35 s，采用濾波值進(jìn)行外推的落點(diǎn)誤差為74.68 km，采用二次平滑值進(jìn)行外推的結(jié)果比之增大5.97%，為79.14 km。

二次平滑外推后，對(duì)最新12次外推落點(diǎn)計(jì)算幾何平均值，其與理論落點(diǎn)的誤差如圖4所示。從圖中可以看到，除了最開始的幾次外推，總體看來，落點(diǎn)平均值誤差相對(duì)不做平均的結(jié)果有大幅減小，減小幅度最大為63.10%，發(fā)生在第395 s，平均減小48.00%；對(duì)于最開始的6次外推，單次外推落點(diǎn)誤差隨外推時(shí)間推進(jìn)急劇減小，從92.55 km減小到47.33 km；對(duì)外推落點(diǎn)計(jì)算幾何平均值后，誤差反而增大，這與前文分析一致，因此對(duì)于這段時(shí)間內(nèi)的外推結(jié)果計(jì)算平均值時(shí)窗口設(shè)為1次比較合適(不計(jì)算平均值)。

綜合以上分析結(jié)果可知，在進(jìn)行常規(guī)濾波后，再進(jìn)行2次方向相反的平滑，以2次平滑后的最新關(guān)聯(lián)點(diǎn)作為初始點(diǎn)進(jìn)行外推可有效提高彈道導(dǎo)彈落點(diǎn)外推的精度；每次外推后，以最近幾次外推落點(diǎn)的幾何平均值作為最終的結(jié)果，又可以顯著提高落點(diǎn)外推精度。

需要注意的是，在本文的彈道仿真試驗(yàn)中，假設(shè)雷達(dá)量測(cè)噪聲是符合正態(tài)分布的，而對(duì)于真實(shí)雷達(dá)系統(tǒng)來說，該假設(shè)并不一定成立。因此，在計(jì)算落點(diǎn)的平均值時(shí)，可以對(duì)每個(gè)落點(diǎn)賦予不同的權(quán)重，比如對(duì)最新外推落點(diǎn)賦予最大的權(quán)重，以降低系統(tǒng)誤差的非正態(tài)性對(duì)最終結(jié)果的影響。

3 結(jié)束語

彈道導(dǎo)彈落點(diǎn)外推的精度在很大程度上取決于初始點(diǎn)信息的準(zhǔn)確性，而初始點(diǎn)信息的誤差來源主要有兩方面：一是雷達(dá)測(cè)量誤差；二是濾波算法誤差。

本文針對(duì)這2類誤差，設(shè)計(jì)并實(shí)現(xiàn)了一種基于平滑優(yōu)化思想的彈道導(dǎo)彈落點(diǎn)外推方法。該方法在常規(guī)UKF的基礎(chǔ)上，設(shè)置合適的平滑窗口寬度，對(duì)窗口內(nèi)的濾波點(diǎn)進(jìn)行一次平滑優(yōu)化；之后再對(duì)優(yōu)化值按時(shí)間反向順序進(jìn)行二次平滑優(yōu)化，選取二次優(yōu)化后的最新關(guān)聯(lián)點(diǎn)作為初始點(diǎn)進(jìn)行外推；最后以最近若干次外推結(jié)果的幾何平均值作為本次外推的最終結(jié)果。利用STK軟件進(jìn)行彈道仿真試驗(yàn)，結(jié)果表明該方法可以有效提高落點(diǎn)外推的精度。

本文方法的核心在于對(duì)雷達(dá)量測(cè)進(jìn)行常規(guī)濾波后，再對(duì)濾波值進(jìn)行2次平滑。這種處理可以顯著提高平滑窗口內(nèi)航跡點(diǎn)信息的準(zhǔn)確性，不僅可以應(yīng)用于彈道外推，還特別適用于雷達(dá)事后數(shù)據(jù)處理。此外，在設(shè)置合適的平滑窗口寬度的情況下，該方法也可以應(yīng)用到雷達(dá)實(shí)時(shí)數(shù)據(jù)處理中。