基于濾波算法的單波束測深數據預處理程序開發

黃瑞貞 星睿

收稿日期：2023-05-16

基金項目：廣東省海洋經濟發展專項（粵自然資合〔2022〕21號）

DOI：10.19850/j.cnki.2096-4706.2023.22.020

摘? 要：在測深過程中單波束因受儀器、海水運動等因素的影響，常出現水深異常值。為解決異常值的問題，利用C#編程語言編寫了中值濾波法和加權濾波法，實現了較好的濾波效果。同時設計開發了單波束測深數據的預處理軟件，其主要功能包括濾波處理、水位改正、精度對比等，具有良好的交互界面，操作便捷，提高了單波束測深數據預處理的效率。

關鍵詞：濾波算法；單波束測深；C#語言；軟件設計

中圖分類號：TP391；P229? ? ? 文獻標識碼：A? ? ? 文章編號：2096-4706（2023）22-0092-05

Development of Single Beam Sounding Data Preprocessing Program Based on Filtering Algorithm

HUANG Ruizhen， XING Rui

（Surveying and Mapping Institute of Lands and Resource Department of Guangdong Province， Guangzhou? 510663， China）

Abstract： In the process of single beam sounding， due to the influence of instruments and sea water movement， the single beam bathymetric data often contains outliers. To solve the problem of outliers， the median filtering method and weighted filtering method were written using C # programming language， achieving good filtering results. And We have designed and developed a preprocessing software for single beam sounding data， its main functions includes filtering processing， water level correction， accuracy comparison， etc. it has a good interactive interface， convenient operation， and improves the efficiency of single beam sounding data preprocessing.

Keywords： filtering algorithm; single beam sounding; C# language; software design

0? 引? 言

水深測量作為海底地形測量的基本手段，是海洋工程建設的基礎工作。當前，聲學手段依舊是獲取海洋地形地貌數據工作中可靠性最高、實用性最好的探測技術。而單波束測深儀作為一種聲學水深測量儀器，其操作簡便、數據精準、性能穩定，在我國水下地形地貌勘察領域應用非常廣泛。但在勘測的過程中，受涌浪、海流等因素的影響，測量船的船姿參數（如航偏角、縱搖、橫搖和動態吃水等）會發生瞬時變化，由此帶來測量誤差，并導致成果失真或地形畸變，所測水深點與實際位置偏差較大。因此，在水深測量中，利用濾波方法消除誤差、剔除虛假信息，對提高水深測量的效率及精度尤為重要。

傳統的粗差探測方法是依靠專業人員手工判別，隨著技術的快速發展，諸多學者對測深數據濾波算法進行了大量的研究，相繼提出了多種通過計算機自動判別異常數據的方法[1]。例如，2007年，李明叁等設計了海洋水深異常數據檢測的順序統計濾波方案，大幅提高了粗差檢測的效率[2]；2009年，孫嵐等比較分析了部分多波束測深數據異常值檢測算法的優缺點[3]；2018年，王君婷等以單波束無驗潮實測水深數據為例，對測深數據中的粗差濾波方法進行對比分析，結果表明于M估計的選權迭代的效果最好[4]；2021年，牛沖以單波束測深數據為例，對中值濾波和加權濾波在不同大小濾波窗口的測深數據處理效果進行了對比分析，給出了不同環境下的濾波方案[5]。有關濾波算法的程序實現，2012年由佘東等設計的數字濾波程序能有效實現濾波[6]；測深數據預處理方面有2016年盧凱樂等基于C#實現多波束測深數據預處理軟件的開發，該軟件具有較高的測量精度和可靠性[7]。但相關文獻有關單波束測深數據濾波處理程序開發的內容極少，而已有測深數據處理軟件中有關濾波處理的功能較為單一，用戶交互性一般，大多需要人工檢查，降低了工作效率。本文基于濾波算法實現了單波束測深數據的濾波處理程序，并開發了測深數據預處理軟件。

1? 測深數據濾波算法

1.1? 中值濾波法

中值濾波法是數字濾波算法中較為常見的一種方法，其基本原理是單波束測深點鄰域范圍內多個測深點排序后的中值代替異常點的測深值，以消除偶然因素引起的波動或脈沖干擾等其他尖銳的噪聲點。設有n個測深點的測深值為h1，h2，…，hn，取窗口長度2M + 1（2M + 1為奇數，且2M + 1≤n）。為檢測測深點i（i ∈ [1，n]）是否為測深異常值，從輸入序列中相繼抽出2M + 1個測深值Hi-M，…，Hi，…，Hi+ M，再將這2M + 1個點按水深值由小到大的順序進行排序，取排序后的中值作為濾波輸出，用計算式表示為：

hi = med{Hk，k ∈ [i - M，i + M]∩[1，n]}，i ∈ [1，n]（1）

其中，2M + 1表示中值濾波的窗口大小，n表示測深點的數量，med{}表示窗口內水深序列取中值的計算。

1.2? 加權濾波法

加權濾波法是對測深點鄰域范圍內多個連續的測深點分別乘上不同的加權系數后再求和，代替該點的測深值。設n個測深點的測深值為h1，h2，…，hn，為了檢測測深點i（i ∈ [1，n]）是否為測深異常值，取濾波窗口的大為2M + 1（2M + 1為奇數，且2M + 1≤n），即從n個點序列中相繼取出2M + 1個測深值H = [Hi-M，…，Hi，…，Hi+M]。所應用的計算式如下：

（2）

對于不在觀察序列中的點，水深值補為零值，其中2M + 1表示濾波窗口大小，A表示水深點對應的權值，約束ΣA = 1。

2? 軟件設計及實現

軟件的設計思路依據單波束測深數據的預處理流程，如圖1所示。在讀取單波束測深數據后，結合參數配置進行濾波處理，處理后可以導出測深數據和異常點數據；可以剔除異常點數據或通過插值方法對其進行改正；讀取驗潮站潮位數據，對濾波及插值后測深數據進行水位改正；測深數據全部預處理后可進行精度統計，獲得主測線和檢查線的精度分析。

2.1? 軟件功能設計

軟件的主要功能包括工程、數據分析、數據導出、參數配置、數據導出、數據分析結果、幫助模塊。軟件功能設計如圖2所示。

工程模塊主要進行文件的讀取，包括讀取dat、EDT等格式的單波束測深數據及驗潮站實測水位數據等，此外還具備數據刷新、圖像刷新等功能。數據分析模塊提供測深數據的濾波處理、異常值的擬合計算、水位改正、精度對比等功能。其中濾波處理方法包括中值濾波法和加權濾波法，異常值的擬合計算包括最小二乘擬合法和反距離加權法，水位改正包括單站改正、雙站改正和三站改正。數據導出包括導出異常點數據、擬合結果數據，以及預處理后測深數據等。參數配置主要為數據處理時的窗口大小提供可選項。數據顯示模塊提供測深數據及水位數據的顯示及對比、濾波處理異常點的顯示查詢、數據處理前后對比等功能。數據分析結果主要顯示數據分析完成后的數量結果和精度統計并繪制圖像顯示。

2.2? 界面設計功能實現

本軟件以.NET為平臺，利用Microsoft公司的Visual Studio 2019編程軟件，采用C#語言進行開發，根據程序功能的需要，添加了相應的組件（菜單和工具欄組件、對話框組件、chart控件、其他公共控件），用以實現程序的文件打開、保存、數據可視化等功能。程序設計主界面如圖3所示。

2.3? 核心代碼

2.3.1? 中值濾波法

讀取單波束測深數據，設置濾波的窗口大小，窗口測深點采樣數量通常為奇數，防止偶數采樣在排序和取值時中值不唯一。對測深數據連續采樣后，按照從小到大的順序對窗口內水深值進行排序，取中間值作為本次采樣值，依據2σ或3σ原則對粗差進行檢驗和標定。

具體C#語言程序實現如下：

int mn = （（int）numericUpDownLBCK.Value-1）/2;//濾波窗口為總共采樣測深點，mn值為取值前后各采樣點數量

int fc = （int）numericfc.Value;? //幾倍均方差，取2或者3

for（int i=mn;i

{

List onewindow = new List（）;? //一個窗體的全部數據

double sum = 0;? //為求平均值而求和

for（int k=i-mn;k

{

sum += depdata[k].z;

onewindow.Add（depdata[k].z）;

}

double aveh = sum / （2 * mn + 1）;? //求水深值的平均

sum = 0;

for （int k = i - mn; k < i + mn + 1; k++）

{

sum += Math.Pow（（depdata[k].z - aveh）， 2）;

}

double sig =Math.Sqrt（sum/ （2 * mn + 1））;//求中誤差

onewindow.Sort（）; //對窗口內水深值進行排序

double medvalue = Math.Abs（depdata[i].z - onewindow[mn]）; //與中值的差值

if（medvalue > fc * sig）

{

depdata[i].limit = true;? //超限數據的標定

}

}

2.3.2? 加權濾波法

加權濾波法同樣要在讀取單波束測深數據后，設置濾波的窗口大小和2σ或3σ作為異常值的判別依據。連續采樣后，對窗口內測深值分別乘以不同的加權系數再求累加和，各權重總合應滿足等于1的約束條件，取加權運算后的和為本次有效采樣值。

具體C#語言程序實現如下：

for （int i = mn; i < n - mn; i++）

{

double[] onewindow = new double[mn*2+1];//一個窗體的全部數據

double sum = 0;? //為求平均值而求和

int m = 0;

for （int k = i - mn; k < i + mn + 1; k++）

{

sum += depdata[k].z;

onewindow[m] = depdata[k].z;

m++;

}

double aveh = sum / （2 * mn + 1）;? //求水深值的平均

sum = 0;

double sumw = 0;? //權重的和

for （int k = i - mn; k < i + mn + 1; k++）

{

sum += Math.Pow（（depdata[k].z - aveh）， 2）;? //差值的平方和

double ww = 1 / （1 + Math.Abs（depdata[k].z - aveh））;? //權重

sumw += ww;

}

double[] AA = new double[mn * 2 + 1]; //權重的矩陣

double AHT = 0;? //A*H的轉置，即水深和權重相乘的總和

for （int k=0;k

{

AA[k] = （1 / （1 + Math.Abs（onewindow[k] - aveh））） / sumw; //每個水深的權重

AHT += AA[k] * onewindow[k];? //水深和權重相乘

}

double sig = Math.Sqrt（sum / （2 * mn + 1））;? //求中誤差

double jq = Math.Abs（depdata[i].z - AHT）; //與加權水深的差值

if （jq > fc * sig）

{

depdata[i].limit = true;? //超限數據的標定

}

3? 軟件測試

以廣東省汕頭市附近海域為試驗區，采用有人船和無人船相結合的方式，運用單波束測深系統采集水下地形數據，定位方式為PPK動態后處理差分模式，獲得測深數據后采用本軟件進行數據預處理測試。

3.1? 測深數據濾波處理測試

為驗證濾波算法的實際應用效果，選取2022年5月5日近岸測深數據進行濾波處理，當日海況為一級風，考慮風浪影響，在航道內布設垂直于岸線的東西走向測深線，測線間隔為100 m，最大水深為22 m，最小水深為2 m，合計23 744個測深點。將單波束測深文件導入軟件，在參數配置框內設置濾波窗口大小、異常值的判別原則即均方差倍數，在數據分析中選擇數據濾波，選擇中值濾波法或加權濾波法進行粗差探測。經過多次試驗，針對中值濾波法選擇濾波窗口大小為5，取2σ為異常值的判別依據，結果顯示異常點數量為462個；針對加權濾波法選擇濾波窗口大小為9，同樣取2σ為異常值的判別依據，結果顯示異常點數量為325個。濾波數據完成后在數據分析結果欄顯示濾波結果，點擊輸出結果，將導出異常點文件和剔除異常點的測深文件，點擊顯示結果，將在右邊圖像框中顯示整條測深數據的曲線和異常點的分布，如圖4所示。在實際生產中，應結合不同的測量條件及干擾情況選用不同的參數和濾波算法，實現最佳濾波效果。

3.2? 測深數據預處理過程

測深數據經過濾波處理后，可對濾波結果進行處理，剔除粗差或選擇異常值的擬合內插，在參數配置中設置內插窗口的大小，在數據分析中選擇合適的異常值內插方法，可選最小二乘擬合法或反距離加權法。

水位改正計算方法復雜煩瑣，耗時費力，為減少內業工作量，為本文編制相應的水位改正程序。根據測深數據的位置及驗潮站分布和控制范圍，判斷所選用的水位改正方法：若整個測區在某站控制范圍內，則采用單站改正；若測區在兩站之間且兩站水位變化相似，則采用雙站改正；若測區范圍大、潮汐變化復雜，則采用三站改正[8]。采用單站改正時，讀取潮位數據并進行抽稀或擬合內插，使時間序列與測深數據時間一致，即可進行水位改正。雙站改正和三站改正均采用時差法，利用相鄰站的同步水位數據求解站間潮時差，通過空間內插潮時差將水位內插至任一測深點的測深時刻。本軟件操作簡便，在測深數據經過粗差處理后，導入驗潮站文件，在數據分析中選擇適宜的水位改正方法即可完成。

對全部測深數據完成預處理后，進行精度統計，包括測深數據擬合前后的對比以及主檢測線的精度對比。測深數據擬合前后對比統計其最近點差值的中誤差、最大值、最小值、平均值及各差值區間的分布比例；主檢測線的精度對比統計其檢查線10 m范圍內與主測線最近點間的差值中誤差、標準偏差、最大值、最小值、超限率、數據量、超限數據及各差值區間的分布比例等。統計完成后在數據分析結果欄中顯示主要分析結果，點擊輸出結果將輸出詳細對比數據及結果。

4? 結? 論

本軟件基于C#語言進行編寫，實現了單波束測深數據的預處理。對于測深數據的粗差探測，采用了中值濾波法和加權濾波法兩種處理算法，可查看并導出異常值；對異常點的處理，采用了最小二乘擬合法或距離加權法進行內插擬合，同樣可繪制曲線查看對比擬合效果。本軟件還可對濾波處理后的測深數據繼續進行水位改正和精度對比，優化自動化處理流程，操作簡單便捷，減少了煩瑣的數據處理工作，降低了作業人員的干預誤差，可有效提升數據質量，提高測深數據內業數據整編的工作效率。

本文總體完成了具有良好交互界面的單波束測深數據預處理應用軟件的設計和實現。但在實際的應用中，還有許多其他的濾波算法可用于單波束測深數據的預處理還有許多其他的濾波算法，及測深數據的處理。今后在現有算法研究的基礎上做進一步的研究，結合實際工作環境及需求狀況，對程序進行持續不斷的完善和改進，以期發揮其良好的應用價值。

參考文獻：

[1] 黃謨濤，翟國君，王瑞，等.海洋測量異常數據的檢測（英文）[J].測繪學報，1999（3）：269-276.

[2] 李明叁，張杰，嚴懷志，等.序統計濾波估計檢測海洋測深異常數據 [J].海洋測繪，2007，114（1）：7-11.

[3] 孫嵐，王海棟，余成道，等.多波束測深數據異常值檢測算法比較[J].海洋測繪，2009，29（5）：57-60.

[4] 佘東.單片機數字濾波程序設計與實現[J].四川兵工學報，2012，33（9）：78-80.

[5] 王君婷，王振杰.單波束無驗潮測深數據濾波方法[J].工程勘察，2018，46（3）：48-52+74.

[6] 牛沖.不同窗口下中值濾波與加權濾波方法效果分析[J].測繪科學，2021，46（8）：21-27+40.

[7] 盧凱樂，王勝平.基于C#的多波束測深數據預處理軟件設計與實現[J].江西科學，2016，34（1）：56-59.

[8] 劉云波，李柏.時差法在測深數據水位改正中的應用[J].測繪地理信息，2014，39（6）：14-16.

作者簡介：黃瑞貞（1986—），女，漢族，廣東揭陽人，工程師，本科，研究方向：海洋測繪；通訊作者：星睿（1997—），女，土族，青海海東人，助理工程師，本科，研究方向：測繪工作及軟件開發。