一種聲學測試彈著點坐標的穩健估計方法

王玉龍 盧小汐 張亞輝 李宏凱 蘇磊

摘? 要： 針對火炮彈著點聲檢靶系統中彈丸著靶坐標精度不穩定的問題，建立抗差最小二乘法。該方法通過等價權函數將抗差估計與最小二乘法結合起來，利用殘差對異常值進行降權處理，從而消除異常值對精度的不良影響。最后，對火炮彈著點聲檢靶系統進行仿真分析。分析和試驗表明，抗差最小二乘法對異常值是穩健的，且精度優于一般最小二乘法，當靶幅為10 m×10 m且風速測試誤差為3 m/s時，火炮彈著點聲檢靶系統針對火炮彈丸彈著點坐標測試精度優于5 cm。

關鍵詞： 彈著點坐標; 穩健估計; 抗差最小二乘法; 抗差估計; 殘差變化; 聲定位

中圖分類號： TN911.1?34; TJ306? ? ? ? ? ? ? ? ? ?文獻標識碼： A? ? ? ? ? ? ? ? ? ? 文章編號： 1004?373X（2020）09?0025?04

A robust estimation method of acoustic test for impact point coordinate

WANG Yulong， LU Xiaoxi， ZHANG Yahui， LI Hongkai， SU Lei

（Huayin Ordnance Test Center， Huayin 714200， China）

Abstract： A robustified least square method is established for the instability coordinate precision of the sound detection target system for the projectile impact point. In the method， the robust estimation is combined with the least square method by the equivalent weight function， and the residual is used to reduce the weight of the outliers to eliminate the adverse effects of the outliers on the accuracy. The simulation analysis of sound detection target system for the projectile impact point was carried out. The analysis and test results show that the robustified least square method is robust to outliers and its accuracy is higher than that obtained with the general least square method. When the target area is 10 m×10 m and the measurement error of wind speed is 3 m/s， the tested coordinate accuracy of the impact point obtained with the sound detection target system for projectile impact point is within 5 cm.

Keywords： impact point coordinate; robust estimation; robustified least square method; robust estimation; residual variation; acoustic positioning

0? 引? 言

聲學立靶彈著點測試法屬于基于到達時間差的無源定位，一般可得到一組非線性雙曲方程組。一般采用最小二乘法（Least Square，LS）進行求解[1?3]。在彈著點聲學靶系統試驗測試過程中發現，有時傳感器測量的時差存在粗差，即異常值或野值。由于測試系統中每個傳感器所處位置不同，根據其測試波形處理得到的波達時刻也不同，很難從表面上區分哪些測試數據屬于異常值或野值。在采用最小二乘法確定著靶坐標時，因為最小二乘法對異常值的敏感性，使得系統測試精度變低。

本文將抗差估計[4?5]與最小二乘法結合，建立了穩健估計方法（Robust Least Square，RLS），消除了異常值對定位結果精度的影響，且計算流程與最小二乘法一致。最后，給出了仿真與試驗結果。

1? 定位原理

聲學靶系統主要采用傳聲器陣列，采集超音速彈丸產生的激波信號，再根據傳聲器坐標與激波時間差，確定彈丸通過空中虛擬平面的坐標。聲學靶及立靶坐標系示意圖如圖1所示。其中，傳聲器采用雙圓環布陣，每個圓環上均勻分布12個傳聲器，[xz]平面為水平面，[z]軸指向炮口方向，[xy]平面為豎直平面，[xy]平面即為虛擬靶平面，雙圓環陣列位于靶平面內。

計算分析和實驗表明，風對定位精度的影響不可忽略。假定風平行地面運動，且在彈丸著靶時刻靶面周圍的局部風場是均勻的。因為風在[z]方向的分量[vz]影響較小[6?8]，所以主要考慮風在[x]方向的分量[vx]。利用矢量疊加原理可得到定位方程組如下：

[（x1-x-vx（t0+t1））2+（y1-y）2=v0（t0+t1）（x2-x-vx（t0+t2））2+（y2-y）2=v0（t0+t2）? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ??（xn-x-vx（t0+tn））2+（yn-y）2=v0（t0+tn）] （1）

式中：[ti，i=1，2，…，n]是激波到達各個傳感器的相對時間（即時延或時差）;[t0]是激波視速度在靶面上傳播至首先觸發傳感器的時間;[v0]是視速度。

方程組（1）中的未知參量為[x，y，v0，t0]。當[n≥4]時方程組（1）有解;當[n=4]時，方程組（1）是非線性方程組;[n>4]時，方程組（1）是超定非線性方程組（矛盾方程組）。求解方程組（1）即得[x，y，v0，t0]。

2? 彈著點穩健估計方法

2.1? 迭代公式

最小二乘法理論完善，算法成熟，得到了廣泛應用。但是，當觀測數據含有粗差（或異常值）時，經典最小二乘法的結果可能因異常值的存在而不可靠。在聲靶試驗過程中，由于傳感器偶然抖動、風吹等不確定因素的影響，傳感器測量的時差有時會出現異常值。此時，應用最小二乘法計算時，導致著靶測試精度低。穩健估計（也叫抗差估計）是指在存在異常值的情況下，通過選擇適當的估計方法，使估計結果盡可能少地受到異常值干擾。穩健估計的原則是充分利用有效信息，限制利用可用信息，排除有害信息。

抗差最小二乘法通過等價權將抗差估計與最小二乘法結合起來，其實質是等價權函數的設計。當實測時差不含異常值時，抗差最小二乘法的結果與最小二乘法結果一致;當實測時差含有異常值時，抗差最小二乘法能夠消除異常值的影響。

下面詳細推導抗差最小二乘法的迭代公式，記：

[fi=（x-xi+vx（t0+ti））2+（y-yi）2-v20（ti+t0）2，? ? ? ? ? ? ? ? ? ? ? ? ? ? ?i=1，2，…，n]? ?（2）

設[f： D?R4→Rn]，[f=（f1，f2，…，fn）T]。求解非線性方程組（1）等價于求解如下非線性最小二乘問題。

[minx∈D12fT（x）f（x）]? ?（3）

利用Taylor級數展開將非線性問題轉化為線性問題，在[x0=（x0，y0，v00，t00）]展開得：

[f（x0）+DfT（x0）（x-x0）=0]? ?（4）

式中：[Df（x0）=Df1DxDf1DyDf1Dv0Df1Dt0Df2DxDf2DyDf2Dv0Df2Dt0????DfnDxDfnDyDfnDv0DfnDt0x0]是Jacobi矩陣;[DfiDx=2（x-xi+vx（t0+ti））];[DfiDy=2（y-yi）];[DfiDv0=][-2v0（ti+t0）2];[DfiDt0=-2v20（ti+t0）];[i=1，2，…，n]。

記[b=-f（x0）]，[A=DfT（x0）]，[X=（x-x0）]，則：

[AX=b] （5）

進而，有：[x1=x0+（ATP0A）-1ATP0b]。

一般地，抗差最小二乘迭代公式如下：

[xk+1=xk+（ATkPkAk）-1ATkPkbk，k=0，1，2，…] （6）

式中[Pk]稱為等價權矩陣或等價權函數。特別地，當取[Pk=I]時，抗差最小二乘法退化為最小二乘法。

就數學表達式而言，抗差最小二乘法與加權最小二乘法基本相同，其區別是權函數的內容不同。加權最小二乘法的權是先驗的，而抗差最小二乘法的權是殘差的函數。

2.2? 等價權矩陣的確定

權函數的作用在于：將所有的觀測數據劃分為正常觀測值、可利用觀測值和粗差觀測值（異常值）三個部分。對于正常觀測值，使其保持原有的權不變，對余差較大的可利用觀測值進行降權處理，而那些粗差觀測值，則使其權為零。經過這樣的處理，就可以根據各個觀測值及“后驗權”，求得最接近于正常情況下的結果。本文中權函數選為IGGⅢ權函數（IGG是中國科學院測量與地球物理研究所的英文縮寫），其表達式如下：

[Pi（vi）=1，? ? vi≤k0k0vik1-vik1-k02，? ? k0k1] （7）

式中：[vi=δiσ]，[δi=fi]是殘差，[σ]是標準化殘差，[σ=i=1nδi（n-1）]，[k0∈（1，1.5）]，[k1∈（2.5，8.0）]，本文取[k0=1.5]，[k1=3]。

IGGⅢ權函數充分考慮了正常階段、可疑階段和淘汰階段三個部分，公式簡單，易于實現，如圖2所示。

2.3? 算法流程圖

抗差最小二乘法的計算流程如圖3所示。

2.4? 等價權和殘差變化的算例

以某型榴彈試驗為例進行分析，圖4給出了迭代過程中不同傳感器對應的殘差和等價權，其中橫坐標代表傳感器編號，縱坐標分別代表殘差、標準化殘差和等價權。

由圖4可知，在迭代的過程中，有多個傳聲器的等價權發生變化，隨著迭代的不斷進行，最后6號傳聲器的等價權減小，并降為0，即最終判斷6號傳感器存在異常。圖5給出了對應的彈著點誤差隨迭代次數的變化曲線。

由圖5可知，隨著迭代次數的增加，水平誤差基本不變，而高低誤差由42.3 cm降低到1.0 cm，即抗差最小二乘法消除了異常值的影響。

3? 仿真分析

影響聲學靶測量誤差的因素主要有：傳感器坐標測試誤差、傳感器時延誤差、溫度測試誤差和風速測試誤差。傳感器坐標測試誤差和時延誤差由設備決定，溫度變化比較穩定。聲靶計算修正時用到的環境風速是激波在靶面附近傳播幾毫秒間的瞬態風速，測試誤差較大。下面以靶幅10 m×10 m，傳感器坐標測試誤差為1 mm，傳感器時延誤差為5 μs，溫度測試誤差為5 ℃，風速測試誤差分別為1 m/s，2 m/s，3 m/s和4 m/s時，給出誤差分布圖，如圖6，圖7所示。

由圖6和圖7知，隨著風速測試誤差的增加，水平誤差和高低誤差也隨之增加。表1給出了在不同風速測試誤差時的均方根誤差。

4? 試驗數據結果與分析

4.1? 鞭炮模擬試驗

采用鞭炮爆炸激發傳聲器，事先測量出鞭炮炸點和傳聲器的相對關系，作為真值，檢測聲學靶的精度。分別采用LS法、RLS法和交會法進行計算，對比不同算法的精度。其中，交會法采用文獻[8]中的交會計算方法。試驗結果如圖8所示。

由圖8可知，交會法的水平均方根誤差為5.0 cm，高低均方根誤差為6.0 cm;LS法的水平均方根誤差為2.8 cm，高低均方根誤差為8.1 cm;RLS法的水平均方根誤差為2.8 cm，高低均方根誤差為2.9 cm。從而可知，LS法容易受到異常值的影響，而RLS法能夠保證測量的穩健性，并且RLS法的精度優于交會法。

4.2? 某型榴彈試驗

某型榴彈試驗同樣分別采用LS法、RLS法和交會法進行計算，對比不同算法的精度。結果如圖9所示。

由圖9知，交會法的水平均方根誤差為5.5 cm，高低均方根誤差為7.6 cm;LS法的水平均方根誤差為3.2 cm，高低均方根誤差為15.4 cm;RLS法的水平均方根誤差為3.3 cm，高低均方根誤差為4.5 cm。與鞭炮模擬試驗結果類似，RLS法精度最高。

5? 結? 論

本文針對火炮彈著點聲學靶系統定位計算方法進行研究，建立了基于抗差最小二乘法的穩健估計方法，提高了彈丸著靶坐標的測試精度。對火炮彈著點聲學靶系統進行仿真分析，風速測試誤差越大，著靶坐標測試精度越低。當風速測試誤差為3 m/s時，著靶坐標測試精度為5 cm。

參考文獻

[1] 王磊，陳昭男.基于被動聲學的超聲速武器脫靶量解算方法[J].應用聲學，2018，37（3）：385?390.

[2] 苗晟，周維，唐浩，等.一種聲源定位系統設計[J].計算機科學，2013，40（11A）：398?400.

[3] 劉哲，陳日林，滕鵬曉，等.基于平面傳聲器陣列的聲源定位系統[J].聲學技術，2011，30（2）：123?128.

[4] 張洋，張志剛，錢棟，等.基于IGGⅢ模型的高精度抗差坐標轉換算法研究[J].全球定位系統，2017，42（5）：29?32.

[5] 李蟬，張士峰，張力軍.應用抗差無跡卡爾曼濾波的再入彈道處理技術[J].國防科技大學學報，2017，39（1）：1?5.

[6] MSI. Strafe scoring target type 590 air?to?ground training [EB/OL]. [2014?07?20]. http：//www.msinstruments.com.

[7] 蔣東東，狄長安.八點線式聲學立靶彈著點檢測系統研究[J].電子測量技術，2010，33（4）：19?21.

[8] 馮斌，石秀華.雙三角陣聲靶測試系統研究[J].應用聲學，2012，31（2）：140?144.

[9] 董杰，高樓.彈丸任意角度入射彈著點聲學檢測模型[J].計算機仿真，2015，32（1）：10?14.