劉 毅, 王 磊, 王勝利, 胡亮亮, 高文龍
(1.山東科技大學測繪與空間信息學院, 青島 266590; 2.自然資源部,海洋測繪重點實驗室, 青島 266590;3.山東科技大學海洋科學與工程學院, 青島 266590)
多波束探測技術已經成為目前水下地形探測的主要方法[1]。近些年隨著近岸淺水調查及淺水水下地形項目需求的日益增長,淺水多波束測深系統得到了廣泛的應用[2-3]。在水下地形測量過程中,淺水多波束安裝于船體側舷容易受到噪聲、復雜的海洋環境及設備參數設置不合理等多種因素的影響,測深數據中難免含有粗差[4]。
為了獲得精確的水下地形,需要剔除多波束數據粗差,現有的濾波方法有人工剔除粗差和自動濾波兩種。人工剔除粗差方法準確度較高,但是需要消耗大量的勞動力,且效率較低,無法滿足大量數據處理[5]。因此許多專家學者提出自動濾波方法。李宏武等[6]利用基于最小二乘支持向量機(least squares support vector machines,LS-SVM)算法探測深度異常值;Calder等[7]提出CUBE(combined uncertainty and bathymetry estimator)算法;趙荻能等[8]提出了一種CUBE、曲面濾波參數聯合選優流程;趙祥鴻等[9]提出利用BP神經網絡剔除多波束測深數據粗差;楊安秀等[10]提出面向多波束測深數據的雙向布料模擬自動濾波方法。
Zhang等[11]提出了基于布料模擬濾濃(cloth simulation filtering,CSF)算法對陸地激光點云濾波,其方法需要調整參數較少,精度高,適用不同復雜地形,給激光點云濾波提供了新的思路。本文中則利用CSF算法剔除多波束點云粗差,如果將CSF算法直接用于多波束點云濾波,會導致對凸起礁石和凹陷溝槽等水下地形過度濾波。因此在CSF算法基礎上,提出一種基于改進的CSF算法多波束數據濾波來解決對于復雜海底地形過度濾波問題。首先,對多波束點云數據進行點云分割;然后,根據點云數據的面積和標準差調整每塊布料濾波算法的參數;其次,構建布料濾波模型剔除粗差點;最后,將分割后的點云拼接生成水下地形點云數據,解決水下地形起伏較大時過度濾波的問題。多波束測量系統是在不斷運動狀態下進行的一種動態測量,由于測量平臺受海浪起伏、風、流等海洋環境效應的干擾,多波束測深數據采集過程難免出現假信號,從而采集的數據與實際地形存在差異,這就是海洋測量信息獲取過程中的粗差問題[12]。因此,粗差的判別非常重要,通常水下地形是連續變化,本文中利用布料模擬連續變化的水下地形,將離散點和異常線狀點云判定為粗差。
Zhang等[10]提出的CSF算法是利用布料模擬技術進行濾波。布料模擬技術可以想象成一塊布放在一個水下地形的上方,布塊在重力作用下落下,如果布料足夠柔軟,則可以貼在地形的表面上,最后布塊的形狀就是地形。基于布料模擬技術,Zhang等[10]開發了CSF算法用于激光雷達提取地面點,首先將原始點云翻轉,模擬一塊布落到翻轉后的地面上,最終確定布塊的形狀,并將其形狀作為點云分類的依據。本文中則利用CSF算法剔除多波束點云數據粗差點(圖1)。

圖1 布料模擬濾波算法原理
在布料模擬過程中,布料是由具有質量和互連的粒子組成的網格(圖2)。網格節點上的粒子沒有大小,質量相同且不變,粒子所在空間中的位置決定了布料的形狀。布料格網中連接兩個粒子的結構稱之為“虛擬彈簧”并且遵循胡克定律。

圖2 布料格網模型
布料格網中粒子的位置和速度由作用于粒子上的力決定,由牛頓第二定律可得,粒子位置與力之間的關系為
(1)
式(1)中:X(t)表示粒子在t時刻的位置;Fe(X,t)表示粒子所受外力,由粒子自身的重力G和與地形接觸的碰撞力Fc組成,粒子在運動方向上碰撞到一些物體則會產生碰撞力;Fi(X,t)表示粒子在X位置和t時間的內力,由粒子間相互作用的“虛擬彈簧”產生(圖3)。

圖3 虛擬彈簧
建立布料格網模型時,首先,粒子的運動被限制在垂直方向,因此可以通過比較粒子和地形的高差來檢測碰撞;然后,當一個粒子到達地形時,這個粒子就被設定為不可移動;最后,將力的分成離散的兩個步驟,先計算重力對粒子的影響,后根據“虛擬彈簧”產生的內力修改粒子的位置(圖4)。

圖4 布料格網節點位移計算
僅從重力計算每個粒子的位置,即內力得零時可得
(2)
式(2)中:m表示粒子質量;g表示重力加速度;X(t)表示t時刻節點粒子位置;Δt表示時間長。當已知節點粒子初始位置與時間步長,則可以解算出格網節點的位置。
為了限制粒子在反轉的點云空隙區域的位移,需要考慮到粒子間的內力。因為內力作用,粒子會停留在網格中,并且回到初始位置。遍歷所有“虛擬彈簧”,對于每個彈簧比較形成這個彈簧的兩個粒子之間的高差。如果兩個粒子都是可以移動,在相反方向上移動相同的量;如果一個粒子不可移動,另一個將被移動。每個粒子的位移計算公式為
(3)
式(3)中:d為粒子的位移量;當粒子可以移動時b=1,不可移動時b=0;pi為p0的相鄰粒子的位置,n為點進行標準化到垂直方向上的單位向量(0,0,1)T。
1.2.1 水下地形點云分割
水下地形是連續變化的,由于海洋水體的運動會不斷地侵蝕水下地形,并且一條測線覆蓋范圍相對較廣,水下會有凹陷溝槽或者凸起礁石等地形。如果對一整條測線都使用相同參數進行CSF算法濾波,則容易在地形起伏較大處剔除水下地形點。因此需要對水下地形點云分割,當分割后點云面積大于S1=10 000 m2,則再對該點云進行分割,分割完成后可以針對不同區域調整參數進行CSF算法濾波(圖5)。

圖5 點云分割
1.2.2 自適應參數選取
水下地形數據濾波時,CSF參數調整是決定濾波質量的關鍵因素。CSF算法需要調整參數為6個分別為:是否有斜坡St、布料格網分辨率GR、布料剛性RI、時間步長dT、距離閾值hcc、迭代次數n。
(4)

基于改進的CSF算法流程(圖6),可以概括為以下步驟。

圖6 算法流程
步驟1基于PCL對點云數據進行分割,使每個點云面積小于S1。
步驟2基于點云數據的標準差設置濾波參數:是否有斜坡St;布料剛性RI;布料格網分辨率GR。
步驟3計算布料在重力作用下格網節點高程值,計算節點之間彈簧力作用下的節點高程值。
步驟4重復步驟(3)當所有迭代次數大于設定n次,或者未達到n次且當所有節點粒子迭代后均未移動則迭代停止。
步驟5計算水下地形點云點與節點的高程差,若高程差小于或等于距離閾值hcc,則被分類為水下地形點,否則為粗差點。
步驟6將所有水下地形點云拼接。
實驗中采用兩處多波束點云數據驗證該算法有效性,分別為溫州測區與青島測區中條帶數據,其中青島測區平均水深41.15 m,數據共包含水深點9 204 749,條帶跨度較大,水下地形起伏較大;溫州測區數據平均水深9.25 m,數據共包含水深點3 551 232,條帶跨度較小,水下地形起伏較小。選擇這兩數數據主要由于分別代表了兩種典型的多波束點云數據粗差,兩次數據粗差點均由噪聲引起且無法避免,其中溫州數據粗差集中于條帶中央,呈線狀,而青島數據粗差主要位于條帶兩側,呈毛刺狀(圖7)。所有數據均利用Reson SeaBat T50-P多波束測深系統采集,表1為Reson SeaBat T50-P多波束技術參數指標。

圖7 實驗區域

表1 T50-P多波束測深系統主要技術參數指標
為了驗證本文算法有效性,分別利用CSF算法、人工剔除粗差與本文改進的CSF算法對這兩處多波束點云數據進行濾波。在測試中,CSF算法參數設定為St=true,GR=1,RI=1,n=500,dT=0.5,hcc=0.65。
圖8、圖9分別采用三種濾波方式對兩處多波束點云得出的結果,通過對比分析可得:

圖8 基于溫州數據多種濾波方法比對

圖9 基于青島數據多種濾波方法比對
(1)三種濾波方式均可對多波束點云數據粗差點剔除,可以看出所有方法對平坦地形的水下粗差剔除效果較好。
(2)CSF算法由于參數固定,雖然對相對平坦地形濾波效果較好,但是對于復雜的水下地形難以達到理想的效果,產生了一定程度的過度濾波,導致部分水下地形數據的丟失,并且當測線面積較大時,所需計算時間和計算機內存成倍增加。
(3)人工剔除粗差時會遺漏少量粗差點。本文算法則彌補了CSF算法的劣勢,針對不同水下地形進行自動參數調整,凸顯了水下地形濾波的優勢。通過對整體數據對比分析,驗證了本文算法不僅能有效剔除多波束點云數據粗差也能更好地保留較大起伏的水下地形特征。
為了更好地量化本文算法優勢,將這兩條測線分別采用人工剔除粗差、CSF算法與本文算法濾波,濾波前后的各項指標如表2所示。與人工方法粗差點剔除比NR是指的利用算法剔除粗差點個數與人工剔除粗差點個數的比,即NR=Nn/Na,其中Na為人工剔除粗差點個數,Nn為算法剔除的粗差點個數[10]。
由表2可以看出,采用本文算法和CSF算法得到的NR,對于相對平坦區域和起伏較大區域分別從1.03下降到1.01和從1.46下降到1.19。從表2中可以看出人工剔除粗差效果雖然很好,但是最小水深均要比兩種算法要小,經檢查發現有遺漏粗差點;CSF算法在平坦地形條件下濾波效果較好,但是在起伏較大地形條件下存在過度濾波問題,導致局部地形數據丟失;本文算法表現較好,平坦和起伏較大地形均有良好的濾波效果,且多波束數據完整。

表2 多種方法濾波性能比對
(1)人工剔除粗差效果較好,但是會遺漏粗差點,表現為最小水深不符合實際水下地形。
(2)CSF濾波算法對相對平坦的水下地形濾波有較好的效果,但是遇到有較大坡度的水下地形時,會出現過度濾波的現象。
(3)基于點云分割和自適應參數調整提出了改進的CSF算法,與CSF算法相比,克服了過度濾波現象,而且對于不同地區與人工方法粗差點剔除比分別從1.03下降到1.01和從1.46下降到1.19,保證了數據的完整性。與人工剔除粗差相比,避免了遺漏粗差點,人工干預很少。但是對于起伏更大地形,需進一步驗證。