王懷計，?？』郏?，王瑞富

（山東科技大學 測繪與空間信息學院，山東 青島 266590）

海底地形地貌的準確性表達對船只航行的安全性有著重要的影響，為了獲得準確的海底地形地貌，為海洋學各種研究提供基礎，需要對海底水深數據進行實測[1-6]，而實測的水深數據往往是離散點的形式，若某個區(qū)域的多波束測深數據缺少，且現有其他來源的水深數據的密度不能有效地反映海底地形地貌，通常需要通過空間內插方法增補水深點來得到足夠多的海底地形地貌表面數據[7-8]。常用的空間插值方法有經典的反距離加權（inverse distance weighted，IDW）、自然鄰域、樣條函數法等，其中，IDW 算法以其原理簡單而被廣泛應用。

近年來，國內外學者對IDW 算法進行了許多研究，在算法的最優(yōu)插值參數選擇方面，張錦明等比較了該算法中權指數[9-10]、搜索點數[11]、搜索方向[12-13]等參數對插值誤差的影響，得出了插值算法中的最優(yōu)插值參數；在算法的應用方面，蔣偉達等利用IDW 對1958—2014 年埕島7 個不同時期的水深點進行數據插值，分析各個時期的地形特征，得出了埕島水下地形演變規(guī)律[14]；在算法的改進方面，何立恒等認為已知樣本點在內插點的全圓方位上分布不均勻會對插值結果產生影響，提出了顧及夾角的改進IDW 算法，該方法能有效解決樣本點分布不均對插值結果的影響[15]；樊子德等認為在不同的地形中人為調節(jié)權指數較為煩瑣，提出了一種改進的IDW 算法，該方法能自適應地調整權指數，提高了算法的自適應性[16]。

通過上述對IDW 的分析可知，無論是在算法的插值參數選擇方面還是在算法的應用和改進方面，IDW 插值的計算量問題都在于：若要求得一個待插值水深點的水深值，需要遍歷所有已知水深點并進行距離的計算。在海量數據下，對數據點進行有效地索引可以減少計算量。另外，在算法改進方面，大多IDW 均假設空間過程具有平穩(wěn)性，導致復雜海底地形的待插值水深點插值精度較低。因此，針對傳統(tǒng)的IDW 方法在復雜海底地形適應性差的局限性，本文提出了一種顧及特征水深點距離重分配的反距離加權插值算法（distance-redistribution of characteristic depth points IDW，DCIDW），該算法首先對待插值水深點建立二階相鄰點，并在此基礎上建立空間索引，利用二階相鄰點的個數自動調節(jié)搜索點數；然后對待插值水深點周圍的樣本水深點進行特征提取，特征水深點能反映樣本水深點之間的相關性；最后根據距離分配因子進行距離重分配。該算法可在顧及樣本水深點內部之間特征性的情況下，提高待插值水深點的精度。

1 特征水深點距離重分配DCIDW算法

經典IDW（classical IDW，CIDW） 算法假設待插值水深點的水深值與各個樣本水深點之間的距離相關，通常采用距離平方的反比來定權，而忽略了海底地形的不平穩(wěn)變化。在復雜海底地形中，CIDW 插值方法將所有的樣本水深點視為同類別的點，不區(qū)分樣本水深點之間的差異，此時相鄰的水深點之間的水深值相差較大的情況如圖1 所示。

圖1 復雜海底地形傳統(tǒng)水深插值誤差示意圖

根據反距離加權插值算法原理，待插值水深點A 的水深值由樣本水深點B、C、D、E 以及其他周圍樣本水深點的水深值共同加權決定，因此A′的水深值應在最深點C 和最淺點D 之間。樣本水深點B、C 與A 點相距較近，根據權函數要求，應賦給點B、C 較大的權，但是實際上，樣本水深點D、E 對待插值水深點A′的影響更為重要，由于距離較遠，所賦的權較小，導致待插值水深點與真實水深點之間存在較大的誤差。圖1 中的A′為忽略空間特征點D 和E 插值導致的結果，這不僅導致水深點不能反映真實的海底地形，而且會給航海安全帶來重大隱患。因此，探索樣本水深點中特征水深點的提取，進一步進行距離重分配，提高待插值水深點的真實性是非常必要的。

特征水深點的距離重分配算法主要包括特征水深點提取、特征點距離重分配、插值3 個步驟，詳細流程如圖2 所示。

圖2 特征水深點距離重分配算法流程圖

本文提出的特征水深點距離重分配算法，在樣本水深點與待插值水深點之間距離之和不變的約束條件下，首先對實驗區(qū)域內的樣本水深點進行特征點提取，來反映樣本水深點之間的空間相關性；然后提出距離重分配指標，對符合指標的待插值水深點進行距離重分配；最后構建一個顧及樣本水深點之間空間屬性的插值模型，提高特征點的權重，從而提高插值的精度。

1.1 特征水深點的提取

海底水深點的離散性造成了其空間關系分析上的復雜性，Delaunay 三角網根據離散點之間的相互關系，對離散點水深值的深淺進行描述。特征水深線是對水深點之間坡向關系的描述，選取局部區(qū)域內的相對淺點來反映相鄰水深點之間的地形變化情況，特征水深點的數量越多，則地形震蕩趨勢越明顯。基于上述原理，本文首先分析離散點之間的空間關系，然后構建特征水深線，并在此基礎上，完成特征水深點的選取。

本文在對離散的水深點構建Delaunay 三角網的基礎上，對離散的水深點構建特征水深線，根據特征水深線上各點之間的水深關系，提取出特征水深線上的特征水深點。特征水深線的提取是根據三角網之間的拓撲關系，遞歸追蹤相鄰點中的最淺點，依次連接各淺點，形成的過渡虛擬連接線。其能為特征水深點的提取提供判斷依據，因此特征水深線是一種用來定義水深點之間關系的概念線。在特征水深線的基礎上，對特征水深線上的水深點進行分析，比較線上各點的左右相鄰點的深度差，進而選取特征水深線上的相對淺點作為特征淺點。

本文算法對特征水深點進行距離重分配的前提是樣本水深點中的特征水深點的提取，該過程的主要任務是將離散的樣本水深點結構化，判斷樣本水深點之間的空間關系，獲取特征水深點[17-18]。離散的樣本水深點不能反映各樣本水深點之間的相關性，對離散的樣本水深點構建Delaunay 三角網是建立點群空間相關性常用的方法之一，它是一系列相鄰、不重疊的三角形構成的三角網，根據三角網之間的空間特性和拓撲關系，能有效地表達海底地形起伏形態(tài)[19]。

根據三角網之間的拓撲關系，對所有已知樣本點遞歸追蹤一階相鄰點中的最淺點，生成多條特征線，記錄每一條特征線上的特征樣本點。特征水深線的提取如圖3 所示。

圖3 水深特征線生成過程

特征水深點的提取過程如下：

（1） 單條特征線的提取。對樣本水深點構建Delaunay 三角網，將三角形存入數組TIN；以全局水深最淺點Pmin為起始點，將Pmin存入數組PL；搜索起始點的相鄰點，找到相鄰點中與起始點深度差最小的點PN，追加存入PL。PL內的點不會再被搜索到；以PN為起始點，重復上述操作，直到找不到相鄰點為止。按序連接PL中的點，得到一條特征線；

（2）所有特征線的提取。找到剩余未被記錄的水深點中的最淺點，作為起始點，重復步驟（1）中的操作，直到所有的點都被搜索到為止，最終得到若干條特征線。

（3）特征水深點的提取。特征水深點是特征水深線上的特征淺點，能反映局部區(qū)域內的水深差值情況，特征水深點的提取過程如圖4 所示。

圖4 基于特征水深線的特征水深點提取過程

具體步驟為：①設特征線上水深點數為n，記為P1，P2，…，Pn，點的水深值記為H1，H2，…，Hn；②若Hi＜Hi-1且Hi＜Hi+1（i=2，3，…，n），則Pi是水深淺點，存入數組Pc；③若連續(xù)多個水深點的深度值相同，則與連續(xù)起點的前一點的水深值和連續(xù)終點的后一點的水深值進行比較，選取Pi和Pi+1中的一點作為特征點；④重復上述步驟，直至特征水深點提取完成。

1.2 基于特征點的距離重分配

CIDW 插值算法不考慮周圍的樣本水深點的特征性，只簡單地考慮距離倒數的平方對權函數的影響。本文在提取特征點之后，在所有樣本水深點到待插值水深點的距離和不變的約束下，量化比較待插值點的二階相鄰點中特征點的個數占比與特征點距離之和與總距離占比，實現對復雜的海底地形特征點的距離重分配。

（1）二階相鄰點的選擇。待插值點的二階相鄰點選擇的實質為空間鄰近分析問題，Delaunay 三角網作為空間鄰近分析工具之一，廣泛應用于空間數據聚類和空間鄰近關系探測[20]，故本文應用Delaunay 獲取待插值點周圍的一階相鄰點[21]，在一階相鄰點的基礎上，獲取待插值點的二階相鄰點，其中包含特征點和非特征點，樣本點均勻地分布在待插值點的周圍。

（2）特征點距離重分配。根據待插值點二階相鄰點中特征點的個數Nfeature和每個特征點到待插值點的距離dfeature_i，計算得到待插值點的二階相鄰點中特征點個數占比Ratio_N 和特征點距離占比Ratio_D，從而可以獲取每個待插值點是否需要距離重分配判斷指標，計算公式如下：

從上式提出特征點距離重分配指標：若Ratio_D ＜Ratio_N 表明特征樣本點距離待插值點較近，特征點對待插值點的權重受地形影響較小，因此無需對特征點的距離進行重分配；若Ratio_D ≥Ratio_N 表明特征點距離待插值點較遠，特征點對待插值點的權重受地形影響較大，按照閾值的比例對特征點的距離重新分配，分配公式如下：

特征樣本點的距離重分配算法具體步驟如下：①獲取待插值水深點的一階相鄰點，在一階相鄰點的基礎上獲取待插水深點的二階相鄰點n；②獲取待插值水深點與二階相鄰點的距離di；③比較閾值Ratio_N 與Ratio_D 大小，判斷是否需要對特征點的距離重新分配；④重復上述步驟，直至待插值水深點距離重分配完成。

樣本點的距離和個數是影響待插值點精度的兩個重要因子，本文引入了特征點距離占比和特征點個數占比兩個指標，若特征點距離占比大于特征點個數占比，表示其相較于其他樣本水深點與待插值水深點之間的距離較遠，因此用特征點個數占比代替特征點重分配之后的距離占比，能使特征水深點與待插值水深點之間的距離變近，在二階相鄰點總距離不變的約束下，其他樣本水深點重分配之后的距離和為總距離與特征水深點重分配距離的差值。

1.3 插值

由于待插值水深點的二階相鄰點中可能存在特征樣本點距離待插值水深點較遠，但對待插值水深點影響較大的情況，特別是距離待插值點的水深值與較近的其他水深點的水深值差別較大會使CIDW插值算法得到的插值結果與真值相差較大。基于此，本文在CIDW 插值算法的基礎上，對每個特征點的距離重新分配，計算公式如下：

式中，dfeature_i為每個特征水深點與待插值水深點之間的原始距離，d′feature_i為每個特征水深點與待插值水深點之間的重分配距離，d′others_i為每個其他水深點與待插值水深點之間的重分配距離。

基于上述公式，本文提出了一種顧及待插值水深點的二階相鄰樣本點中特征樣本點的IDW 插值模型：

2 實驗對比分析

2.1 特征水深點的選取

本文實驗數據來源于夏威夷茂伊島附近的機載激光雷達測深數據，特征水深點的提取實驗過程如下。

（1）本實驗選取區(qū)域內部分水深點作為實驗數據，包括待插值水深點以及其周圍的已知水深點，用于構建Delaunay 三角網，構建結果如圖5 所示。

圖5 三角網構建

（2）提取特征水深點。根據水深特征線構建方法，對整體水深點構建水深特征線，構建特征水深線的目的是為了分析已知樣本水深點之間的關系，提取樣本水深點中的特征水深點。所選取的部分水深點構建結果如圖6 所示。

圖6 特征水深線和特征水深點提取結果

2.2 基于特征水深點的樣本水深點的距離重分配

（1）在已提取的特征水深點基礎上，首先對待插值水深點構建二階相鄰線，得到待插值水深點的二階相鄰點，結果如圖7 所示。

圖7 待插值水深點的二階相鄰點

（2）按照待插值水深點的二階相鄰點中特征水深點與所有二階樣本水深點的個數比值和特征水深點與所有二階樣本點的距離比值進行分析，待插值水深點與樣本水深點之間的距離需要進行重分配，待插值水深點與樣本水深點之間的距離在總距離之和不變的約束下，對各個樣本水深點的距離進行重新分配，分配過程示意圖如圖8 所示。

圖8 樣本水深點距離重分配示意圖

（3）對待插值水深點與樣本水深點之間的距離進行重分配之后，按照式（7）對待插值水深點進行水深計算，待插值水深點與周圍二階相鄰點的實際地形如圖9 所示。

圖9 待插值水深點實際水下地形圖

圖中待插值水深點的實際水深值為13.9 m，藍色點為待插值水深點的二階相鄰點中其他水深點，與待插值水深點之間的平面距離分別為15.6 m、4.9 m、5.9 m、18.0 m、21.5 m、17.5 m、25.6 m、23.3 m、22.6 m、22.8 m，其他樣本點與待插值水深點之間重分配的平面距離分別為16.9 m、5.3 m、6.4 m、19.5 m、23.3 m、19 m、27.8 m、25.3 m、24.5 m、24.7 m，為更好地表現出分配之后的距離變化，其他水深點的距離重分配結果如表1 所示。

表1 其他水深點的距離重分配表

紅色點為經過特征提取之后的特征水深點，與待插值水深點之間的平面距離分別為13.5 m、30.0 m、26.3 m、22.0 m、15.6 m、23.4 m；根據本文提取的樣本水深點距離重分配方法，特征樣本點與待插值水深點之間重分配的平面距離分別為11.9 m、26.5 m、23.3 m、19.5 m、13.8 m、20.7 m，為更好地表現出分配之后的距離變化，特征水深點的距離重分配結果如表2 所示。

表2 特征水深點的距離重分配表

依據CIDW 算法，權函數的冪指數為2，周圍水深點數為12，P 點的內插值為14.05 m。依據本文提出的DCIDW 插值算法，利用式（7） 計算P點的內插值為13.95 m，此值與實測水深值更接近，減弱了樣本水深點變化不均勻對內插水深值的影響。

為進一步驗證特征采樣點的距離重分配算法的可行性，對本文選取的整體海區(qū)的樣本水深數據進行實驗，研究區(qū)域樣本水深數據共1077 組，其中862 組作為訓練集樣本進行建模，215 組作為測試集樣本進行驗證。研究區(qū)域平均水深為16.37 m，最大水深為29.48 m，最小水深為8.49 m。

本文選取的215 組待插值水深點中，有105 組進行了距離重分配，110 組未進行距離重分配，插值結果有兩類情況，如圖10 所示。圖10（a）中由于實際水深周圍地形的復雜性，本文提出的DCIDW插值方法能有效接近實際水深值，相較于CIDW 方法能有效逼近真實海底地形，圖10（b）中本文提出的DCIDW 插值算法相較于實際水深更淺，而CIDW 算法相較于實際水深更深，給航行安全帶來隱患。將本文DCIDW 算法和CIDW 算法的插值結果與真實水深作差（圖11），分析可知，本文提出的方法在復雜海底地形情況下相較于傳統(tǒng)方法能保證待插值水深點的內插質量。

圖10 插值水深結果分類圖

圖11 插值水深與實際水深誤差示意圖

將本文提出的算法與CIDW、自然鄰近點插值（natural neighbor interpolation，NNI）、樣條函數插值（spline function interpolation） 算法進行對比分析，其中CIDW 參考點數量選取12 個，本文提出的算法待插值點的參考數量為二階相鄰點，兩種IDW 插值算法冪指數均取值為2，待插值點從樣本點中隨機篩選，實驗采用交叉驗證方法驗證算法的誤差精度，選取平均絕對誤差（mean absolute error，MAE）、均方根誤差（root mean square error，RMSE）來表達待插值點插值結果的絕對誤差，相對誤差（relative error，RE）來表達待插值點插值結果的相對誤差。

實驗結果如表3 所示。為了直觀地顯示4 種水深插值方法的實驗結果誤差，將4 種水深插值方法對應的精度評價指標MAE 和RMSE 繪制成柱形圖，如圖12 所示。

表3 四種水深插值方法的實驗結果

圖12 四種水深插值方法實驗結果

分析圖12 可以發(fā)現，CIDW、NNI 和Spline 插值算法的MAE 分別為0.46 m、0.44 m、0.48 m，NNI 插值精度比CIDW、Spline 插值算法都要高，而本文提出的顧及特征水深點距離重分配DCIDW插值算法的MAE 在4 種算法中最低，為0.43 m；CIDW、NNI 和Spline 插值算法的RMSE 分別為0.59 m、0.59 m、0.63 m，Spline 插值算法的RMSE為0.63 m，相比較于CIDW、NNI 兩種插值算法，插值精度較低，而本文提出的顧及特征水深點距離重分配DCIDW 插值算法的RMSE 在4 種算法中最低，為0.53 m。通過對圖12 分析可知，本文方法在精度上優(yōu)于其他三種方法。相比于其他方法，本文方法進行待插值水深點插值時進行了特征提取和距離重分配，降低了空間非平穩(wěn)變化的影響，插值精度有所提高。

3 結語

本文針對CIDW 在復雜海底地形適應性差的局限性，提出了一種顧及特征水深點距離重分配的反距離加權插值算法。該算法顧及樣本水深點的周圍特征點對插值的影響，根據特征點距離分配因子對距離重新分配，從而提高了反距離加權插值算法在復雜海底地形中的精度。最后通過實測水深數據實驗驗證了本文方法的優(yōu)越性和可行性，該方法能有效地減弱在待插值水深點插值過程中樣本水深點鄰近點水深值分布不均勻對內插值的影響，與已有方法相比提高了插值結果的精度。

當然，本文僅針對一些典型區(qū)域進行了驗證，對于其他地形變化復雜的區(qū)域，本方法還需要進一步的研究。另外，本文僅單方面地考慮了特征水深點對待插值水深點水深值的約束，對于待插值水深點與等深線、地形線不同因素相互間的協(xié)調約束，有待進一步地探索。