冷兆龍，劉高峰，王子齊

(1.海軍工程大學， 武漢 430033；2.中國人民解放軍91213部隊，山東 青島 266000)

1 引言

潛標也被稱為水下浮標。通過浮體和錨系設(shè)備，潛標可以停留在預設(shè)的作業(yè)位置和作業(yè)深度，持續(xù)工作數(shù)月甚至更長的時間[1-2]。以往，潛標主要應(yīng)用于海洋環(huán)境監(jiān)測等民用領(lǐng)域。為創(chuàng)新反潛作戰(zhàn)樣式，文獻[3]提出了一種將潛標陣應(yīng)用于反潛作戰(zhàn)的構(gòu)想，即在重要的海峽、航道或其他海域，前置部署由若干被動聲探測潛標組成的潛標陣，用于反潛目標的偵察和指示。在潛標陣中，各潛標使用被動聲探測方式，僅能獲得目標方位，再通過方位對目標運動要素進行估計，這類問題可歸屬于純方位目標運動分析(BOTMA)問題。文獻[4]提出了一種戰(zhàn)場威脅約束下的觀測站運動軌跡優(yōu)化算法；文獻[5-7]分別從計算復雜度、計算精度、消除估計有偏性的角度對BOTMA算法進行了改進；文獻[8]針對未知分布的噪音，構(gòu)建了一種改進的自適應(yīng)濾波器；文獻[9]改進了伯努利濾波器，用于多目標情況下的純方位目標跟蹤；文獻[10-12]針對濾波器初值選擇問題進行了分析。已有的BOTMA研究，從研究內(nèi)涵來看，大多聚焦于觀測站運動軌跡、濾波器改良、濾波器初值選擇等算法優(yōu)化方面；從應(yīng)用場景來看，主要有陸基監(jiān)聽站、艦艇、潛艇、反潛機等裝備搭載平臺。

由于潛標陣應(yīng)用的特殊性，已有的BOTMA算法大多無法直接運用。本文建立了一種“單潛標數(shù)據(jù)獨立處理(single submarine buoy data independent processing，SSBDIP)和多潛標數(shù)據(jù)集中融合(multi submarine buoy data centralized fusion，MSBDCF)”的潛標陣目標運動分析描述模型，采用最小二乘法，開展非機動/機動目標運動要素解算研究。

2 潛標陣目標運動分析問題描述

2.1 傳統(tǒng)多站BOTMA問題算法特性

根據(jù)觀測站數(shù)量，BOTMA問題可分為單站BOTMA和多站BOTMA。相關(guān)研究已表明，單站BOTMA問題在觀測站不進行有效機動情況下，目標運動狀態(tài)不具備完全可觀測性。多站BOTMA問題中，觀測數(shù)據(jù)維度增加了，因此對觀測站是否機動沒有要求。傳統(tǒng)多站BOTMA問題的解決思路如圖1所示。

圖1 傳統(tǒng)多站BOTMA問題的解決思路框圖

先獲取目標位置的粗略估計，然后利用粗略估計的位置優(yōu)化目標運動模型參數(shù)，同時利用運動模型修正位置估計，通過循環(huán)反饋的方式，得到較高精度的估計結(jié)果。

傳統(tǒng)多站BOTMA問題中，對位置的粗略估計通常采用方位線交匯定位的方式，如圖2所示。圖2(a)中，觀測站數(shù)量為2時，取觀測所獲得方位線的交點作為位置的粗略估計；圖2(b)中，觀測站數(shù)量大于2時，由于觀測誤差影響，方位線往往不能交于一點，而是形成了一個多邊形，這時可以按一定準則從多邊形中選擇某個點(如選擇多邊形的重心)作為位置的粗略估計。顯然，方位線交匯定位需要不同觀測站對目標同一時刻的觀測數(shù)據(jù)，或者說需要各觀測站的實時觀測數(shù)據(jù)，這也是絕大多數(shù)多站BOTMA算法對觀測站的要求。

圖2 傳統(tǒng)多站BOTMA問題中位置的粗略估計示意圖

相對于傳統(tǒng)以平臺為中心的探潛方式，潛標陣在反潛目標探測應(yīng)用中具有一定的特殊性，主要表現(xiàn)在以下方面：

1) 潛標陣遠離本土部署在重要海域，與指揮所之間難以建立有線通信。為將潛標陣探測的目標信息回傳到指揮所，潛標需要上浮，通過衛(wèi)星鏈路與指揮所進行數(shù)據(jù)通信。

2) 各潛標觀測的目標角度數(shù)據(jù)多，通過衛(wèi)星鏈路回傳時，傳輸數(shù)據(jù)量需進行一定的限制，否則會占用大量通信資源，進而影響通信的可靠性和穩(wěn)定性。

3) 指揮所得到的目標信息數(shù)據(jù)是潛標回傳之前的數(shù)據(jù)，存在一定時延。

4) 為保持潛標的隱蔽性，各潛標相互獨立且彼此之間不通信。

2.2 潛標陣目標運動分析描述模型

潛標陣目標運動分析問題屬于多站BOTMA問題,但不具備直接運用傳統(tǒng)多站BOTMA問題相關(guān)算法的條件?；蛘哒f，基于傳統(tǒng)多站BOTMA問題算法，目標對于潛標陣“部分可觀測”而非 “完全可觀測”。已有研究表明：對于靜止單站被動觀測，可以解算出目標運動的部分要素，該結(jié)論可以推廣應(yīng)用于潛標陣目標運動分析之中?；舅悸肥牵菏紫?，各潛標獨立解算，得到目標運動的部分要素；其次，通過多潛標數(shù)據(jù)集中融合，解算出目標運動的全部要素。提出的潛標陣目標運動分析描述模型如圖3所示。

圖3 潛標陣目標運動分析描述模型框圖

當目標經(jīng)過潛標陣部署海域時，各潛標被動感知目標方位，并對方位數(shù)據(jù)自行處理，得到目標運動的部分要素。當潛標判斷目標已離開該潛標探測范圍后，自行上浮，將處理結(jié)果和潛標位置、潛標觀測目標的起止時刻、運動分段起始時刻的角度觀測值等數(shù)據(jù)，通過衛(wèi)星鏈路回傳至指揮所。指揮所的情報中心獲得多個(2個以上)潛標數(shù)據(jù)，進行時空對準，然后解算出目標運動的全部要素。由于潛標在目標離開探測范圍后才上浮回傳數(shù)據(jù)，上浮過程可能需要耗時數(shù)分鐘左右，因此解算結(jié)果是非實時的，存在時延。針對這一問題，可以通過預估目標實時位置(利用目標航向、航速)的方式來解決。盡管這樣處理存在誤差，考慮到潛艇的運動速度和機動能力都是不很強，預估的目標運動要素可認為基本足以支撐反潛作戰(zhàn)指揮決策，也可以為反潛武器運用提供一定參考。

3 單潛標數(shù)據(jù)獨立處理算法設(shè)計

單潛標數(shù)據(jù)獨立處理(SSBDIP)算法目的在于利用最小二乘法，估計得到目標運動的部分要素。

3.1 非機動目標的數(shù)據(jù)處理算法

假設(shè)目標做勻速直線運動，如圖4所示。

圖4 單個潛標對勻速直線運動目標的觀測示意圖

任選一點作為坐標系原點，以正東為軸正方向，正北為x軸正方向，建立平面直角坐標系。目標進行勻速直線運動，航向為Km，速度在x軸和y軸方向上的分量為vx、vy，在任意時刻t目標的位置坐標為[rx(t),ry(t)]。潛標位置坐標為(x,y)，可以對周圍半徑為Dr的圓形區(qū)域范圍進行被動探測。記潛標首次探測到目標的時刻為t0，此時目標相對潛標的方位為β(t0)，以y軸正方向為0°，順時針方向為正。潛標持續(xù)探測目標，直到目標離開潛標探測區(qū)域，把這些目標方位數(shù)據(jù)分別記為：β(t1)，β(t2)，…，β(tn)，其中，t1，t2，…,tn為目標方位數(shù)據(jù)所對應(yīng)的時刻。

根據(jù)三角函數(shù)關(guān)系，對于任意時刻ti，均有：

(1)

假設(shè)目標做勻速直線運動，有：

(2)

將式(2)代入式(1)，整理后可得：

cosβ(ti){[ry(t0)-y]tanβ(t0)-x+vx(ti-t0)}=

sin(ti)[ry(t0)-y+vy(ti-t0)],i∈[0,1,…,n]

(3)

已有研究表明：靜止單站被動觀測條件下，如果觀測站(潛標)不在目標運動軌跡及其延長線上，目標就是部分可觀測的。對于式(3)，只要獲得了至少3個不同時刻的觀測數(shù)據(jù)，便可利用最小二乘法對目標運動的2個要素[ry(t0)-y]/vy和vx/vy進行估計，估計結(jié)果為：

(4)

其中，

(5)

通過上述處理過程，得到了目標運動的2個要素：[ry(t0)-y]/vy和vx/vy。這樣，潛標在利用衛(wèi)星鏈路通信時，不需要將目標的全部觀測方位信息回傳，大大降低了衛(wèi)星鏈路通信的數(shù)據(jù)量。同時，[ry(t0)-y]/vy和vx/vy均為非時敏參數(shù)，因此對實時性沒有嚴格要求。

3.2 機動目標的識別處理

描述機動目標的運動模型有很多，在潛標陣機動目標識別過程中，選擇了一種分段勻速直線運動模型，將目標運動描述為若干段相連的勻速直線運動。該模型適用范圍廣，可以描述轉(zhuǎn)向、變速等不同類型的機動目標，而且模型參數(shù)較為簡單，便于計算。

使用分段勻速直線運動模型對機動目標觀測數(shù)據(jù)處理的過程中，關(guān)鍵問題在于識別和判斷各段勻速直線運動的起止。潛標僅能獲取目標方位角度信息，無法直接識別目標的運動狀態(tài)是否發(fā)生變化。文獻[13-14]中提出了純方位觀測條件下通過角度預測判斷目標是否發(fā)生機動的方法。本文在這一方法基礎(chǔ)上進行了改進，識別的基本過程為：

過程1) 假定目標處于勻速直線運動狀態(tài)，利用式(4) 可獲得[ry(t0)-y]/vy和vx/vy的估計結(jié)果。

(6)

過程3) 計算目標角度預測數(shù)據(jù)與目標角度觀測數(shù)據(jù)之間的差Δβ(ta+Δtc),即:

(7)

過程4) 設(shè)置檢測閾值εsingle，用于判斷識別目標是否進入了機動段。

若Δβ(ta+Δtc)>εsingle，可認為在ta+Δtc時刻，目標已經(jīng)進入了新的機動段。記錄下[t0,ta]時間里單潛標數(shù)據(jù)處理得到的目標運動部分要素，然后以ta+Δtc作為新的運動段起點t0′，返回過程1)。

過程5) 重復上述循環(huán)處理過程，直至完成了所有目標觀測數(shù)據(jù)的處理。

需要特別指出的是，過程2)中的時間窗Δtc和檢測閾值εsingle不能隨意選取，具體取值需根據(jù)潛標觀測的精度等因素確定。如圖5所示，目標在tZ時刻進入了新的運動狀態(tài)?？梢钥闯?，在tZ之前，預測值與觀測數(shù)據(jù)之間偏差值較??；在tZ之后，由于目標機動，目標真實位置與機動前的軌跡偏離越來越大，預測值與真實值之間的偏差也逐漸積累增大。如果Δtc設(shè)置過小，可能無法積累較大的偏差，進而導致無法檢測出目標運動狀態(tài)的變化；如果Δtc設(shè)置過大，會增大[ry(t0)-y]/vy和vx/vy的估計誤差，并將其誤判為目標運動狀態(tài)的變化。同理，選擇檢測閾值εsingle時，過大可能導致無法識別出目標運動狀態(tài)的變化，過小會將[ry(t0)-y]/vy和vx/vy的估計誤差誤判為目標運動狀態(tài)的變化。

從圖5中還可以看出，tZ與ta、(ta+Δtc)之間沒有確定性關(guān)系。因此，單個潛標僅能判斷目標在某一時間段里是否進入了新的運動段，無法計算目標機動時刻tZ，需要在多潛標數(shù)據(jù)集中融合過程中對tZ進行估計。

圖5 潛標對機動目標的識別處理過程曲線

通過SSBDIP算法，每個潛標便可利用觀測數(shù)據(jù)得到數(shù)組[ry(t0)-y]/vy和vx/vy，每組對應(yīng)一段勻速直線運動。

4 多潛標數(shù)據(jù)集中融合算法設(shè)計

單潛標的目標觀測數(shù)據(jù)通過SSBDIP算法處理后，得到了目標運動的2個非時敏要素：[ry(t0)-y]/vy和vx/vy，并通過衛(wèi)星鏈路將目標運動的2個要素，以及運動段起始時刻t0、運動段結(jié)束時刻tend、運動分段起始時刻的角度觀測值β(t0)、潛標位置坐標(x,y)等數(shù)據(jù)回傳至指揮所。將單個潛標傳回的一個運動分段的數(shù)據(jù)記為Xi(i用于區(qū)分不同潛標或者同一潛標觀測到的多個運動段)，有：

(8)

4.1 非機動目標的多潛標數(shù)據(jù)集中融合算法

對于非機動目標，MSBDCF算法計算過程如下。

不同潛標估計的vx/vy之間存在偏差。一般來說，潛標對目標持續(xù)觀測時間越長，對目標運動部分要素的估計精度就越高，可以通過時間加權(quán)平均的方式對vx/vy進行融合處理。設(shè)指揮所收到了N個潛標的數(shù)據(jù)，有：

(9)

式(9)中，Vx/Vy表示融合結(jié)果，可解算出目標航向Km：

Km=arctan(Vx/Vy)

(10)

2) 目標位置和速度的計算

由于各潛標相互獨立且彼此之間不通信，各個潛標回傳至指揮所的數(shù)據(jù)基于自身觀測數(shù)據(jù)和其首次發(fā)現(xiàn)目標的時間，因而需要進行時空對準。

對準方式為：在接收到多個潛標回傳數(shù)據(jù)Xi后，定義最先感知到目標信息的時刻作為統(tǒng)一的時間基點，并對各潛標首次觀測到目標的時刻校準。校準后，用T1,T2,…,TN分別表示第1個、第2個，…，第N個潛標發(fā)現(xiàn)目標時刻，其中，T1=0。用x1,x2,…,xN和y1,y2,…,yN分別代表這些潛標位置的橫坐標和縱坐標。

對準后，各潛標的[ry(t0)-y]/vy應(yīng)改寫為[ry(Ti)-yi]/vy。

利用式(2)對[ry(Ti)-yi]/vy轉(zhuǎn)化，有：

(11)

用ai表示[ry(Ti)-yi]/vy，則式(11)可看作二元一次方程：

就知識主題最晚結(jié)束年級而言，中國在4個知識主題設(shè)置時間均明顯前置．相較于六國平均水平，“未知數(shù)、變量的使用”“代數(shù)式概念”“代數(shù)式的運算”“代數(shù)式的證明”設(shè)置時間較之六國均前置兩年．從一定程度上可以說明中國代數(shù)思維課程內(nèi)容設(shè)置時間相對較短．

ry(T1)+vy(Ti-T1-ai)=yi

(12)

式(12)中僅有ry(T1)和vy2個未知數(shù)。當指揮所接收到2組以上的數(shù)據(jù)Xi，便可以利用最小二乘法進行求解。

另外，可能出現(xiàn)一種特殊情況，觀測到目標的所有潛標正好處于一行，即y1=y2=…=yN。此時式(12)退化為：

ry(T1)+vy(Ti-T1-ai)=y1

(13)

這種情況下，即使獲取了多組Xi，式(13)也無法求解。此時，需要利用式(1)對式(12)進行變換，整理后得：

(14)

式(14)中，僅有rx(T1)和vy等2個未知數(shù)。而且，對于不同位置的若干潛標，如果y1=y2=…=yN，則x1,x2,…=,xN必然互不相等。指揮所獲得2個或以上Xi后，便可進行求解。

獲得ry(T1)或rx(T1)和vy后，利用式(1)以及Vx/Vy，可以計算得到rx(T1)或ry(T1)和vx，再通過式(2)可以計算任意時刻目標位置。

4.2 機動目標的多潛標數(shù)據(jù)集中融合算法

對于機動目標，潛標對目標觀測數(shù)據(jù)通過SSBDIP算法，可得到若干組數(shù)據(jù)Xi，每組Xi對應(yīng)一段勻速直線運動。多潛標數(shù)據(jù)集中融合過程中，需要先剔除掉可能存在的嚴重偏差，再將不同潛標對同一運動段的Xi關(guān)聯(lián)起來，最后利用4.1節(jié)公式解算各勻速直線運動段的全部運動要素。

嚴重偏差是一種特殊情況，如圖6所示。單潛標對機動目標識別處理時，假定目標處于勻速直線運動，利用觀測數(shù)據(jù)處理獲得目標運動的2個要素：[ry(t0)-y]/vy和vx/vy，待運動要素估計收斂到一定精度范圍后，進行目標角度預測。如果目標剛進入潛標探測范圍后不久，在[ry(t0)-y]/vy和vx/vy收斂前就進行了機動(圖6中實線)，那么潛標不僅無法識別此變化，反而會將其目標運動狀態(tài)變化前后的軌跡擬合成一段勻速直線運動，導致了嚴重偏差估計的出現(xiàn)(圖6中虛線)。

圖6 嚴重偏差的示意圖

顯然，此時單潛標數(shù)據(jù)處理得到的Xi與目標運動是不匹配的，該組數(shù)據(jù)會影響整體的解算效果，應(yīng)進行剔除。

考慮到潛標陣的主要探測對象為潛艇，通常不會進行頻繁的機動。因此，可以采用以下辦法剔除嚴重偏差：

指揮所將收到的各潛標傳回的目標各個運動段的數(shù)據(jù)Xi，按運動段起始時間Ti進行排序。依次選取時間上相鄰的2組數(shù)據(jù)，按照4.1節(jié)中的計算方法進行融合解算。將這些解算結(jié)果中的vy分別記為：vy(1,2),vy(2,3),…,vy(i-1,i),，括號中的數(shù)字下標表示參與解算數(shù)據(jù)Xi的序號。

過程1) 取vy(1,2),vy(2,3),…,vy(i-1,i)中相鄰的2個數(shù)值，計算其相對平均偏差δk為：

(15)

過程2) 設(shè)置判斷閾值εmulti用于評價的大小。如果δk≤εmulti，說明vy(k-1,k)和vy(k,k+1)大小接近，參與解算的3組數(shù)據(jù)Xk-1、Xk和Xk+1屬于同一段勻速直線運動。返回到過程1)，計算δk+1。

過程3) 如果δk>εmulti，說明vy(k-1,k)和vy(k,k+1)大小相差較大，參與解算的3組數(shù)據(jù)Xk-1、Xk和Xk+1不屬于同一段勻速直線運動?？紤]到過程2)中循環(huán)處理的條件，不妨假設(shè)δk-1≤εmulti，那么Xk-1和Xk屬于一段勻速直線運動分段。因此，Xk+1有2種可能情況：①Xk+1屬于一段新的勻速直線運動；②Xk+1是嚴重偏差估計。

過程4) 計算δk+2并利用εmulti評價其大小。

如果δk+2≤εmulti，說明Xk+1、Xk+2和Xk+3屬于一段新的勻速直線運動；

如果δk+2>εmulti，說明Xk+1、Xk+2和Xk+3不屬于同一段勻速直線運動。對于Xk+1來說，沒有與其屬于同一段勻速直線的其他數(shù)據(jù)?？紤]到潛艇目標不會進行頻繁的機動，可認為Xk+1就是需要剔除的嚴重偏差估計。

過程5) 完成對Xk+1的判斷后，回到過程1)，計算δk+3。

通過上述過程，完成了剔除嚴重偏差估計和數(shù)據(jù)關(guān)聯(lián)工作。關(guān)聯(lián)后的數(shù)據(jù)，按照4.1節(jié)的方法便可解算各段勻速直線運動的要素。

另外，單潛標目標觀測數(shù)據(jù)通過SSBDIP算法，無法給出目標機動時刻tZ。解算出各個勻速直線運動分段的要素后，取相鄰2段運動的軌跡延長線交點作為前一運動段的終點(同時也是后一運動段的起點)，計算目標到達該位置的時刻，作為tZ的估計結(jié)果。

5 仿真實驗與分析

5.1 仿真實驗想定

假設(shè)，某特定海域部署了5行×5列，共計25枚潛標組成的潛標陣，相鄰兩行、相鄰兩列間的距離均為10 km，如圖7所示。潛標的探測能力，參考文獻[15]，取探測范圍半徑Dr=8 km，觀測周期T為1 s，觀測誤差假定為高斯分布，均方差取1.5°，均值為0。生成4條模擬潛艇航跡，航跡參數(shù)見表1。SSBDIP算法中，時間窗Δtc取200 s，檢測閾值εsingle取2.5°；MSBDCF算法中，判斷閾值εmulti取0.1。

表1 模擬潛艇航跡參數(shù)

5.2 仿真實驗分析

利用Matlab軟件，對每條航跡各進行50次蒙特卡洛仿真實驗。圖8～圖11給出了算法對4條航跡的仿真航跡和位置估計誤差變化曲線。

圖8 航跡1的實驗結(jié)果

圖9 航跡2的實驗結(jié)果

圖10 航跡3的實驗結(jié)果

圖11 航跡4的實驗結(jié)果

表2給出了算法對目標機動的識別結(jié)果以及誤差，表2中數(shù)據(jù)四舍五入省略了小數(shù)部分，潛艇目標實際機動位置為(19 995 m,11 311 m)。

表2 目標機動識別結(jié)果以及誤差

在大幅度轉(zhuǎn)向運動目標(航跡2)情況下，算法的位置估計誤差較大；在小幅度轉(zhuǎn)向運動目標(航跡3)情況下，算法的位置估計誤差較小；變速不變向運動目標(航跡4)情況下，算法的位置估計誤差介于前面2種情況之間。具體分析如下：① 從圖9可見，對大幅度轉(zhuǎn)向運動目標，算法在目標機動時刻(2 000 s)附近的位置估計誤差很大，是因為算法對目標機動時刻tZ的估計不夠準確引起，之后呈逐步下降趨勢。② 從圖10可見，目標機動后，算法的位置估計誤差稍遜于勻速直線運動狀態(tài)下的目標，但位置估計誤差的變化幅度沒有超過400 m。③ 從圖11可見，對于變速不變向運動目標，位置估計誤差在機動時刻附近很大，但由于目標航向未發(fā)生變化，因此位置估計誤差在一段時間后迅速降低。

6 結(jié)論

1) 對勻速直線運動目標的位置估計。在設(shè)定的仿真條件下，在潛標陣觀測范圍內(nèi)，算法的位置估計誤差在200 m以內(nèi)。

2) 對機動目標的識別判斷。3種不同運動方式的機動目標，其算法均能夠識別出目標機動，但機動時刻的估計誤差較大。實際上，減小時間窗Δtc可以降低tZ的估計誤差，但對于小幅度轉(zhuǎn)向機動目標等情況，可能無法識別出目標機動。

3) 通過仿真實驗可見：建立的SSBDIP-MSBDCF的綜合算法，實現(xiàn)了潛標陣純方位目標運動分析從“部分可觀測”到“完全可觀測”的突破，綜合算法的穩(wěn)定性和精度良好。