周天，沈嘉俊，陳寶偉，周曹韻，李海森，王楠

(1. 哈爾濱工程大學(xué) 水聲技術(shù)重點(diǎn)實(shí)驗(yàn)室，黑龍江 哈爾濱 150001；2. 海洋信息獲取與安全工業(yè)和信息化部重點(diǎn)實(shí)驗(yàn)室(哈爾濱工程大學(xué))，黑龍江 哈爾濱 150001；3. 哈爾濱工程大學(xué) 水聲工程學(xué)院，黑龍江 哈爾濱 150001；4. 中國(guó)船舶重工集團(tuán)公司 第七一九研究所，湖北 武漢 430064)

近年來，三維前視聲吶作為一種重要的水下探測(cè)系統(tǒng)，在海底地形測(cè)繪、石油管道檢測(cè)、沉船打撈、水下航行器避障、水下機(jī)器視覺等領(lǐng)域具有廣泛的應(yīng)用[1-3]。目前，國(guó)外多個(gè)公司已研制出實(shí)用化的三維前視聲吶設(shè)備。其中典型代表如TriTech公司Eclipse[4]、CodaOctopus公司的Echoscope4G[5]，Reson公司SeaBat 7130[6]等。Eclipse聲吶采用“T”型陣列結(jié)構(gòu)，以水平向波束形成、垂直向掃描的方式完成三維探測(cè)，然而垂直掃描方式對(duì)于三維探測(cè)來說仍然效率偏低；Echoscope4G聲吶則采用了二維接收面陣設(shè)計(jì)，采用二維方位搜索技術(shù)完成三維探測(cè)，其陣元總數(shù)為M2=482，系統(tǒng)復(fù)雜度可近似描述為O(M2)。顯然地，該聲吶雖然具備目標(biāo)三維探測(cè)能力，但系統(tǒng)實(shí)現(xiàn)復(fù)雜度非常高。相比之下，SeaBat 7130聲吶采用二維稀疏矩形面陣設(shè)計(jì)，其陣元總數(shù)為MN,系統(tǒng)復(fù)雜度可近似描述為O(MN)。因而，相比Echoscope聲吶，SeaBat 7130聲吶系統(tǒng)具有更低的系統(tǒng)復(fù)雜度。

現(xiàn)階段公開報(bào)道國(guó)內(nèi)三維前視聲吶研制工作中，文獻(xiàn)[7]采用“T”型陣列結(jié)構(gòu)，文獻(xiàn)[8]采用覆蓋不同探測(cè)區(qū)域的多條水平線陣結(jié)構(gòu)，垂直向均采用相控掃描的方式完成垂直向目標(biāo)方位估計(jì)，也即多次發(fā)射接收探測(cè)信號(hào)，進(jìn)行多Ping聯(lián)合三維目標(biāo)探測(cè)，探測(cè)效率較低。本文對(duì)比分析國(guó)內(nèi)外前視聲吶系統(tǒng)方案的局限性，在文獻(xiàn)[6,9]基礎(chǔ)上，提出具備單Ping高效率三維探測(cè)能力的前視聲吶系統(tǒng)解決方案，并進(jìn)行了理論推導(dǎo)、計(jì)算機(jī)仿真和工程實(shí)現(xiàn)驗(yàn)證。該系統(tǒng)方案涉及了一種二維稀疏陣布陣方案以及基于此基陣模型的一維和Vernier聯(lián)合方位估計(jì)(one-dimensional and Vernier combined DOA estimation, ODVCDE)算法；仿真實(shí)驗(yàn)，驗(yàn)證了所提出方法在估計(jì)精度和計(jì)算復(fù)雜度方面優(yōu)于2D波束形成算法；通過水池主被動(dòng)目標(biāo)探測(cè)實(shí)驗(yàn)和外場(chǎng)試驗(yàn)，驗(yàn)證了聲吶系統(tǒng)和相關(guān)算法的工程實(shí)用性。

1 稀疏陣列模型

本文設(shè)計(jì)的前視聲吶接收陣列結(jié)構(gòu)如圖1所示。該陣列為二維面陣設(shè)計(jì)，由3個(gè)均勻線性子陣構(gòu)成，3個(gè)子陣之間陣元間距由一對(duì)互質(zhì)的整數(shù)(5和3)構(gòu)成，其中每一個(gè)子陣由M個(gè)陣元構(gòu)成，為了避免相位模糊陣元間距為半波長(zhǎng)d=λ/2[10-11]。因此，該陣列由MN個(gè)陣元組成，其中，M=48、N=3。

圖1 前視聲吶接收陣列結(jié)構(gòu)Fig.1 Receiving array configuration in forward-looking sonar

信號(hào)模型如圖2所示，在快拍號(hào)t，考慮K個(gè)遠(yuǎn)場(chǎng)窄帶相干回波信號(hào)同時(shí)入射二維接收陣列，第n行第m列陣元接收回波信號(hào)模型可表述為：

〈n(t)〉m,n,t=1,2,…,L

(1)

2πd(mcosηkcosθk+Snsinηk)/λ

(2)

式中：s(t)=exp(j2πft/fs)為主動(dòng)探測(cè)信號(hào)；f和fs分別為工作頻率和采樣頻率；n(t)為加性復(fù)高斯白噪聲且與探測(cè)信號(hào)s(t)統(tǒng)計(jì)獨(dú)立；θk和ηk分別為方位角和俯仰角；λ=c/f為波長(zhǎng)；d為陣元間距。本文中，第1行第1列陣元設(shè)置為相位參考點(diǎn)。

(3)

將式(1)改寫為矩陣形式：

(4)

式中⊙為Hadamard乘積。其中，Ak∈CN×M，

(5)

X(t)∈CN×M，S(t)∈CN×M，N(t)∈CN×M與導(dǎo)向矩陣Ak結(jié)構(gòu)類似，這里不再贅述。

圖2 三維笛卡爾坐標(biāo)系下的目標(biāo)方位Fig.2 Orientation representation of the target in three-dimensional Cartesian coordinates

2 空間方位估計(jì)算法

2.1 空間角Θ估計(jì)

根據(jù)線陣幾何性質(zhì)，空間域的導(dǎo)向波束具有錐形形狀[12]，也即沿同一錐面入射的信號(hào)具有相同的空間角，如圖3所示。

圖3 導(dǎo)向波束的三維波束Fig.3 Three-dimensional beam pattern of steered beams

對(duì)于空間角Θ估計(jì)問題，基于每一個(gè)均勻線性子陣，本文分別設(shè)置第1號(hào)陣元為相位參考點(diǎn)，進(jìn)而本文有Sn=0，也即各子陣在空間角估計(jì)問題上的信號(hào)模型是一致的。第m號(hào)陣元相對(duì)參考陣元延遲可寫作：

〈τk〉m,n=πmdcosΘk/λ

(6)

進(jìn)而，可直接使用一維方法完成空間方位估計(jì)。

依據(jù)空間角分析俯仰角可行解范圍。根據(jù)如圖3所示的三維錐形波束圖案，有：

RY=RcosΘ

(7)

RZ=Rsinη

(8)

式中：R表示陣列參考點(diǎn)與目標(biāo)的距離；RY和RZ為R投影到Y(jié)軸和Z軸方向的距離，滿足：

(9)

sin2Θ≥sin2η

(10)

由于Θ∈[0,π]，η∈[-π/2,π/2]，進(jìn)而有：

(11)

顯然，俯仰角估計(jì)范圍受空間角約束。例如，當(dāng)Θ=140°時(shí)，|η|≤50°，此為下一節(jié)Vernier法解俯仰角的約束條件。

2.2 Vernier法[13-14]俯仰角估計(jì)

由于子陣間陣元間距大于半波長(zhǎng)，目標(biāo)回波相位差可以在[-π,π]范圍內(nèi)展開：

sinηp=(2πp+Δφ12)/[π(S1-S2)]

(12)

sinηq=(2πq+Δφ23)/[π(S2-S3)]

(13)

式中：

(14)

(15)

式中：angle(·)為計(jì)算復(fù)數(shù)相位角的函數(shù);〈·〉*為取復(fù)數(shù)共軛的函數(shù)，〈X〉n表示取矩陣X中的第n行。進(jìn)而由式(12)和式(13)得到目標(biāo)回波方位的候選子集為：

Λp={η12|p∈[?-(S2-S1)/2」,?(S2-S1)/2」]}

(16)

Λq={η23|q∈[?-(S3-S2)/2」,?(S3-S2)/2」]}

(17)

根據(jù)互質(zhì)特性，由候選子集S12和S23可得到目標(biāo)回波方位估計(jì)：

η=S12∩S23

(18)

然而，上述等式成立建立在相位差Δφ12和Δφ23估計(jì)無偏的前提下，為了提高對(duì)于估計(jì)誤差的寬容性，本文定義下面Vernier法估計(jì)的正則化表達(dá)式：

(19)

根據(jù)式(19)可得到俯仰角η估計(jì)。

2.3 基于最小角定理的方位角估計(jì)

根據(jù)式(9)，目標(biāo)距離基陣參考點(diǎn)沿X軸分量為：

(20)

進(jìn)而，方位角可寫作：

(21)

簡(jiǎn)化后：

sec2θ=cos2η/cos2Θ

(22)

由于聲吶功能定位為前視聲吶，目標(biāo)位于一、四、五、八卦限，有：

(23)

整理得cosΘ=cosθcosη，也即最小角定理。對(duì)于方位角θk估計(jì)，由2.1節(jié)得到空間角估計(jì)Θk，由2.2節(jié)得到俯仰角估計(jì)ηk，再由最小角定理可得到方位角估計(jì)：

θk=arccos(cosΘk/cosηk)

(24)

至此，完成目標(biāo)二維方位估計(jì)(包括方位角θk和俯仰角ηk)的相關(guān)算法，總結(jié)相關(guān)算法流程如下：1)對(duì)MN個(gè)通道進(jìn)行一維空間角搜索，得到Θ；2)由式(14)、(15)計(jì)算子陣間相位差Δφ12和Δφ23；3)由式(12)、(13)展開得到俯仰角候選子集Λp和Λq；4)由式(19)解算得到俯仰角估計(jì)η；5)結(jié)合步驟1)和步驟4)結(jié)果，由式(24)解算目標(biāo)方位角θ。

系統(tǒng)復(fù)雜度分析：1)硬件復(fù)雜度。本文所述聲吶硬件復(fù)雜度為MN個(gè)陣元，而Echoscope4G二維面陣聲吶硬件復(fù)雜度為M2個(gè)陣元，這里N<

3 仿真及分析

國(guó)外單Ping三維前視探測(cè)方案中，Echoscope4G聲吶采用二維波束形成算法[5,7]，本節(jié)擬將所提出的算法與波束形成算法進(jìn)行仿真比較，定量分析方位估計(jì)精度和算法計(jì)算復(fù)雜度。定義均方根誤差(root mean square error，RMSE)為：

(25)

考慮一個(gè)方位角85°，俯仰角7°的回波信號(hào)，信噪比以5 dB為步長(zhǎng)從-5 dB變化到20 dB，采用2D波束形成算法和本文提出的算法對(duì)回波空間方位進(jìn)行估計(jì)，B=1 024，MC=1 000次蒙特卡洛仿真實(shí)驗(yàn)結(jié)果如圖4所示。計(jì)算機(jī)配置為Intel(R) Core(TM) i5-8250 u CPU，16 GB RAM。

圖4 仿真結(jié)果對(duì)比信噪比Fig.4 Simulated results versus SNR

如圖4所示為2種算法估計(jì)方差隨信噪比變化曲線，2種算法估計(jì)均方根誤差均隨信噪比增加而單調(diào)遞減，然而低信噪比條件下2種算法估計(jì)均方根誤差稍差于CRLB，且隨著信噪比的增加2種方法逐漸接近CRLB[15]，所提出的方法估計(jì)性能方面優(yōu)于2D波束形成算法。此外，本文統(tǒng)計(jì)2D波束形成算法和所提出的方法的平均計(jì)算時(shí)間分別為32.864 1 s和0.009 8 s，顯然，所提出算法計(jì)算復(fù)雜度方面更優(yōu)，這與2.2節(jié)的理論分析一致。

4 試驗(yàn)及結(jié)果分析

為進(jìn)一步檢驗(yàn)上述算法的有效性，本文依托哈爾濱工程大學(xué)研制的三維前視聲吶分別于哈爾濱工程大學(xué)信道水池、吉林市松花湖地區(qū)開展試驗(yàn)研究，以下試驗(yàn)中三維空間原點(diǎn)均設(shè)為接收基陣第1行第1列陣元。三維前視聲吶關(guān)鍵參數(shù)如表1所示。

4.1 水池聲源被動(dòng)探測(cè)實(shí)驗(yàn)

本節(jié)實(shí)驗(yàn)于哈爾濱工程大學(xué)水聲技術(shù)重點(diǎn)實(shí)驗(yàn)室信道水池中開展，信道水池寬6 m，水深5 m，距離聲吶基陣10 m處放置一個(gè)工作頻率為150 kHz的聲源，聲吶基陣和聲源均置于水面以下2.5 m，聲源發(fā)射和聲吶基陣接收機(jī)使用同步觸發(fā)，聲吶基陣分別豎直安裝從-37°旋轉(zhuǎn)至43°、水平安裝從-71°旋轉(zhuǎn)至72°，每隔1°旋轉(zhuǎn)一次并采集聲源發(fā)射信號(hào)，依托二維波束形成算法和所提出的算法進(jìn)行實(shí)驗(yàn)數(shù)據(jù)處理。水池實(shí)驗(yàn)場(chǎng)景和三維前視聲吶安裝布放實(shí)景如圖5所示。

表1 三維前視聲吶關(guān)鍵參數(shù)

圖5 水池實(shí)驗(yàn)場(chǎng)景和前視聲吶布放示意Fig.5 Deployment and rotating of the forward-looking sonar

如圖6所示為旋轉(zhuǎn)至14°時(shí)各通道相位輸出曲線，由圖可見，隨著陣元編號(hào)增加，相位連續(xù)變化。為了更清晰地描述相位差曲線，基陣從-37°旋轉(zhuǎn)至43°，子陣間相位差曲線在圖7中給出。由圖7可見，子陣間相位差隨基陣旋轉(zhuǎn)角度連續(xù)變化，基本符合理論上的線性關(guān)系。

圖6 子陣輸出相位箱線Fig.6 Phase boxplot of subarray output

圖7 子陣相位差箱線Fig.7 Phase differences boxplot between two subarrays

解析聲源方位后的DOA估計(jì)如圖8所示。由圖可見，估計(jì)的DOA隨基陣旋轉(zhuǎn)角度連續(xù)變化，基本符合理論上的線性關(guān)系。水池實(shí)驗(yàn)驗(yàn)證了第2節(jié)相關(guān)算法的有效性。值得注意的是，由于接收陣元存在幅度、相位等誤差，估計(jì)結(jié)果與理論值存在一定誤差。

圖8 方位估計(jì)結(jié)果Fig.8 Estimated results of DOA

4.2 水池目標(biāo)主動(dòng)探測(cè)實(shí)驗(yàn)

本實(shí)驗(yàn)同樣于信道水池中開展，聲吶基陣置于水面以下2.5 m，距離聲吶基陣2.5、5.2、6.1 m處分別放置3個(gè)直徑280 mm的塑料球目標(biāo)，各目標(biāo)距離水面深度分別為2.0、2.3、2.5 m。如圖9所示為應(yīng)用波束形成技術(shù)形成的聲吶圖像。以聲吶基陣原點(diǎn)所在位置為參考原點(diǎn)，目標(biāo)估計(jì)結(jié)果圖如圖10所示。

圖9 水池實(shí)驗(yàn)場(chǎng)景聲圖像Fig.9 Acoustic image of tank experiment scene

由圖10可見，依據(jù)本文所提出的算法能夠?qū)λ刂心繕?biāo)方位角、俯仰角進(jìn)行估計(jì)，聯(lián)系回波到達(dá)時(shí)間，可完成目標(biāo)的三維定位。

圖10 水池實(shí)驗(yàn)估計(jì)結(jié)果Fig.10 Estimated results in tank experiment

4.3 外場(chǎng)試驗(yàn)

外場(chǎng)試驗(yàn)于2019年8月在吉林市松花湖進(jìn)行。3個(gè)目標(biāo)分別呈三角形布放，每個(gè)目標(biāo)懸浮于水中；水面懸掛一個(gè)浮標(biāo)；水底設(shè)置一個(gè)重物用于定位目標(biāo)。3個(gè)目標(biāo)均置于水面以下5 m位置處，聲吶系統(tǒng)固定安裝在舷側(cè)并置于水面以下1 m位置處，采用光學(xué)激光測(cè)距儀測(cè)量3個(gè)目標(biāo)距離聲吶分別為157.1、160.7、165.7 m。通過聲吶發(fā)射探測(cè)信號(hào)、接收回波數(shù)據(jù)可對(duì)目標(biāo)方位進(jìn)行估計(jì)。外場(chǎng)湖試場(chǎng)景及前視聲吶、目標(biāo)安裝布放實(shí)景如圖11所示。

圖11 外場(chǎng)湖試場(chǎng)景及三維前視聲吶、目標(biāo)安裝布放實(shí)景Fig.11 The scene of field lake experiment and deployment of the three-dimensional forward-looking sonar and targets

如圖12為應(yīng)用波束形成技術(shù)形成聲吶圖像，可以顯著地看出有3個(gè)目標(biāo)。以聲吶基陣原點(diǎn)所在位置為參考原點(diǎn)，應(yīng)用本文中所提方法，目標(biāo)三維空間位置估計(jì)結(jié)果如圖13所示。

圖12 外場(chǎng)試驗(yàn)場(chǎng)景聲圖像Fig.12 Acoustic image of field experiment scene

圖13 外場(chǎng)湖試估計(jì)結(jié)果Fig.13 Estimated results in field lake experiment

由圖13可見，依據(jù)本文所提出的算法能夠?qū)ν鈭?chǎng)湖試布放的目標(biāo)進(jìn)行方位角、俯仰角估計(jì)，依據(jù)回波到達(dá)時(shí)間，可完成目標(biāo)的三維定位。

5 結(jié)論

1)相比Echoscope4G三維前視聲吶基于二維面陣的方案，所提出的方法具有更低的系統(tǒng)復(fù)雜度和實(shí)現(xiàn)成本，僅僅利用MN個(gè)陣元即可完成目標(biāo)的三維探測(cè)，陣元數(shù)僅為二維方案的N/M，計(jì)算復(fù)雜度為二維方案的N/MB。

2)仿真實(shí)驗(yàn)表明，相比2D波束形成算法，所提出的JSODVM算法具有更高的方位估計(jì)精度、更低的計(jì)算復(fù)雜度。結(jié)合三維前視聲吶的水池、湖上試驗(yàn)進(jìn)一步驗(yàn)證了上述算法能夠準(zhǔn)確估計(jì)目標(biāo)在三維空間中的方位，工程實(shí)用性強(qiáng)。

本文更關(guān)注于相關(guān)算法的研究，在未來的研究中，將針對(duì)陣列校準(zhǔn)這一問題開展深入研究。