曹樹新,李 巖,李帥孝,徐樹威

(1.南京理工大學(xué) 能源與動力工程學(xué)院,江蘇 南京 210094;2.陸軍裝備部 沈陽地區(qū)軍事代表局駐沈陽地區(qū)第二軍事代表室,遼寧 沈陽 110004;3.齊齊哈爾建華機(jī)械有限公司,黑龍江 齊齊哈爾 161006)

隨著快速機(jī)動能力的不斷增強(qiáng)和偵察校射技術(shù)的不斷發(fā)展,現(xiàn)代戰(zhàn)爭對火炮的快速反應(yīng)能力提出了更高要求,氣象保障作為火炮射擊保障的核心環(huán)節(jié)之一,有必要在提升氣象保障的速度和精度方面加以研究[1]。近年來,數(shù)值天氣預(yù)報在火炮氣象保障中的應(yīng)用極大縮短了氣象保障的準(zhǔn)備時間[2-5],但是,初始場測量處理誤差和數(shù)值模式建模誤差,使得數(shù)值預(yù)報結(jié)果中含有較大的預(yù)報誤差,且尤以風(fēng)的預(yù)報誤差最為顯著[6]?，F(xiàn)行有效的數(shù)值修正方法以數(shù)值天氣預(yù)報的預(yù)報誤差及其變化規(guī)律為依據(jù),對臨近天氣預(yù)報的預(yù)報誤差進(jìn)行數(shù)值修正[7],提升氣象要素的預(yù)報精度。因此,如何快速準(zhǔn)確獲取預(yù)報誤差已成為提升彈道解算用數(shù)值天氣預(yù)報的預(yù)報精度需首要解決的問題之一。

目前,利用成熟的火炮外彈道跟蹤觀測技術(shù)獲取的實測彈道數(shù)據(jù),為快速準(zhǔn)確獲取彈道解算用數(shù)值天氣預(yù)報的預(yù)報誤差提供了新途徑。本文為獲取風(fēng)的預(yù)報誤差,使用風(fēng)修正系數(shù)描述風(fēng)的預(yù)報誤差,提出一種采用無跡卡爾曼濾波算法,直接從實測彈道數(shù)據(jù)和原始預(yù)報氣象數(shù)據(jù)中辨識風(fēng)修正系數(shù)的方法,并通過數(shù)值比對和彈道重構(gòu),對風(fēng)修正系數(shù)的辨識精度進(jìn)行分析。

1 系統(tǒng)模型

1.1 WRF模式簡介及參數(shù)設(shè)置

WRF(weather research and forecast,WRF)模式是由美國多家科研機(jī)構(gòu)于2000年開發(fā)的新一代中尺度數(shù)值模式和數(shù)據(jù)同化系統(tǒng),其廣泛應(yīng)用于中小尺度到全球尺度的數(shù)值預(yù)報和大氣模擬,并可提供指定區(qū)域在未來時段內(nèi)的時空網(wǎng)格化氣象數(shù)據(jù)。

本文中WRF模式預(yù)報區(qū)域位于我國東北地區(qū)的松嫩平原,該區(qū)域為溫帶大陸性季風(fēng)氣候,太陽輻射充足。模式水平方向上采用單重網(wǎng)格結(jié)構(gòu),格點數(shù)為121×121,格點間距為4 km,垂直方向上分為81層,模式頂層氣壓為1 000 Pa,初始場使用NCEP發(fā)布的分辨率為0.25°×0.25°的再分析資料,積分步長為15 s,預(yù)報時長為4 h。

1.2 假設(shè)條件

影響彈丸飛行狀態(tài)的氣象要素包括氣壓、虛溫、縱風(fēng)和橫風(fēng),為了具體分析彈道解算用數(shù)值天氣預(yù)報的氣象保障效果,本文以WRF模式的短時預(yù)報為例,一方面使用均方根誤差、平均絕對誤差和誤差標(biāo)準(zhǔn)差[8]對WRF模式預(yù)報的氣象要素進(jìn)行檢驗,另一方面,基于外彈道仿真,分析WRF模式中單一氣象要素的預(yù)報誤差對彈丸飛行狀態(tài)的影響。其中,預(yù)報氣象數(shù)據(jù)和實測氣象數(shù)據(jù)描述的均為該地2019年12月14日12時的高空大氣狀態(tài),外彈道仿真則使用某型榴彈的六自由度剛體彈道模型,射角為51°,射向為100.487°,射程約為30 km。WRF模式預(yù)報氣象要素的檢驗結(jié)果如表1所示,WRF模式單一氣象要素的預(yù)報誤差對彈丸飛行狀態(tài)影響的分析結(jié)果如表2所示。

表1 WRF模式預(yù)報氣象要素檢驗結(jié)果

表2 彈丸飛行狀態(tài)偏差均值表

由表1和表2可得,相比于風(fēng)的預(yù)報誤差,氣壓和虛溫的預(yù)報誤差及其離散程度以及預(yù)報誤差引起的彈道偏差均相對較小,因此,本文只將風(fēng)修正系數(shù)作為待辨識參數(shù),并忽略氣壓和虛溫的預(yù)報誤差對辨識風(fēng)修正系數(shù)的影響,故本文做出如下假設(shè):

①忽略氣壓和虛溫的預(yù)報誤差對辨識風(fēng)修正系數(shù)的影響;

②風(fēng)修正系數(shù)在辨識區(qū)間內(nèi)為常值。

1.3 彈道模型

為了從實測彈道數(shù)據(jù)和原始預(yù)報氣象數(shù)據(jù)中辨識風(fēng)修正系數(shù),需構(gòu)建彈丸運動模型,建立實測彈道數(shù)據(jù)與氣象數(shù)據(jù)之間的聯(lián)系?？紤]到目前彈道跟蹤觀測數(shù)據(jù)主要以位置和速度為主,本文選用質(zhì)點彈道模型,此彈道模型在有風(fēng)條件下的彈道方程[9]為

(1)

式中:

x、y、z,vx、vy和vz分別為彈丸在發(fā)射坐標(biāo)系下的坐標(biāo)分量和速度分量;m,d,l和JC分別為彈丸質(zhì)量、彈丸直徑、彈丸長度和彈丸極轉(zhuǎn)動慣量;Cx,C′y,m′z和m′xz分別為空氣阻力系數(shù)、升力系數(shù)導(dǎo)數(shù)、靜力矩系數(shù)導(dǎo)數(shù)和極阻尼力矩系數(shù)導(dǎo)數(shù);p,Tv,wx,wz,vs,Rd和k分別為當(dāng)?shù)貧鈮?、?dāng)?shù)靥摐?、?dāng)?shù)乜v風(fēng)、當(dāng)?shù)貦M風(fēng)、當(dāng)?shù)芈曀?、干空氣氣體常數(shù)和比熱比;kL和kC分別為縱風(fēng)修正系數(shù)和橫風(fēng)修正系數(shù);g,v0,s和η分別為當(dāng)?shù)刂亓铀俣?、?biāo)定初速、彈道弧長和膛線纏度。

2 無跡卡爾曼濾波算法辨識風(fēng)修正系數(shù)

為了實現(xiàn)從非線性彈道模型中辨識風(fēng)修正系數(shù),本文將無跡卡爾曼濾波算法作為辨識算法。由卡爾曼濾波穩(wěn)定性原理[10]可知,如果濾波穩(wěn)定,估計值和估計均方誤差將隨濾波時間的增長逐漸不受濾波初值的影響,并最終收斂,實現(xiàn)無偏估計。但是,每個估計值和每個估計均方誤差的收斂速度不同,若濾波時間較短,可能存在估計值和估計均方誤差沒有完成收斂的情況。因此,本文將原始預(yù)報氣象數(shù)據(jù)作為彈道模型的基礎(chǔ)氣象條件,采用無跡卡爾曼濾波算法對實測彈道數(shù)據(jù)進(jìn)行迭代處理,并將每次迭代所得風(fēng)修正系數(shù)和狀態(tài)變量的估計均方誤差作為下次迭代的初值,直至滿足迭代終止條件,其流程如圖1所示。

圖1 無跡卡爾曼濾波算法辨識風(fēng)修正系數(shù)流程圖

由式(1)和實測彈道數(shù)據(jù)的種類可知,無跡卡爾曼濾波算法中的觀測變量為x、y、z,vx、vy、vz,狀態(tài)變量為x、y、z,vx、vy、vz、kL和kC,即:

Z=(xyzvxvyvz)T

(2)

X=(x1x2x3x4x5x6x7x8)T=

(xyzvxvyvzkLkC)T

(3)

將式(3)代入式(1),得:

(4)

結(jié)合圖1,無跡卡爾曼濾波算法辨識風(fēng)修正系數(shù)的主要步驟如下。

①第1步。選定初值。

②第2步。無跡卡爾曼濾波計算,詳細(xì)計算步驟可參考文獻(xiàn)[10]。

③第3步。判斷是否滿足迭代終止條件,若沒有滿足,則返回第2步。為了提高迭代效率,本文設(shè)置2個迭代終止條件,條件一為達(dá)到最大迭代次數(shù)nmax,條件二為連續(xù)2次迭代所得縱風(fēng)修正系數(shù)之差εL和橫風(fēng)修正系數(shù)之差εC均小于給定閾值。

為了提升辨識結(jié)果精度,減小包括氣動參數(shù)誤差和起始擾動誤差在內(nèi)的其他誤差對辨識結(jié)果的影響,本文將彈丸風(fēng)洞吹風(fēng)數(shù)據(jù)作為彈道模型中氣動參數(shù)的真值,并將多組辨識結(jié)果的平均值作為最終的辨識結(jié)果。

3 數(shù)值計算與分析

為了檢驗無跡卡爾曼濾波算法辨識風(fēng)修正系數(shù)的效果,本文基于某實測彈道數(shù)據(jù)和WRF模式原始預(yù)報氣象數(shù)據(jù)對其進(jìn)行數(shù)值分析,相關(guān)參數(shù)設(shè)計如下。

①實測彈道數(shù)據(jù)的位置誤差(DRMS,68%):水平10 m,垂直10 m;速度精度(DRMS,68%):水平0.2 m/s,垂直0.5 m/s;

②彈道數(shù)據(jù)選擇距地面5 200～6 800 m的降弧段彈道數(shù)據(jù);

③WRF模式參數(shù)設(shè)置與第1.1節(jié)中的參數(shù)設(shè)置相同;

④無跡卡爾曼濾波算法迭代終止條件一為最大迭代次數(shù)nmax=200,迭代終止條件二為連續(xù)2次迭代差值εL<0.000 1且εC<0.000 1。

基于上述參數(shù),本文以其中一組實測彈道數(shù)據(jù)為例,在辨識過程中x、y、z的估計均方誤差隨迭代次數(shù)n(取前10次,下同)的變化如圖2所示,vx、vy、vz的估計均方誤差隨迭代次數(shù)的變化如圖3所示,kL、kC的估計均方誤差隨迭代次數(shù)的變化如圖4所示,kL、kC隨迭代次數(shù)的變化如圖5所示。

圖2 迭代過程中的位置估計均方誤差

圖3 迭代過程中的速度估計均方誤差

圖4 迭代過程中的風(fēng)修正系數(shù)估計均方誤差

將3組辨識結(jié)果取平均值作為該辨識區(qū)段內(nèi)的辨識結(jié)果,即辨識的縱風(fēng)修正系數(shù)為2.084 056,橫風(fēng)修正系數(shù)為1.254 748。實測的風(fēng)修正系數(shù)如表3所示,使用加權(quán)平均值(權(quán)重由彈丸在每個高度段內(nèi)的飛行時間確定,下同)代表實測的風(fēng)修正系數(shù),即實測的縱風(fēng)修正系數(shù)為1.944 045,橫風(fēng)修正系數(shù)為1.206 837。

表3 實測的風(fēng)修正系數(shù)

將實測的風(fēng)修正系數(shù)和辨識的風(fēng)修正系數(shù)進(jìn)行對比可得,辨識的縱風(fēng)修正系數(shù)和橫風(fēng)修正系數(shù)的相對誤差分別為7.20%和3.97%,產(chǎn)生上述誤差的原因可能如下:一是辨識過程中未考慮氣壓和虛溫的預(yù)報誤差對辨識風(fēng)修正系數(shù)的影響,將氣壓和虛溫的預(yù)報誤差折算入風(fēng)修正系數(shù)當(dāng)中;二是彈道模型誤差,結(jié)合六自由度剛體彈道模型可知,狀態(tài)變量還受彈丸姿態(tài)角等因素的影響,而質(zhì)點彈道模型未考慮這些因素;三是實測氣象數(shù)據(jù)與實測彈道數(shù)據(jù)時空不匹配。

原始預(yù)報氣象數(shù)據(jù)和實測氣象數(shù)據(jù)描述的均為同一時刻且同一地點的高空大氣狀態(tài),利用前文辨識所得的縱風(fēng)和橫風(fēng)的修正系數(shù)對原始預(yù)報氣象數(shù)據(jù)中的風(fēng)進(jìn)行修正,得到修正預(yù)報氣象數(shù)據(jù),3種氣象數(shù)據(jù)中風(fēng)的加權(quán)平均值如表4所示。

表4 彈道區(qū)段內(nèi)風(fēng)的加權(quán)平均值

由表4可得,原始預(yù)報氣象數(shù)據(jù)中的縱風(fēng)和橫風(fēng)經(jīng)修正后,其誤差分別減小85.19%和79.27%。

分別使用原始預(yù)報氣象數(shù)據(jù)和修正預(yù)報氣象數(shù)據(jù)對位于距地面5 200～6 800 m的降弧段彈道進(jìn)行重構(gòu),2種重構(gòu)的彈道段和濾波處理的實測彈道段在x-y平面的投影如圖6所示,在x-z平面的投影如圖7所示,在z-y平面的投影如圖8所示,2種重構(gòu)的彈道段與濾波處理的實測彈道段間的狀態(tài)均方根誤差如表5所示。

圖6 重構(gòu)彈道段在x-y平面的投影

圖7 重構(gòu)彈道段在x-z平面的投影

圖8 重構(gòu)彈道段在z-y平面的投影

表5 重構(gòu)彈道段的狀態(tài)均方根誤差

由圖6～圖8及表5可得,修正預(yù)報氣象重構(gòu)的彈道段相比原始預(yù)報氣象重構(gòu)的彈道段而言,更加接近于濾波處理的實測彈道段,因此,修正預(yù)報氣象比原始預(yù)報氣象對真實氣象的描述更為準(zhǔn)確,從而說明辨識所得風(fēng)修正系數(shù)的精度較高以及修正彈道解算用數(shù)值天氣預(yù)報的預(yù)報誤差是必要的。

4 結(jié)束語

本文將原始預(yù)報氣象數(shù)據(jù)作為基礎(chǔ)氣象條件,采用無跡卡爾曼濾波算法從實測彈道數(shù)據(jù)中辨識風(fēng)修正系數(shù),通過對實際辨識結(jié)果進(jìn)行分析,得出以下結(jié)論:該方法辨識所得風(fēng)修正系數(shù)具有較高的精度,可作為彈道解算用數(shù)值天氣預(yù)報的預(yù)報誤差,實現(xiàn)對臨近天氣預(yù)報的數(shù)值修正,滿足工程應(yīng)用需求;無跡卡爾曼濾波算法在辨識風(fēng)修正系數(shù)的過程中兼有濾波功能,能夠在一定程度上減小噪聲對辨識結(jié)果的影響。本文研究結(jié)果可以為提升彈道解算用數(shù)值天氣預(yù)報的預(yù)報質(zhì)量提供一定的理論參考。