基于序列輪廓線的多波束水體影像目標(biāo)三維重建

陳建兵，劉洪霞,2*，李東輝，梅 賽，馮成凱，陽凡林,6

（1.山東科技大學(xué) 測繪與空間信息學(xué)院，山東 青島 266590；2.山東科技大學(xué) 數(shù)學(xué)與系統(tǒng)科學(xué)學(xué)院，山東 青島 266590；3.青島環(huán)海海洋工程勘察設(shè)計院，山東 青島 266061；4.中國地質(zhì)調(diào)查局青島海洋地質(zhì)研究所，山東 青島 266071；5.江蘇海洋大學(xué) 海洋技術(shù)與測繪學(xué)院，江蘇 連云港 222005；6.自然資源部海洋測繪重點(diǎn)實(shí)驗室，山東 青島 266590）

水體數(shù)據(jù)是多波束聲吶采集的最原始的數(shù)據(jù)，包含了波束從換能器至海底完整的聲學(xué)信息[1]，將其按一定方式成像處理可生成多視角的水體影像（Water Column Image，WCI）。水體影像廣泛應(yīng)用于探測魚群分布、水下沉船、海底熱液噴口、冷泉羽狀流以及海洋內(nèi)波等[2-5]物理目標(biāo)物，在水下地形測量、海事搜救和打撈、軍事及地質(zhì)活動等方面具有極其重要的應(yīng)用價值[1]。隨著多波束水體影像在諸多領(lǐng)域的應(yīng)用，相關(guān)學(xué)者對水體數(shù)據(jù)的精細(xì)化處理展開了進(jìn)一步研究。Werf A[6]通過多種成圖方式分析水體數(shù)據(jù)，并精確測量出水下沉船的最淺點(diǎn)，這種方式對操作人員的專業(yè)知識背景要求較高。Liu H等[7]提出一種水體影像自適應(yīng)軟閾值去噪方法，可有效提高影像信噪比，降低水體影像識別難度。然而，目前對水體目標(biāo)的識別與分析通常在二維水體影像中進(jìn)行，難以滿足目標(biāo)物三維成像的需求，因此一些學(xué)者開始將三維可視化技術(shù)應(yīng)用到多波束水體數(shù)據(jù)的處理過程中。Urban P等[8]根據(jù)羽狀流水體影像構(gòu)建三維體素網(wǎng)格，通過體積反向散射強(qiáng)度轉(zhuǎn)換來計算羽狀流氣體通量，極大地提高了羽狀流定量分析的精度。李東輝等[9]對水體影像目標(biāo)識別進(jìn)行分析研究，提出一種目標(biāo)物自動提取方法，并完整的提取出沉船采樣點(diǎn)，為三維重建提供了點(diǎn)云數(shù)據(jù)獲取方法。

三維重建技術(shù)可快速獲取物體的數(shù)字化模型，從而實(shí)現(xiàn)基于現(xiàn)實(shí)的3D虛擬顯示效果，在醫(yī)學(xué)仿生、逆向工程、地礦建模、游戲影視等領(lǐng)域產(chǎn)生巨大應(yīng)用價值。用于三維重建的數(shù)據(jù)形式通常有點(diǎn)云、序列影像、平行輪廓線，此三種數(shù)據(jù)皆可由水體數(shù)據(jù)進(jìn)行相應(yīng)的處理得到，然而關(guān)于水體數(shù)據(jù)三維重建的相關(guān)文獻(xiàn)資料卻鮮有公開發(fā)表，一方面由于多波束聲吶采集的目標(biāo)點(diǎn)云空間分布具有嚴(yán)重的斷層性、復(fù)雜性和非均勻性，直接利用點(diǎn)云重建目標(biāo)物效果不佳；另一方面，水體影像中包含有未改正的載體姿態(tài)、吃水、聲速剖面等方面的誤差，根據(jù)序列影像重建易導(dǎo)致目標(biāo)物變形扭曲。近些年，隨著逆向工程行業(yè)的迅速發(fā)展，基于斷層數(shù)據(jù)（如醫(yī)學(xué)上CT、MRI掃描數(shù)據(jù)）的三維重建方法逐漸成為研究熱點(diǎn)，輪廓線法是一種利用空間斷層數(shù)據(jù)重建的經(jīng)典方法，由Keppel E[10]在1975年提出，后來的研究者在輪廓線法基礎(chǔ)上改進(jìn)了約束條件的判斷形式，形成最大體積法、最小表面法、最短對角線法、同步前進(jìn)法和切開縫合法等[10-14]，這些方法在醫(yī)學(xué)、地質(zhì)等行業(yè)類似問題中展現(xiàn)出良好的適用性。

本文提出一種基于序列輪廓線的多波束水體目標(biāo)物三維重建算法，通過點(diǎn)云分割對多波束單Ping（一個發(fā)射接收周期）數(shù)據(jù)做進(jìn)一步細(xì)分處理并提取輪廓，在提取的序列輪廓線基礎(chǔ)上，提出一種空間格網(wǎng)劃分方法，通過局部優(yōu)化算法建立頂點(diǎn)間拓?fù)潢P(guān)系，對空間格網(wǎng)形成面片填充，從而實(shí)現(xiàn)按“點(diǎn)—線—面”的順序逐步重建水體目標(biāo)。

1 多波束水體目標(biāo)成像特征

多波束聲吶通過連續(xù)接收每條波束的回波信號，可在聲照射扇區(qū)內(nèi)形成高密度采樣，從而繪制水體影像。波束攜帶的能量會穿透整個水體到達(dá)海底，但聲波在水下傳播速度遠(yuǎn)比不上光速，聲信號旅行時不可忽略，換能器連續(xù)發(fā)射聲波的時間須限制在此前發(fā)射的聲信號全部返回之后，因此水體數(shù)據(jù)是多波束聲吶在航線上對水體的間斷式掃描的結(jié)果。均勻、密集、穩(wěn)定性好的采樣數(shù)據(jù)是最理想的三維重建數(shù)據(jù)，而水體數(shù)據(jù)空間分布特征會對三維重建算法產(chǎn)生一定的影響。

1.1 水體目標(biāo)成像機(jī)理

多波束聲吶在工作時，以換能器基陣為中心，在一個平面內(nèi)同時形成數(shù)百甚至上千條聲學(xué)波束，這些波束攜帶的能量從換能器發(fā)射后，經(jīng)水體傳播至海底，波束在水下的傳播過程中（包括水體、海底以及其它物體）會產(chǎn)生聲散射現(xiàn)象，這些后向散射信號再由接收換能器接收并記錄在水體數(shù)據(jù)中，從而完成水體數(shù)據(jù)一Ping觀測。換能器發(fā)射基元會持續(xù)發(fā)射聲波，伴隨載體的移動，可得到一個類似三棱柱形狀的水體掃描區(qū)域（圖1-a）。

垂向扇形圖與航向堆疊圖是多波束水體影像最常見的兩種形式。垂向扇形圖以Ping為單位進(jìn)行繪制，水體回波信號以極坐標(biāo)的形式存儲在水體數(shù)據(jù)文件中，計算采樣點(diǎn)像素坐標(biāo)，并將強(qiáng)度區(qū)間轉(zhuǎn)換為灰度區(qū)間，即可得到一幅垂向扇形圖（圖1-b）。航向堆疊圖可整體上觀察水體情況，將一Ping水體數(shù)據(jù)按反向散射強(qiáng)度堆疊為一列，對當(dāng)前測線所有Ping進(jìn)行堆疊顯示即可獲得航向堆疊圖（圖1-c），圖1-c中藍(lán)線標(biāo)識位置為圖1-b堆疊結(jié)果。

圖1 多波束水體目標(biāo)成圖原理

1.2 水體數(shù)據(jù)空間分布特征

在不改變多波束聲吶采樣率的前提下，水體采樣點(diǎn)在每個波束空間上均勻分布。由于波束在聲照射扇面內(nèi)是以換能器為中心，按固定的角度接收各個方向的聲信號，因此越靠近換能器采樣點(diǎn)越集中，隨著水深的增加波束之間的幾何縫隙逐漸擴(kuò)大（圖2-a），在一Ping水體影像中，邊緣波束海底附近的采樣點(diǎn)最為稀疏。在三維空間中水體數(shù)據(jù)的復(fù)雜性進(jìn)一步擴(kuò)大，通常水體數(shù)據(jù)在Ping平面內(nèi)采樣率高于航線方向，相鄰Ping之間的空白區(qū)域較大，整體上呈斷層式分布，另外由于載體姿態(tài)、吃水以及聲速剖面等諸多不可控因素的影響，多波束沿航線采集的多Ping水體影像之間并非完全平行分布（圖2-b）。因此多波束水體數(shù)據(jù)在空間中表現(xiàn)出嚴(yán)重的斷層性、復(fù)雜性與非均勻性，對三維重建算法造成不利影響。

圖2 多波束水體影像數(shù)據(jù)分布示意圖

2 目標(biāo)序列輪廓線三維重建算法

多波束聲吶獲取的水體數(shù)據(jù)空間分布情況較復(fù)雜，易導(dǎo)致算法重建失敗，為提高水體目標(biāo)三維模型的質(zhì)量，本文基于最短對角線法局部優(yōu)化準(zhǔn)則，通過Ping內(nèi)輪廓提取與Ping間面片填充實(shí)現(xiàn)按照“點(diǎn)—線—面”的順序逐步分析處理數(shù)據(jù)，完成水體目標(biāo)重建。首先根據(jù)點(diǎn)云數(shù)據(jù)生成目標(biāo)整體的序列輪廓線，Ping是水體數(shù)據(jù)的基本組織單位，在Ping的基礎(chǔ)上對點(diǎn)云做進(jìn)一步細(xì)分處理并提取輪廓線，以避免輪廓線失真；在輪廓線面片填充時，由于相鄰Ping輪廓線數(shù)量一致性難以保證，可能會出現(xiàn)多條輪廓線都與同一條輪廓線進(jìn)行面片填充，即“分支”問題，這種情況可通過拆分輪廓線轉(zhuǎn)換成多個一一對應(yīng)的輪廓線填充形式，再使用輪廓線中心平移與劃分空間格網(wǎng)的方法，完成“一對一”的輪廓線面片填充，最終形成目標(biāo)物表面模型。

2.1 最短對角線法原理

輪廓線法是以一系列大致平行的輪廓線作為建模數(shù)據(jù)，相鄰輪廓線之間以一系列細(xì)小的三角形面片填充，從而生成整個表面模型。在生成三角形面片的過程中，需要不斷進(jìn)行判斷與調(diào)整，使得模型整體勻稱、視覺效果好。最短對角線法作為一種局部優(yōu)化準(zhǔn)則，具有計算量小、構(gòu)網(wǎng)均勻、簡便快速的優(yōu)點(diǎn)，應(yīng)用廣泛。該算法原理如圖3所示，假設(shè)多波束相鄰兩Ping的輪廓線分別為P和Q，在四邊形pipi+1qi+1qi中，兩對角線長度分別為r1和r2，若r1＜r2，則將三角形piqiqi+1加入三角網(wǎng)，否則，將三角形pipi+1qi加入三角網(wǎng)，同時，四邊形前進(jìn)構(gòu)建下一個三角形，逐步擴(kuò)展三角網(wǎng)直到重新返回起點(diǎn)。

圖3 最短對角線法示意圖

2.2 多波束Ping內(nèi)點(diǎn)云輪廓提取

多波束水體數(shù)據(jù)是按Ping進(jìn)行采集、存儲和處理的，但在提取水下目標(biāo)序列輪廓線時，以Ping為單位提取可能會導(dǎo)致輪廓線失真，因為每Ping代表水體的一個斷面測量，考慮到復(fù)雜結(jié)構(gòu)的水下目標(biāo)的情況，一Ping數(shù)據(jù)中可能同時包含一個物體多個子結(jié)構(gòu)的剖面（如圖4-a中水下沉船的纜繩與船身），對多個不連續(xù)的剖面生成同一條輪廓線，會降低輪廓線的精度，因此需在輪廓線提取之前對Ping內(nèi)點(diǎn)云進(jìn)行分割。本文根據(jù)點(diǎn)云之間的歐式距離進(jìn)行聚類，將距離較近的點(diǎn)劃分為一類，然后再以類為單位提取輪廓（圖4-b）。依次處理所有Ping點(diǎn)云從而獲得整個目標(biāo)物的序列輪廓線。

圖4 Ping內(nèi)輪廓線生成

2.3 多波束相鄰Ping間填充

相鄰Ping輪廓線對應(yīng)關(guān)系分為“分支”和“一對一”兩種形式，而“分支”問題會導(dǎo)致模型面片交叉、數(shù)據(jù)冗余，難以控制模型的質(zhì)量，需對“分支”問題進(jìn)行合理轉(zhuǎn)化，使面片填充統(tǒng)一在“一對一”的形式下進(jìn)行。

2.3.1 相鄰Ping間“分支”問題“分支”問題是由于水下目標(biāo)在航線方向的幾何結(jié)構(gòu)發(fā)生變化所導(dǎo)致的輪廓線數(shù)量不一致，合理解決“分支”問題是檢驗算法有效性的關(guān)鍵。目前現(xiàn)有文獻(xiàn)提供的解決方案主要有兩種，一是利用插值算法在分支輪廓線與母輪廓線之間添加過渡輪廓線[15-16]；二是創(chuàng)造輔助線將母輪廓線分割成幾塊，然后再對分割出來的輪廓線進(jìn)行填充[17-18]。對于多波束水體數(shù)據(jù)，由于數(shù)據(jù)量較大，添加過渡輪廓線將進(jìn)一步增加計算量，降低算法效率，且插值算法的精度不易驗證；創(chuàng)造輔助線切分凹輪廓線時，切分結(jié)果與“分支”輪廓個數(shù)可能出現(xiàn)差異，且分支數(shù)目越多，處理起來越困難。本文通過在母輪廓線上添加斷點(diǎn)對輪廓線進(jìn)行拆分，將“分支”問題轉(zhuǎn)化為多個“一對一”填充，優(yōu)點(diǎn)是以點(diǎn)拆分輪廓線比輔助線分割更靈活，且不會增加額外的數(shù)據(jù)量。

圖5中第iPing的兩條輪廓線都與第i+1 Ping輪廓線形成填充關(guān)系，通過插入兩個斷點(diǎn)可將第i+1 Ping輪廓線拆成兩條完整的輪廓線，從而將“二對一”填充問題轉(zhuǎn)換成兩個“一對一”填充。為使得輪廓線拆分合理、視覺效果更好，在第i+1 Ping輪廓線上，按照最小距離差的方法插入斷點(diǎn)，具體方法如下：遍歷第i+1 Ping輪廓線上所有的點(diǎn)，分別計算到第iPing兩條輪廓線的最短距離r1和r2，在|r1-r2|達(dá)到最小處插入斷點(diǎn)，利用這些斷點(diǎn)將母輪廓線拆分成兩段。然后再分別計算拆分的兩條輪廓線中心點(diǎn)與第iPing中兩條輪廓線中心點(diǎn)的距離，對這兩Ping輪廓線進(jìn)行匹配，形成兩個“一對一”的填充關(guān)系。

圖5 拆分輪廓線

2.3.2 “一對一”填充 對于相似性較好的輪廓線，可以直接以最短對角線法作為約束準(zhǔn)則進(jìn)行填充。但是，當(dāng)兩輪廓線的位置和形狀差異過大時，可能會出現(xiàn)圖6中的錯誤，從而使算法失效。

圖6 最短對角線法失效的兩種情況

輪廓線錯位易導(dǎo)致圖6-a中的圓錐面，需在面片填充之前進(jìn)行輪廓線中心平移，使兩輪廓線的中心在Ping平面法線方向重合，然后按最短對角線法進(jìn)行面片填充（圖7-b），最后再將點(diǎn)云的拓?fù)潢P(guān)系映射回平移之前的輪廓線上（圖7-c）。

圖7 輪廓中心平移示意圖

圖6-b中的錯誤是由于上下輪廓線相似性不好，且拐點(diǎn)較多，最短對角線法在拐點(diǎn)出現(xiàn)的輪廓線上容易出錯，在兩條輪廓線上不能同步前進(jìn)，算法繼續(xù)進(jìn)行則出現(xiàn)圖8-a中的面片自相交情況。在拐點(diǎn)較多的輪廓線上生成面片時，先在輪廓線之間劃分出一系列空間格網(wǎng)，如圖8-b所示，對每一塊空間格網(wǎng)可以看作是一個上下邊為曲線的空間四邊形，把復(fù)雜輪廓線的填充轉(zhuǎn)化為對每一個空間四邊形的填充，這些空間四邊形上下邊的曲率變化相對簡單，即使出現(xiàn)圖6-b中不同步構(gòu)網(wǎng)的情況，也能將其限制在該格網(wǎng)內(nèi)，避免繼續(xù)傳播影響相鄰的格網(wǎng)。

圖8 劃分空間格網(wǎng)示意圖

劃分格網(wǎng)的具體方法：在下層的輪廓線中均勻插入n個點(diǎn)P={pi|i∈[1,n]}，將該輪廓線等分成n份，在上層輪廓線上插入n個點(diǎn)Q={qi|i∈[1,n]}，且qi是上層輪廓線上與pi距離最近的點(diǎn)，連接這n對斷點(diǎn)便在兩輪廓線間形成n個空間格網(wǎng)，如圖8-b，最后再采用最短對角線法在每一個格網(wǎng)上生成三角形面片，從而在兩輪廓線之間形成完整表面，如圖8-c。

綜上所述，首先對多波束水體點(diǎn)云數(shù)據(jù)提取序列輪廓線，然后按航向順序在相鄰Ping輪廓線之間填充面片，從而形成水體目標(biāo)表面模型，最后通過模型渲染提升視覺效果。具體流程如圖9所示。

圖9 算法流程

3 實(shí)驗與分析

本次實(shí)驗采用Kongsberg公司提供的實(shí)測數(shù)據(jù)。數(shù)據(jù)采集地點(diǎn)為加拿大悉尼市瑪格麗特公主海濱公園附近海域（48°43.323'N，123°21.339' W），使用Kongsberg EM3002型多波束聲吶系統(tǒng)對G.B.Church號水下沉船進(jìn)行掃描，此船全長53 m，寬11 m，沉沒于平坦砂質(zhì)海底（圖10），周圍水深24～27 m。數(shù)據(jù)采集時，多波束聲吶采樣率設(shè)置為7 146.5 Hz，波束數(shù)和波束開角分別為160°和130°，表層聲速為1 477 m/s，船體橫跨145Ping數(shù)據(jù)。沉船形狀復(fù)雜、細(xì)節(jié)保留豐富，具有較高的代表性。

圖10 多波束水深圖和海底反向散射強(qiáng)度圖像

通過解析多波束水體數(shù)據(jù)并繪制成圖，可得到多視角下的水體影像（圖11），采用李東輝等[9]介紹的水體影像目標(biāo)自動提取算法從原始水體數(shù)據(jù)中提取沉船點(diǎn)云，逐Ping對采樣點(diǎn)進(jìn)行聚類分割，含有沉船的145Ping數(shù)據(jù)被細(xì)分為416類，從而提取出416條輪廓線，沉船整體的序列輪廓線如圖11左下圖所示。

圖11 沉船水體影像及序列輪廓線

根據(jù)圖11中的序列輪廓線，分別采用最短對角線法和本文算法進(jìn)行重建，結(jié)果如圖12所示，直接使用最短對角線法所建的模型中，圓錐面和自相交面片情況導(dǎo)致模型中出現(xiàn)多處不合理結(jié)構(gòu)，如桅桿和船尾部位，這體現(xiàn)了傳統(tǒng)輪廓線算法對于復(fù)雜輪廓線處理的局限性。本文算法針對模型中大量出現(xiàn)的圓錐面和自相交面片進(jìn)行改進(jìn)，在可控的小范圍內(nèi)生成首尾相連的三角網(wǎng)，降低了處理復(fù)雜輪廓線的難度，同時，模型結(jié)構(gòu)得到改善，形狀更加合理。本文算法在Visual studio 2013平臺上使用C++語言實(shí)現(xiàn)，并借助VTK（Visualization Toolkit）對結(jié)果進(jìn)行可視化。

圖12 輪廓線算法建模結(jié)果對比圖

為進(jìn)一步驗證本文算法對水體數(shù)據(jù)的適用性，使用α-shape[19]算法與拋雪球算法[20]對沉船點(diǎn)云進(jìn)行重建。α-shape算法與拋雪球算法是經(jīng)典的三維重建算法，其中α-shape算法是在Delaunay三角剖分基礎(chǔ)上，根據(jù)指定的α值選取剖分結(jié)果中的α簇作為表面模型；拋雪球法需要將輸入的點(diǎn)云數(shù)據(jù)轉(zhuǎn)換成體數(shù)據(jù)，按照分布函數(shù)計算每個體素對投影平面像素點(diǎn)的貢獻(xiàn)值，從而實(shí)現(xiàn)點(diǎn)云數(shù)據(jù)的三維重建。對上述兩種算法進(jìn)行多次試驗以挑選最佳參數(shù)值，圖13和表1展示上述兩種算法與本文算法重建效果及統(tǒng)計結(jié)果。

圖13 不同算法重建結(jié)果對比

表1 不同算法重建結(jié)果統(tǒng)計表

從圖11的沉船照片中可以看到，桅桿尺寸較小，水體數(shù)據(jù)在此位置采樣點(diǎn)稀疏，α-shape算法將桅桿還原成螺旋形狀，而拋雪球算法重建結(jié)果中，桅桿粗糙不平，有明顯凸起。在纜繩與船身連接位置，由于二者距離較近，α-shape算法與拋雪球算法均將纜繩與船身連為一體。在船艙位置有一排狹長開口，開口排列方向大致與Ping平面方向垂直，這些開口在α-shape算法與拋雪球算法的重建結(jié)果中均沒有體現(xiàn)。從圖13的3種算法的對比中可看出，本文算法更符合目標(biāo)物形狀，對桅桿、纜繩、船艙開口等細(xì)節(jié)還原準(zhǔn)確，模型中的面片大小均勻，整體視覺效果較好。對三種模型的數(shù)據(jù)量和運(yùn)行時間進(jìn)行統(tǒng)計，從表1中可以看出α-shape總面片和非流形邊數(shù)量最多，非流形邊是3個或以上面片的交線，非流形邊大量出現(xiàn)表明模型存在面片交叉問題；而拋雪球算法面片數(shù)量和運(yùn)行時間最少，但模型特征還原不夠精準(zhǔn)，表面結(jié)構(gòu)冗余。

另外，考慮到多波束系統(tǒng)在連續(xù)工作的過程中，受海況或系統(tǒng)自身問題影響，時有出現(xiàn)某Ping水深和水體數(shù)據(jù)異常，這些異常數(shù)據(jù)不符合目標(biāo)形狀，可事先檢查異常Ping并剔除，在剔除異常Ping之后水體數(shù)據(jù)整體仍為斷層式分布。本文算法針對空間斷層數(shù)據(jù)進(jìn)行組織與優(yōu)化，算法的正確運(yùn)行不依賴于數(shù)據(jù)的空間密度和分布均勻性，具有較強(qiáng)的魯棒性，本文算法對剔除異常Ping之后的數(shù)據(jù)仍然有效。因此，本文算法對多波束水體數(shù)據(jù)三維重建是合理且實(shí)用的。

4 結(jié)論

三維重建技術(shù)可解決多波束水體影像對目標(biāo)成像不夠全面、直觀的弊端，本文充分考慮數(shù)據(jù)的斷層性、復(fù)雜性與非均勻性對重建算法的影響，基于傳統(tǒng)輪廓線法建模原理，提出了按“點(diǎn)—線—面”順序?qū)λw目標(biāo)逐步重建的方法，針對幾何結(jié)構(gòu)復(fù)雜的目標(biāo)，先進(jìn)行Ping內(nèi)點(diǎn)云聚類再提取輪廓線，可生成與真實(shí)情況契合度較好的輪廓，并通過輪廓線拆分、輪廓中心平移與劃分空間格網(wǎng)的步驟能有效解決復(fù)雜輪廓線間面片填充的問題，避免面片交叉、連接錯誤的現(xiàn)象。實(shí)驗證明該算法穩(wěn)健性強(qiáng)，模型準(zhǔn)確度高，相比于傳統(tǒng)輪廓線算法以及目前常用的點(diǎn)云三維重建方法，本文算法整體上無錯誤結(jié)構(gòu)，且對小尺度上的細(xì)節(jié)特征也有精準(zhǔn)還原，在處理復(fù)雜的水體目標(biāo)物時展現(xiàn)了良好的適用性，可在水下目標(biāo)精準(zhǔn)識別與定量分析工作中產(chǎn)生較大的實(shí)用價值。