張明希，許捍衛(wèi),于艷超，俞禮彬

（1. 河海大學 地球科學與工程學院，江蘇 南京 210098）

基于Python的多波束測深數(shù)據(jù)壓縮研究

張明希1，許捍衛(wèi)1,于艷超1，俞禮彬1

（1. 河海大學 地球科學與工程學院，江蘇 南京 210098）

在曲線數(shù)據(jù)壓縮的垂距限差法基礎(chǔ)上，引入總體最小二乘算法對多波束測深的ping條帶數(shù)據(jù)進行壓縮，利用Python工具庫實現(xiàn)總體最小二乘壓縮算法并集成到ArcToolbox中。與傳統(tǒng)道格拉斯-普克壓縮算法具有接近的壓縮比，但精度更高，能夠在整體上更好地表示原始數(shù)據(jù)。

多波束數(shù)據(jù);抽稀;總體最小二乘；Python

多波束測深系統(tǒng)是一種高效、高精度的水下地形測量系統(tǒng)，數(shù)據(jù)繁雜，即使單個ping也能產(chǎn)生幾兆數(shù)據(jù)，對多波束測深數(shù)據(jù)的壓縮處理非常重要。多波束測深系統(tǒng)沿著測線逐ping觀測，并以ping為順序存儲，可以從Ping數(shù)據(jù)的壓縮入手。關(guān)于數(shù)據(jù)壓縮問題，國內(nèi)外學者發(fā)展了許多算法。道格拉斯-普克算法（RDP）從整體上對曲線段進行考慮，通過保留特征點、刪除冗余點來達到壓縮目的[1]?？傮w最小二乘（TLS）方法是近年發(fā)展起來的能同時顧及觀測值誤差和模型系數(shù)矩陣誤差的數(shù)學方法[2]。與RDP不同的是,TLS算法沒有簡單地刪除距離誤差小于限差的點, 而是略微放松對端點的精度要求, 通過TLS作直線擬合,從而提高曲線壓縮的精度，在整體上最優(yōu)逼近原始曲線。盧銀宏等[3]綜合利用RDP算法和TLS算法來壓縮多波束測深數(shù)據(jù)，采用Matlab仿真二維ping數(shù)據(jù)進行壓縮處理，忽略了當曲線的曲率很大時，兩條擬合直線的交叉點可能會遠離原始曲線的情況。本文基于ArcPy站點包和開源Python工具庫，以Python作為腳本語言，開發(fā)了多波束測深三維ping數(shù)據(jù)點抽稀算法，具有易用性和重用性特點。

1 Python空間數(shù)據(jù)處理特點

Python語言作為穩(wěn)定的動態(tài)腳本語言與ArcGIS高度集成，其靈活性和重用性可以方便實現(xiàn)地理空間數(shù)據(jù)的處理流程[4]。目前，ArcGIS10．0已集成基于Python的ArcPy站點包，可理解為ArcGIS的Python API。ArcPy提供了大量類和函數(shù)，可以輕松執(zhí)行ArcGIS工具箱中的工具，并創(chuàng)建原生對象，如幾何、柵格、空間參考等[5]。

另一方面，Python具有高級的內(nèi)建數(shù)據(jù)結(jié)構(gòu)，優(yōu)秀的科學計算擴展庫SciPy、NumPy（包括科學計算、統(tǒng)計分析、可視化等模塊），高效的空間分析開源庫GDAL、Shapely等。其中SciPy包含的Python模塊有最優(yōu)化、線性代數(shù)、積分、插值、特殊函數(shù)、快速傅里葉變換、信號處理和圖像處理、常微分方程等函數(shù)[6]。Shapely庫是針對平面幾何對象的空間分析工具。基于點集理論，其將幾何對象分為點、線、面3種類型，每種數(shù)據(jù)類型又有3個屬性：內(nèi)部集、邊界集、外部集。Shapely以Python風格的代碼形式訪問GEOS庫，將PostGIS的一般操作用Python語言加以實現(xiàn)，適合處理一些非常規(guī)的地理空間數(shù)據(jù)問題。

2 基于Python的總體最小二乘點抽稀

2.1 總體最小二乘直線擬合模型

傳統(tǒng)最小二乘法擬合直線時，一般令直線方程為：

其中，x，y是測點坐標；a，b為待估參數(shù)。最小二乘法認為測量誤差都存在于y中，而x沒有誤差或在進行平差時不予考慮。則誤差方程可以表示為：

按最小二乘準則，得出其最小二乘解為：

一般認為，觀測值y的獲取方法相同，因此P為單位陣。

最小二乘只認為觀測值y含有測量誤差，事實上x、y均是采用測量方法獲取的，按照總體最小二乘思想，應(yīng)該對觀測值x和y都進行改正?？紤]到觀測值x的誤差，則條件方程為：

寫成矩陣形式為：

其中， Vy為y的改正數(shù)；EA為系數(shù)陣A的改正矩陣。由于系數(shù)陣里含有固定列，運用混合總體最小二乘的迭代解法[7],求出參數(shù)a與b。

總體最小二乘直線擬合算法需要在SciPy的線性代數(shù)庫linalg基礎(chǔ)上，利用NumPy豐富的數(shù)組處理函數(shù)，按照總體最小二乘的迭代解法，重新編寫Python腳本計算。

2.2 基于總體最小二乘的多波束數(shù)據(jù)點抽稀算法

RDP算法在給定限差情況下的逼近精度不高。如圖1，當采用RDP算法對所示曲線數(shù)據(jù)壓縮時，會得到線ABC。在保留點AB之間，被刪除的點全部位于保留直線的一側(cè)；在保留點BC之間，被刪除點也全部位于保留直線的一側(cè)。顯然，用ABC表示原始數(shù)據(jù)具有明顯的偏差，且均方誤差較大。如圖1a，采用基于TLS的分段直線擬合算法，對AB和BC之間的原始點分別作直線擬合，得到兩條直線A′B′、B′C′，它們相交于B′。此時，利用A′-B′-C′代替A-B-C，可以提高對A、B和C點的壓縮精度。但在個別情況下，如圖1b，曲線的曲率很大時，兩條擬合直線的交叉點可能遠離原始曲線，明顯偏離了原曲線的形狀。

圖1 直線擬合和交叉點遠離原始曲線的處理

利用多波束測深的ping點數(shù)據(jù)，每個ping的相鄰數(shù)據(jù)點由相鄰波束腳印的中心確定，每個ping描述了某個橫向的地形剖面，其x、y、z坐標可映射到二維平面?；诘栏窭?普克的總體最小二乘三維點抽稀算法的基本步驟如下：

Step1：將多波束測深ping點三維坐標（x，y，z）變換為二維曲線坐標（dx，dy）。計算ping條帶上某點O到空間直線AB（ping條帶上的首尾點）的垂足坐標N（xn，yn，zn），然后求得N到點A和C的距離dx、dy，其中A（x1，y1，z1）、B（x2，y2，z2）為：

Step2：利用RDP算法對曲線坐標作壓縮抽稀，選取抽稀后的特征點（Mx，My），并標記特征點所對應(yīng)的索引值index。

Step3：搜索index和index+1中所包含的二維曲線點坐標對。考慮到TLS至少需要3個點坐標，當點坐標對數(shù)大于2時，利用TLS求出系數(shù)a、b。遍歷所有index，求出直線擬合的交點坐標集（Hi，Hj）（i，j=1，2，…），記錄其序號。

Step4：計算直線擬合的交點坐標（Hi，Hj）和特征點（Mx，My）間距離d。如果d大于限差ε，將直線擬合交點坐標替換成特征點，再運用基于Shapely庫函數(shù)中的曲線intersection方法求出擬合直線段和原始曲線的交點。若某段出現(xiàn)多個交點，排序得出距特征點最近的交點，并補充到直線擬合交點坐標集中；否則只保留直線擬合交點坐標。

Step5：用步驟4得到的二維點集，解算空間平面坐標方程，反算得到抽稀后的三維ping點數(shù)據(jù)。

算法用Python腳本實現(xiàn)后，添加到ArcToolBox中，如圖2。

圖2 總體最小二乘算法點抽稀工具界面

2.3 實 驗

選取某江水下地形多波束測深數(shù)據(jù)ping條帶數(shù)據(jù)做測試，驗證算法的有效性和正確性。該ping帶含有62個測深點。用ArcToolbox中的自定義點抽稀工具腳本進行實例計算。

根據(jù)國際海道測量規(guī)范s44 特等測量的精度要求,抽稀閾值設(shè)為0．26 m[8]。圖3顯示了垂距限差取0．26 m時，分別采用道格拉斯-普克算法和分段總體最小二乘的壓縮效果。其中藍色實線代表原始ping數(shù)據(jù)，黑色點表示RDP算法抽稀效果，紅色點表示TLS點抽稀效果。

圖3 點抽稀效果圖

為評定算法的抽稀精度，計算RDP算法和TLS算法的均方誤差，公式為其中di表示原始曲線上第i個點到壓縮后曲線的垂直距離，n表示總點數(shù)。垂距限差分別取0．26、0．5、0．1進行測試。TLS和RDP的壓縮比近似相等，均方根誤差為當垂距限差較大時,TLS點抽稀算法比RDP算法精度提高更明顯，能更好地逼近原始數(shù)據(jù)。

信息科學版,2013,38(7):850-856

[3] 盧銀宏,岳東杰，宋飛鳳．基于總體最小二乘的Douglas-Peucker算法在多波束測深數(shù)據(jù)抽稀中的應(yīng)用[J]．水利與建筑工程學報,2012,10(2):4-5

[4] 彭海波,向洪普．基于Python的空間數(shù)據(jù)批量處理方法[J]．測繪與空間地理息,2011,34(4):81-82

[5] Dobesova Z．Programming Language Python for Data Processing[C]．International Conference on Electrical and Control Engineering (ICECE), 2011

[6] 張若愚．Python科學計算[M]． 北京：清華大學出版社，2012

[7] 魯鐵定,周世．總體最小二乘的迭代解法[J]．武漢大學學報：信息科學版,2010,35(11):1 351-1 354

[8] 夏偉,黃謨濤,劉雁春,等．Douglas-Peucker算法在多波束測深數(shù)據(jù)抽稀中的應(yīng)用[J]．測繪科學,2009,34(3):159-160

[1] 楊云,孫群,朱長青．曲線數(shù)據(jù)壓縮的總體最小二乘算法[J]．西安電子科技大學學報：自然科學版,2008,35(5):946-950

[2] 王樂洋,許才軍． 總體最小二乘研究進展[J]．武漢大學學報：

P229.3

B

1672-4623（2014）04-0117-03

10.11709/j.issn.1672-4623.2014.04.041

張明希，碩士，主要研究方向為地理信息系統(tǒng)開發(fā)與應(yīng)用。

2013-10-10。

項目來源：國家自然科學基金資助項目（41101374）。