一種估計(jì)艦空導(dǎo)彈脫靶量的方法

裴 喆，史 震

（1.哈爾濱工程大學(xué) 自動(dòng)化學(xué)院，哈爾濱 150001；2.92941 部隊(duì)，遼寧 葫蘆島 125001）

在靶場進(jìn)行艦空導(dǎo)彈攔截靶標(biāo)（靶機(jī)或靶彈）飛行試驗(yàn)時(shí)，脫靶量是必須測量的參數(shù)，因?yàn)樗呛饬繉?dǎo)彈在彈道末端制導(dǎo)與控制性能的關(guān)鍵指標(biāo)。脫靶量測量一般采用兩種方法：一種是設(shè)法測量導(dǎo)彈和靶標(biāo)的外彈道數(shù)據(jù)，利用遭遇段的空間位置參數(shù)重建彈靶運(yùn)動(dòng)軌跡，外推計(jì)算脫靶量[1]。目前主要采用地面雷達(dá)或光電經(jīng)緯儀進(jìn)行外彈道測量，雷達(dá)測量精度較低，光電經(jīng)緯儀測量精度高但易受氣象條件影響。另一種則是設(shè)法測量遭遇段彈靶之間的多普勒頻率，根據(jù)脫靶量與多普勒頻率之間的數(shù)學(xué)關(guān)系估計(jì)脫靶量[2-3]。目前主要采用遙測或靶機(jī)專用脫靶量測量設(shè)備進(jìn)行多普勒頻率測量，靶機(jī)專用脫靶量測量設(shè)備精度較高，但重復(fù)測量成本高，且無法裝在靶彈上。遙測是靶場主要的內(nèi)彈道測量方法，技術(shù)成熟、工作可靠，遙測數(shù)據(jù)中包含引信探測的多普勒頻率。實(shí)際彈靶遭遇時(shí)，由于靶標(biāo)體目標(biāo)效應(yīng)的影響，多普勒頻率經(jīng)常包含許多孤立野值和成片野值，而如何去除成片野值在公開發(fā)表的文獻(xiàn)中還沒有見到通用有效的方法[4]。本文首先討論利用無野值的多普勒頻率估計(jì)脫靶量的方法，然后研究利用小波分解識別并去除野值的方法。

1 脫靶量數(shù)學(xué)模型

導(dǎo)彈和靶標(biāo)遭遇過程中，由于相對速度大，遭遇時(shí)間短，可以認(rèn)為相對速度保持不變。脫靶量數(shù)學(xué)模型如圖1所示，圖中rV為導(dǎo)彈相對靶標(biāo)的運(yùn)動(dòng)速度，ρ為脫靶量（ρ垂直于 Vr），ti為多普勒頻率采樣時(shí)間，tρ為脫靶時(shí)間。

圖1 脫靶量數(shù)學(xué)模型

對于不同的ti(i=1,2,3,…,n)，多普勒頻率Ed(i)與脫靶量ρ的數(shù)學(xué)關(guān)系可由式(1)表示[2]：

式中：λ為引信天線波長；Vr可由引信開機(jī)前的導(dǎo)引頭遙測數(shù)據(jù)得出。由于引戰(zhàn)配合需要，引信在 tρ到來前引爆戰(zhàn)斗部，tρ時(shí)沒有Ed(i)數(shù)據(jù)，tρ和ρ只能由已測的Ed(i)數(shù)據(jù)根據(jù)式(1)計(jì)算得到。實(shí)際遙測測量的多普勒頻率fd(i)可由式(2)表示：

式中：ε (i)表示隨機(jī)誤差；εy(i)表示野值。

利用fd(i)估計(jì)ρ 需要盡量消除 ε (i)和 εy(i)，從而得到與Ed(i)非常接近的曲線。

2 多普勒頻率數(shù)據(jù)僅含隨機(jī)誤差時(shí)的脫靶量估計(jì)

先討論fd(i)僅含隨機(jī)誤差ε (i)的情況，即εy(i)=0。多數(shù)情況下 ε (i)近似于零均值高斯白噪聲，一般采用擬合或?yàn)V波的方法去除 ε (i)。小波變換可同時(shí)進(jìn)行時(shí)域和頻域分析，具有可變的時(shí)頻分辨率，通過小波多尺度分解可將信號分解為反映整體趨勢的低頻分量和反映噪聲細(xì)節(jié)的高頻分量。

對fd(i)利用db4 小波進(jìn)行4層分解，得到低頻系數(shù)a1、a2、a3、a4和高頻系數(shù)d1、d2、d3、d4，再對各層系數(shù)進(jìn)行重構(gòu)。分析各層系數(shù)發(fā)現(xiàn)：從第1層到第4層，高頻系數(shù)將噪聲按頻率降低的順序分別提取出來；低頻系數(shù)所含噪聲越來越少，多普勒頻率趨勢越來越明顯。a1與fd(i)曲線如圖2所示，由于d1提取了白噪聲最尖銳的部分，a1不但沒有丟失Ed(i)的有用信息，而且比fd(i)光滑許多。對a1進(jìn)行最小二乘擬合，使其與擬合曲線fd' (i)殘差序列的平方和 S (i)最小，表達(dá)式如式(3)所示：

最后由fd' (i)根據(jù)式(1)計(jì)算出ρ。

圖2 無野值多普勒頻率a1系數(shù)的擬和結(jié)果

為驗(yàn)證方法的可行性和有效性，現(xiàn)利用蒙特卡洛方法進(jìn)行仿真。給定 Vr、tρ、ρ、ti的一組設(shè)定值，取 Vr=1 200 m/s、tρ=30.000 s、ρ=10 m、ti=29.992 s時(shí)引信起爆。由式(1)得到Fd(i)，再模擬生成 ε (i)～N(0，1 0002)，按式(2)合成fd(i)序列。然后對fd(i)序列的小波系數(shù)a1進(jìn)行最小二乘擬合，得到ρ和 tρ的估計(jì)值。由于Fd(i)都是ρ和 tρ的非線性函數(shù)，須采用非線性最小二乘擬合[5]。仿真1 000次的統(tǒng)計(jì)結(jié)果如表1所示，可以看出：

a）利用a1進(jìn)行擬合精度很高，其中一次的擬合曲線fd' (i)如圖2所示，fd' (i)與Fd(i)幾乎重合。ρ的估計(jì)均值與設(shè)定值僅差約0.01 m，估計(jì)標(biāo)準(zhǔn)差僅約為0.24 m。統(tǒng)計(jì)還表明ρ的估計(jì)值近似服從正態(tài)分布，因此ρ的估計(jì)精度能夠控制在（?0.75 m～+0.75 m）范圍內(nèi)，完全滿足測量要求。

b）由于對實(shí)際數(shù)據(jù)進(jìn)行擬合時(shí)事先并不知道真值，仿真中設(shè)定了3組不同的初值，結(jié)果表明擬合算法對不同初值有很好的收斂性。

因此，當(dāng)多普勒頻率數(shù)據(jù)僅含隨機(jī)誤差時(shí)，利用該擬合算法估計(jì)脫靶量是可行而且有效的。

表1 a1系數(shù)擬合統(tǒng)計(jì)結(jié)果

3 多普勒頻率數(shù)據(jù)有野值時(shí)的脫靶量估計(jì)

3.1 野值識別

式(2)中的εy(i)表示野值，一批數(shù)據(jù)中有部分?jǐn)?shù)據(jù)與其余數(shù)據(jù)相比明顯不一致的稱為野值，或稱離群值、異常值[6]，一般分為孤立型和成片型。由于彈靶遭遇時(shí)靶標(biāo)的體目標(biāo)效應(yīng)影響，實(shí)際的引信多普勒頻率數(shù)據(jù)經(jīng)常含有野值。直接對有野值fd(i)的a1系數(shù)進(jìn)行擬合，結(jié)果的局部圖如圖3所示。擬合曲線fd' (i)與Fd(i)偏離較大，ρ估計(jì)值也與設(shè)定值有較大偏差，因此必須對野值進(jìn)行識別、去除。

圖3 有野值多普勒頻率a1系數(shù)的擬和結(jié)果

通過隨機(jī)加入不同類型、不同數(shù)量、不同位置、不同幅值的野值 εy(i)對fd(i)進(jìn)行模擬，利用db4小波進(jìn)行4層分解，將重構(gòu)后的各層系數(shù)與無野值多普勒頻率小波分解后的各層系數(shù)對比分析發(fā)現(xiàn)：

a）孤立野值和白噪聲主要影響d1的幅值，因此a1能反映成片野值的輪廓，定位比較準(zhǔn)確，a1局部圖如圖3所示。而且a1與擬合曲線fd' (i)的殘差序列能夠反映成片野值的大小。

b）成片野值主要影響d2、d3，它們局部幅值的增大說明成片野值的存在。如果設(shè)置一定門限可用于成片野值的識別，且這種識別方法與擬合曲線fd' (i)無關(guān)。

由此可總結(jié)出fd(i)有無野值的識別準(zhǔn)則：

a）孤立野值和兩點(diǎn)成片野值可直接利用fd(i)與擬合曲線fd' (i)的殘差序列識別，門限值一般取

±3σ。[7]

b）成片野值可利用d2、d3的局部幅值識別，門限值可通過與無野值多普勒頻率的d2、d3幅值統(tǒng)計(jì)比較得到。

c）成片野值也可利用a1與擬合曲線fd' (i)的殘差序列識別，多次仿真發(fā)現(xiàn)無野值多普勒頻率的a1與擬和曲線fd' (i)的殘差序列中幾乎沒有連續(xù)三點(diǎn)數(shù)據(jù)幅值都超過1.5σ的。因此，識別兩點(diǎn)以上成片野值的門限值可取±1.5σ。

3.2 野值去除

3.2.1 野值去除順序

確定野值去除順序應(yīng)首先分析野值的特殊性，通過對不同類型、數(shù)量、位置、幅值的野值進(jìn)行模擬發(fā)現(xiàn)：若幅值相同、分布位置相同，則成片野值相比孤立野值對擬合影響大；若幅值相同、類型相同，則位于后半?yún)^(qū)的野值比前半?yún)^(qū)的對擬合影響大，這是因?yàn)槎嗥绽疹l率是非線性的，離靶標(biāo)較遠(yuǎn)時(shí)變化較緩，較近時(shí)急劇下降。另外，由于靶標(biāo)體目標(biāo)效應(yīng)逐漸增強(qiáng)，通常后半?yún)^(qū)野值數(shù)量更多，幅值更大，且起伏劇烈。由于這些特性，擬合曲線fd' (i)對Fd(i)有時(shí)有較大偏離（尤其后半?yún)^(qū)），如圖3所示。

野值有正有負(fù)、有大有小，導(dǎo)致fd' (i)對Fd(i)有不同的偏離方向和不同的偏離程度，因此野值去除順序必須保證fd' (i)向Fd(i)逐步收斂。偏離的方向有全區(qū)上升、全區(qū)下降、前降后升、前升后降4種情況。全區(qū)上升時(shí)必須先去除正野值，后去除負(fù)野值；全區(qū)下降時(shí)與之相反。由野值的特殊性決定，前降后升時(shí)必須先去除后半?yún)^(qū)正野值和前半?yún)^(qū)負(fù)野值，最后去除其他野值；前升后降時(shí)與之相反。為表征fd' (i)偏離Fd(i)的方向和程度，將a1與fd' (i)的殘差序列定義為偏離因子，因?yàn)樵撔蛄心軌蜉^準(zhǔn)確地反映成片野值的位置和幅值。偏離因子初始值為0，它的正負(fù)對應(yīng)野值的正負(fù)，它的大小對應(yīng)偏離的程度。將偏離因子大小根據(jù)成片野值不同的長度、幅值、位置進(jìn)行分檔，長度越長、幅值越大，則偏離因子越大，且后半?yún)^(qū)的偏離因子比前半?yún)^(qū)的大。再將殘差序列按時(shí)間等分為幾個(gè)（一般6個(gè)即可）分區(qū)，計(jì)算各分區(qū)偏離因子以及全區(qū)偏離因子總和，這樣就可以識別偏離方向及程度了。全區(qū)上升或下降表現(xiàn)為全區(qū)偏離因子總和較大，且各分區(qū)偏離因子大部分為正或大部分為負(fù)。前降后升或前升后降表現(xiàn)為全區(qū)偏離因子總和較小，且前半?yún)^(qū)的偏離因子與后半?yún)^(qū)的偏離因子正負(fù)相反。

由此可總結(jié)出野值去除的順序：先確定擬合曲線的偏離方向和程度；然后對比前后半?yún)^(qū)，若后半?yún)^(qū)野值多則先去除，少則全區(qū)同時(shí)去除；無論哪個(gè)區(qū)域都須先去除較長成片野值，再去除較短成片野值，最后去除孤立野值。

3.2.2 野值替換方法

野值替換是指在野值點(diǎn)fd(j)處的擬合曲線fd' (j)數(shù)據(jù)增加適當(dāng)?shù)恼{(diào)整值（一般為k?σ），生成新的fd(j)數(shù)據(jù)：

孤立野值替換直接使用式(4)。成片野值替換需調(diào)整為式(5)：

因?yàn)槿绻善爸档奶鎿Q數(shù)據(jù)相同會(huì)影響再次擬和的結(jié)果，所以應(yīng)在式(4)基礎(chǔ)上增加與成片野值相同長度的零均值高斯白噪聲ε (j,j+m)～N(0,σ2)，m+1為成片野值長度。

k的正負(fù)與野值的正負(fù)相反，野值為正時(shí)，調(diào)整值為負(fù)；野值為負(fù)時(shí)，調(diào)整值為正。k的大小取決于擬和曲線fd' (i)偏離Fd(i)的程度，可用仿真統(tǒng)計(jì)的方法確定k的取值范圍。

3.3 算法仿真

為驗(yàn)證野值識別和去除方法的有效性，并確定k的取值范圍，編制了仿真程序，程序流程分為：

第1步：模擬生成如式(2)所示的fd(i)數(shù)據(jù)。其中野值 εy(i)由人工隨機(jī)設(shè)置，成片野值長度最長為7點(diǎn)，幅值最大5倍σ，數(shù)量不大于全部數(shù)據(jù)的二分之一。

第2步：利用db4 小波對fd(i)進(jìn)行4層分解，并重構(gòu)各層系數(shù)。對a1進(jìn)行非線性最小二乘擬合，得到擬合曲線fd' (i)、擬合殘差序列及ρ的估計(jì)值。

第3步：按3.1節(jié)中的準(zhǔn)則對fd(i)進(jìn)行野值識別，若有野值，執(zhí)行第4步；若無野值，則程序結(jié)束，ρ即為所求脫靶量。

第4步：按3.2節(jié)中的方法去除部分野值，對野值替換后生成的新fd(i)數(shù)據(jù)循環(huán)執(zhí)行第2步。

通過設(shè)置不同野值進(jìn)行仿真后結(jié)果表明，k的取值范圍主要與野值分布有關(guān)。后半?yún)^(qū)野值較少時(shí)，全區(qū)替換時(shí)k 一般在0＜k≤1范圍內(nèi)取值；后半?yún)^(qū)野值較多時(shí)，后半?yún)^(qū)k的取值范圍能達(dá)到 0＜k＜2，將后半?yún)^(qū)大部分野值替換后，再在0＜k≤1范圍內(nèi)取值對全區(qū)替換即可。

為提高替換精度，可將 k?σ的大小按0.25·σ分檔，即0.25·σ、0.5·σ、0.75·σ、σ，具體選檔需根據(jù)擬合曲線偏離程度確定。

仿真還表明，如果野值數(shù)量小于全部數(shù)據(jù)的三分之一，后半?yún)^(qū)野值較少時(shí)，野值去除后的脫靶量估計(jì)精度能控制在（?1 m～+1 m）范圍內(nèi)；后半?yún)^(qū)野值較多時(shí)，也能控制在（?1.5 m～+1.5 m）范圍內(nèi)，仍能滿足測量要求，表明算法是合理有效的。但如果野值數(shù)量大于全部數(shù)據(jù)的三分之一，擬合曲線有時(shí)會(huì)發(fā)散。

4 結(jié)束語

本文研究了一種利用小波分解和非線性最小二乘擬和相結(jié)合估計(jì)艦空導(dǎo)彈脫靶量的方法。通過仿真驗(yàn)證可知，對于同時(shí)含有成片野值和孤立野值的引信多普勒頻率數(shù)據(jù)，如果野值數(shù)量小于全部數(shù)據(jù)的三分之一，則利用該方法估計(jì)的脫靶量誤差在（?1.5 m～+1.5 m）范圍內(nèi)，滿足測量要求，但野值較多時(shí)的去除方法仍需加以研究。

[1]袁俊.由雷達(dá)測量數(shù)據(jù)求解脫靶量算法研究[D].南京:南京理工大學(xué),2003.

[2]林松.用導(dǎo)彈“視頻遙測”信號推算導(dǎo)彈脫靶量[J].上海航天,2002(3):35-38.

[3]魏國華,吳嗣亮,王菊,等.脫靶量測量技術(shù)綜述[J].系統(tǒng)工程與電子技術(shù),2004,26(6):768-772.

[4]胡紹林,孫國基.靶場外測數(shù)據(jù)野值點(diǎn)的統(tǒng)計(jì)診斷技術(shù)[J].宇航學(xué)報(bào),1999,20(2):68-74.

[5]張光澄.非線性最優(yōu)化計(jì)算方法[M].北京:高等教育出版社,2005:159-168.

[6]張德然.統(tǒng)計(jì)數(shù)據(jù)中異常值的檢驗(yàn)方法[J].統(tǒng)計(jì)研究,2003(5):53-55.

[7]費(fèi)業(yè)泰.誤差理論與數(shù)據(jù)處理[M].5 版.北京:機(jī)械工業(yè)出版社,2008:43-49.