李志明，范錦彪

(中北大學 儀器科學與動態(tài)測試教育部重點實驗室，山西 太原 030051)

0 引言

彈著點是評價武器著靶精度的重要指標，準確獲取彈著點位置對于武器威力評估具有重要意義[1]。彈丸落地爆炸時，會形成破壞區(qū)、塑性帶以及彈性形變區(qū)三個區(qū)域，并最終在彈性形變區(qū)形成在地下傳播的彈性波——地震波。對動能彈而言，當彈丸質(zhì)量、落地時速度不同時，其對應(yīng)動能不同，產(chǎn)生地震波的幅值和頻率也不同，動能越大，產(chǎn)生地震波的幅值越大，傳播的距離也會越遠，因此可為彈丸的落點定位提供有用信息。隨著微震定位技術(shù)發(fā)展及在許多領(lǐng)域的應(yīng)用[2-4]，彈著點定位也得到了越來越多學者的重視。張躍華等[5]改進了Geiger方法和震源掃描算法，并進行聯(lián)合計算，大幅提高了定位精度并縮短了計算時間；李鵬宇等[6]采用米字形陣列，預先通過來波方向(DOA)算法估計波速，然后將波速估計值代入到達時間差(TDOA)算法中計算初始位置，最后通過泰勒級數(shù)算法收斂定位；文獻[7-8]通過構(gòu)建方程組的方式求解彈著點坐標。

目前的彈著點定位方法多數(shù)都需要提前測量地震波波速或?qū)ζ溥M行反演，但由于地震波在地表介質(zhì)中傳播時在各個方向上的速度都不一致，即便是在同一方向上也會隨著距離的增大而衰減，因此會導致定位誤差增加。而彈著點定位的關(guān)鍵是得到彈丸落點坐標，無需求解波速，因此本文提出一種基于波速方差的目標函數(shù)，并利用單純形法搜索定位點坐標。然而利用單純形法進行搜索時，當選擇的初始單純形不同時，定位結(jié)果會存在一定的偏差，個別條件下會使目標函數(shù)收斂于局部最小值，出現(xiàn)異常的定位結(jié)果。本文針對這種情況，采用帶噪聲的密度聚類(DBSCAN)算法對隨機生成不同的初始單純形進行重復定位的結(jié)果進行判別，剔除異常值，對剩余結(jié)果求平均作為最終的彈著點坐標。

1 無速度參數(shù)測量模型

在定位系統(tǒng)中，可以直接獲得的數(shù)據(jù)只有各個測點的坐標以及測點拾取到的地震波到達時間，彈丸落地的時刻、地震波走時和地震波傳播速度等參數(shù)難以在第一時間獲得，而對于定位來講，這些恰恰是獲得彈著點位置的重要參數(shù)，傳統(tǒng)采用基于地震波TDOA算法可以有效避免對彈丸落地時刻和地震波走時的求解，但是仍然需要提前測量波速。董隴軍等[9]提出基于TDOA的新方法雖然不需要提前測量波速，但是該方法將波速作為未知參數(shù)與震源坐標一起反演求解，然而在實際應(yīng)用中，地震波的傳播速度是不均勻的，在地表介質(zhì)中傳播時，各個方向的速度都不一致，即便是在同一方向上也會隨著距離的增加而衰減，因此在求解時會導致定位誤差增加。李健等[10]提出無需測速的目標函數(shù)，雖然不需要提前測量波速以及對波速進行反演，但是，其目標函數(shù)較為復雜，當測點數(shù)量增加時，運算量也會隨之增加。因此，本文提出了如下基于波速方差的目標函數(shù)。

假設(shè)彈丸撞擊地面時，地表發(fā)生的是彈性形變，且地面無大型溝渠、無外界干擾源。在待測試區(qū)域構(gòu)建平面直角坐標系，由于坐標系原點所在位置不會對定位結(jié)果產(chǎn)生影響，所以可任意指定原點位置。假設(shè)彈著點坐標為(x0,y0)，彈丸落地時刻，即彈丸撞擊地面的時刻為t0，第i、第j個測點的坐標分別為(xi,yi)、(xj,yj)，第i、第j個測點拾取到的地震波初至時間，即彈丸落地時刻與地震波走時之和分別為ti、tj，地震波傳播速度為v，則有

(1)

(2)

兩式相減，可得

ti-tj=(Li-Lj)/v，

(3)

式中：

(4)

(5)

易得

(6)

(7)

因此，速度的方差為

(8)

(8)式即為無速度參數(shù)的目標函數(shù)。

測點的位置坐標是一定的，在一次彈丸落地事件中，測點拾取到地震波的初至時間也是一定的。如果(x0,y0)不是真實的彈著點坐標，那么計算得到的各個速度值之間偏差會較大，因此計算得到的速度方差也會較大。只有真實的彈著點坐標才會使速度方差最小，所以，尋找到可以使速度方差最小的點即為真實的彈著點。

2 模型仿真驗證

2.1 模型參數(shù)

目標函數(shù)(8)式中，未知參數(shù)僅為彈著點坐標，無需求解波速，減少了誤差來源，具有更高的可靠性。用單純形法對其尋優(yōu)即可實現(xiàn)定位，為驗證算法是否有效，構(gòu)建如圖1所示的仿真定位模型。假定監(jiān)測范圍為100 m×100 m的區(qū)域，共布置8個測點，編號分別為1～8. 在實際應(yīng)用中，每個測點由一個檢波器和一個波形記錄儀組成。A、B兩點為假定的彈著點。測點及假定彈著點的坐標如表1所示。

表1 測點及假定彈著點坐標

圖1 仿真定位模型

假設(shè)彈丸落地時刻t0=0 s，地震波的傳播速度v0=400 m/s，但由于地震波在地表介質(zhì)中傳播的速度是不均勻的，因此為了盡量模擬實際情況，在彈著點至各個測點的方向上引入±3%的速度誤差，引入方式如(9)式所示：

v=400±400×3%×rand，

(9)

式中：rand為0～1之間的隨機數(shù)。

v代入(1)式即可得到引入速度誤差后測點拾取到的地震波初至時間，8個測點拾取到彈著點A和B的地震波初至時間如表2所示。

表2 各測點拾取彈著點A和B地震波初至時間

2.2 模型驗證

表1所示測點坐標以及表2所示地震波初至時間代入目標函數(shù)(8)式，并通過單純形法搜索其最小值。本文的彈著點定位是平面定位問題，因此單純形法的初始迭代單純形為三角形，分別以三角形的各頂點坐標作為彈著點計算3種情況下的目標函數(shù)值，根據(jù)單純形法的迭代規(guī)則，通過拉伸、收縮、對稱等變換，使其不斷向著目標函數(shù)值最小的方向移動，直到滿足迭代終止條件即可得到準確的彈著點坐標，本文中單純形法的迭代終止條件設(shè)為三角形各邊長均值不大于5 cm.當選取的初始迭代三角形不同時，單純形法的定位結(jié)果相互之間也會有一定的偏差，甚至會使目標函數(shù)收斂于局部最小值，導致出現(xiàn)定位異常的情況。鑒于上述情況，在待監(jiān)測區(qū)域內(nèi)隨機生成不同的初始迭代三角形進行重復定位，觀察其定位結(jié)果的分布情況。模型中彈著點A和B的200次定位結(jié)果分布情況如圖2所示。

圖2 彈著點A和B的200次定位結(jié)果分布圖

由圖2可以看出，200次的定位結(jié)果在大多數(shù)情況下都比較集中，相互之間的偏差較小，只有少數(shù)定位出現(xiàn)異常的情況。

為了最終獲得彈著點坐標，需首先剔除定位異常的結(jié)果。DBSCAN算法是一維或多維特征空間中的非參數(shù)、基于密度的聚類算法。與劃分聚類算法不同，DBSCAN算法不需要預先聲明聚類的數(shù)量，而是通過評估樣本的緊密程度來劃分對應(yīng)的類別，理論上可以找出任何形狀的聚類。

DBSCAN算法將簇定義為密度相連點的最大集合，可以將數(shù)據(jù)集中不包含在任何簇中的對象視為噪聲，通過這樣的方法，能夠?qū)⒍ㄎ唤Y(jié)果中的異常值剔除。本文中鄰域半徑ε設(shè)置為0.05，密度閾值Minpts設(shè)置為4，DBSCAN算法運行結(jié)果如圖3所示。

圖3 彈著點 A和B的DBSCAN算法運行結(jié)果

剔除圖3中的異常值后，分別對其余有效值的橫、縱坐標分量求均值作為最終的定位結(jié)果。彈著點A的最終定位結(jié)果為(18.498 6 m, 61.191 6 m)，與A點的距離誤差為1.916 8 m；彈著點B的最終定位結(jié)果為(38.222 3 m,73.410 3 m)，與B點的距離誤差為3.845 9 m.測點至彈著點距離計算誤差絕對值平均值[11]為

(10)

式中：n為測點數(shù)量；li為計算彈著點位置到各測點的距離；di為實際彈著點位置到各測點的距離。

根據(jù)(10)式可得，彈著點A的定位誤差為1.97%，彈著點B的定位誤差為3.97%，定位誤差較小，所以該彈著點定位方法可行。

3 驗證試驗

本文通過靜爆試驗對該方法進行驗證，靜爆試驗與彈丸侵徹產(chǎn)生的地震信號一樣，都是作用時間短、信號能量高度集中的信號，區(qū)別在于靜爆是由于超壓作用于地面所產(chǎn)生的地震波，而彈丸侵徹是通過撞擊產(chǎn)生的地震波。爆炸或者彈丸侵徹地面時都會形成具有一定能量的沖擊波，沖擊波的壓力遠遠大于土壤的抗壓強度，因此靠近震源區(qū)域的土壤結(jié)構(gòu)被破壞從而形成破壞區(qū)，之后隨著能量的消耗，沖擊波會衰減為應(yīng)力波，此時應(yīng)力波依然會對土壤造成塑性形變，最后隨著能量再次衰減，塑性波變?yōu)閺椥圆?，即為地震波?/p>

2019年11月在華北某靶場進行了地面靜爆試驗，藥柱當量在50～60 kg梯恩梯(TNT)之間，共布置6個測點，選用西安石油儀器勘探廠生產(chǎn)的動圈式速度型DFJ3A檢波器，其主要技術(shù)指標[12]如表3所示。

表3 DFJ3A檢波器主要技術(shù)指標

測點和爆點分布情況如圖4所示。圖4中，爆點所在位置坐標為(0 m, 0 m)，為保證各個測點保持時鐘同步，采用了斷線觸發(fā)的方式，斷線的時刻記為0 s時刻，地面靜爆試驗中所測得的典型地震波信號如圖5所示。

圖4 爆點和測點分布圖

圖5 實測典型地震波信號

測點坐標以及各個測點拾取到的地震波初至時間如表4所示。

表4 測點坐標及地震波初至時間

同樣，表4中數(shù)據(jù)代入(8)式并通過單純形法進行200次定位，定位結(jié)果的分布情況如圖6所示。

圖6 實測200次定位結(jié)果分布圖

由圖6可以看出，與仿真定位模型的結(jié)果相似，也出現(xiàn)了少量定位異常的情況，DBSCAN算法的聚類結(jié)果如圖7所示。剔除圖7中異常結(jié)果后，分別對其余結(jié)果的橫、縱坐標分量求均值，即可得到爆點的定位結(jié)果為(-0.166 1 m, 0.239 0 m)，與真實爆點位置的距離誤差為0.291 1 m，測點至爆點距離誤差絕對值的平均值為0.38%.

圖7 實測DBSCAN算法運行結(jié)果

4 結(jié)論

本文提出一種基于波速方差目標函數(shù)的彈著點定位方法，通過仿真定位模型和地面靜爆試驗對該方法進行驗證。得出以下主要結(jié)論：

1)用單純形法對基于波速方差的目標函數(shù)尋優(yōu)并結(jié)合DBSCAN算法可以有效實現(xiàn)彈著點定位。

2)該彈著點定位方法無需提前測量地震波波速或?qū)Σㄋ龠M行反演，可以減小地震波波速不均勻性給定位結(jié)果帶來的誤差，提高定位精度。

然而，靜爆與彈丸侵徹產(chǎn)生的地震波存在較大差異，所以不宜采用靜爆試驗作驗證，但在本文的目標函數(shù)中已經(jīng)排除了波速的影響，且沒有用到幅度和頻率信息。所以仍有一定的借鑒意義，后續(xù)將進行大當量戰(zhàn)斗部試驗進一步驗證。