數(shù)字水深模型(DBM)是在一定范圍內(nèi)描述水深信息的數(shù)據(jù)集，用于實現(xiàn)對水下地形地貌形態(tài)的數(shù)字化表達。數(shù)字水深模型是海洋基礎(chǔ)地理信息數(shù)字成果的重要內(nèi)容，是海洋經(jīng)濟建設(shè)、海洋災(zāi)害防治、海洋生態(tài)環(huán)境以及海洋科學(xué)研究等的基礎(chǔ)地理信息資源。我國近海的水深測量均為不同時間采用不同的單波束和多波束設(shè)備進行，測量比例尺不盡相同，導(dǎo)致數(shù)據(jù)源密度不同，所以給近海DBM的建立造成較大的困擾。為了最大程度地利用這些數(shù)據(jù)，有必要尋求合適的方法把這些多源多時相的水深數(shù)據(jù)進行集成和融合建立DBM，以滿足現(xiàn)階段對DBM的需求。

多源多時相水深數(shù)據(jù)融合建立DBM是目前國際上水深模型建立研究的熱點，國內(nèi)外學(xué)者對多源多時相水深數(shù)據(jù)建立DBM有著諸多研究。目前國際上對多源數(shù)據(jù)建立DBM的插值方法主要采用連續(xù)曲率張力樣條法、三角網(wǎng)法、克里金法以及最近鄰法。GEBCO組織全球相關(guān)專家研制了全球海陸一體數(shù)字高程和水深模型，目前的版本為2020年發(fā)布的GEBCO_2020網(wǎng)格，是全球海洋工作者廣泛采用的模型，它的建模方法也受到了較廣泛的關(guān)注。該模型中陸地采用SRTM的數(shù)據(jù)，海洋部分采用全球公開的多波束、單波束測深數(shù)據(jù)，以及衛(wèi)星反演的水深數(shù)據(jù)，用連續(xù)曲率張力樣條法對數(shù)據(jù)進行插值，采用“移去-恢復(fù)”法進行模型融合，從而構(gòu)建全球的數(shù)字高程-水深模型。

B?CKSTR?M基于單波束、多波束、海圖數(shù)據(jù)和衛(wèi)星反演數(shù)據(jù)建立了瑞士韋特恩湖的數(shù)字水深模型，認(rèn)為連續(xù)曲率張力樣條法是最合適的插值方法。Amante評估了不同插值技術(shù)構(gòu)建DBM的精度，認(rèn)為連續(xù)曲率張力樣條法是最精確的插值技術(shù)。我國的研究者也采用張力樣條插值算法建立DBM。可見，學(xué)者們對多源數(shù)據(jù)的DBM的建立采用最廣泛的插值方法是連續(xù)曲率張力樣條法，但是對建模的精度以及影響因素鮮少給出定量的客觀評價和系統(tǒng)分析。目前多源水深數(shù)據(jù)建立DBM的研究主要集中在深海大洋，沒有明確的精度要求，數(shù)據(jù)分辨率及精度和我國近海的情況不盡相同，它采用的方法是否一定適合我國近海需求，值得探討。作為國家地理信息數(shù)字成果的近海淺水區(qū)DBM建立是我們關(guān)注的重點，其對分辨率和精度要求有明確的規(guī)定，所以建模精度評價是不能缺少的，而對建模影響因素的分析對提高建模精度有很大的意義。

多源水深數(shù)據(jù)最大的特點是水深數(shù)據(jù)密度不均勻。為提高DBM的建模精度，建立DBM需要分為兩步：首先，采用合適的插值方法建立初步的不同分辨率的DBM；然后，采用合適的模型融合方法建立最終的DBM。在數(shù)據(jù)測量精度一定的情況下，插值方法的選擇對DBM的建立尤為重要。為探討淺海區(qū)多源水深數(shù)據(jù)建立DBM的問題，本文以渤海海區(qū)的多波束、單波束數(shù)據(jù)為例，選取目前常用的克里金(Kriging)法、連續(xù)曲率張力樣條插值法和狄洛尼(Delaunay)三角網(wǎng)三種插值方法分別建立多分辨率的DBM，針對數(shù)據(jù)密度不同采用疊加多分辨率網(wǎng)格模型融合方法建立最終的DBM，并評價不同插值方法在不同條件下建模的精度，分析影響DBM不確定性因素，從而給出適合我國近海海區(qū)多源水深數(shù)據(jù)建立DBM的合理化技術(shù)建議。

⒈插值方法

⑴克里金(Kriging)法

Kriging數(shù)學(xué)模型是基于地統(tǒng)計學(xué)的網(wǎng)格插值方法，在許多領(lǐng)域有廣泛應(yīng)用。這種方法與傳統(tǒng)插值方法的不同之處在于估計元觀測樣本數(shù)值時，不僅考慮待插值點與鄰近有觀測數(shù)據(jù)點的空間位置，還考慮了已知點之間的位置關(guān)系，而且利用已有觀測值空間分布的結(jié)構(gòu)特點，使其估計比傳統(tǒng)方法更精確，更符合實際。

Kriging方法用協(xié)方差函數(shù)和變異函數(shù)來確定高程變量隨空間距離變化的規(guī)律，在有限區(qū)域內(nèi)對區(qū)域化變量進行無偏最優(yōu)估計，使內(nèi)插函數(shù)處于最優(yōu)狀態(tài)。

 (1)

 式中：z₀為待插值點；z_i為第i個采樣點的實測值；w_i為第i個采樣點的權(quán)重系數(shù)，取決于測量點、預(yù)測位置的距離和預(yù)測位置周圍的測量值之間空間關(guān)系的擬合模型。本文采用線性模型為變異函數(shù)擬合模型。克里金法有較大的計算要求，需要復(fù)雜的編程，在處理較大數(shù)據(jù)量時效率相對較低。該方法被廣泛用于數(shù)字地形模型的建立，并被認(rèn)為是精度較高的模型建立方法。

⑵連續(xù)曲率張力樣條插值法

Smith和Wesse通過在最小曲率網(wǎng)格算法中引進張力參數(shù)，將最小曲率格網(wǎng)插值算法概括為一個更具普遍適用性的算法，即連續(xù)曲率張力樣條法，該方法可以消除最小曲率法存在較大的波動和無關(guān)變形點缺點：

式中：T_I為張力參數(shù)，下標(biāo)I表示擬合區(qū)域內(nèi)部，T_I可在[0，1)區(qū)間取值。當(dāng)T_I≒0時，式(2)為最小曲率格網(wǎng)化方程。邊界條件為：

 ⑶

式中：T_B為邊界張力參數(shù)，在[0，1)區(qū)間取值。

這種算法最主要的是確定張力參數(shù)，參數(shù)越大，格網(wǎng)化結(jié)果越平滑。該算法用來計算了全球廣泛使用的DBM，例如GEBCO發(fā)布的GEBCO_2019網(wǎng)格模型以及南北極的網(wǎng)格模型。Akkala等認(rèn)為最適合張力樣條插值的是不規(guī)則間隔的數(shù)據(jù)。根據(jù)眾多學(xué)者的經(jīng)驗[5-9，20]，張力參數(shù)為0.25～0.45，格網(wǎng)化結(jié)果既能反映局部的變化，又比較平滑，精度較高，本次的張力參數(shù)取0.35。

⑵狄洛尼三角網(wǎng)

狄洛尼三角網(wǎng)法將三角不規(guī)則網(wǎng)中的已知數(shù)據(jù)點連接起來。“空外接圓法則”是狄洛尼三角網(wǎng)生成原理，這個法則指如果三角網(wǎng)中每個三角形的外接圓只包含建立三角網(wǎng)的3個點，不包含點集中的任何其他點，則被認(rèn)為是狄洛尼三角網(wǎng)。狄洛尼三角網(wǎng)是目前受廣泛應(yīng)用的三角網(wǎng)構(gòu)網(wǎng)形式，它結(jié)構(gòu)良好、數(shù)據(jù)結(jié)構(gòu)簡單、存儲效率高，適應(yīng)各種分布密度的數(shù)據(jù)。該方法是建立數(shù)字水深模型的常用插值方法，在許多研究中取得了較好的結(jié)果。

⒉模型融合方法

無論采用哪種插值方法構(gòu)建模型，都會出現(xiàn)相同的問題，即由于數(shù)據(jù)密度的不一致，若只采用統(tǒng)一分辨率的網(wǎng)格大小會導(dǎo)致：網(wǎng)格過大，體現(xiàn)不出高密度數(shù)據(jù)表達細(xì)節(jié)的優(yōu)勢；網(wǎng)格過小，低密度數(shù)據(jù)區(qū)內(nèi)插會產(chǎn)生空白或者突變的偽值。為解決這個矛盾，最好的辦法是設(shè)計一個可變分辨率網(wǎng)格的模型，根據(jù)數(shù)據(jù)密度調(diào)整網(wǎng)格大小，使網(wǎng)格的大小和數(shù)據(jù)密度匹配，但是目前大多數(shù)的地形分析和可視化軟件不支持可變分辨率網(wǎng)格模型。所以，需要將不同分辨率的模型進行融合，形成一種恒定分辨率的模型。多數(shù)文獻選擇采用“移去-恢復(fù)”法對模型進行融合來解決上述問題，取得了較好的效果。“移去-恢復(fù)”法首先根據(jù)數(shù)據(jù)密度分別建立低分辨率和高分辨率的兩種網(wǎng)格模型，低分辨率網(wǎng)格保證稀疏數(shù)據(jù)有正確的內(nèi)插值且填充空白，但是移去了地形細(xì)節(jié)，高分辨率的網(wǎng)格保證高密度數(shù)據(jù)的細(xì)節(jié)，但是在稀疏區(qū)會出現(xiàn)空白；然后將低分辨率的網(wǎng)格重采樣內(nèi)插成和高分辨率網(wǎng)格分辨率一致的網(wǎng)格作為基礎(chǔ)網(wǎng)格，將基礎(chǔ)網(wǎng)格和高分辨率網(wǎng)格比較得出二者差值網(wǎng)格，用基礎(chǔ)網(wǎng)格加上差值網(wǎng)格得到最終的恢復(fù)細(xì)節(jié)的網(wǎng)格。這種方法針對數(shù)據(jù)密度只分為2種類型，如果研究區(qū)的數(shù)據(jù)密度有多種，則可以采用疊加融合方法。本文采用一種疊加多分辨率網(wǎng)格的方法建立恒定網(wǎng)格大小的模型，該方法是在“移去-恢復(fù)”法的基礎(chǔ)上發(fā)展而來。

疊加多分辨率網(wǎng)格方法是根據(jù)研究區(qū)數(shù)據(jù)密度分成多個不同分辨率的網(wǎng)格，將不同分辨率網(wǎng)格模型進行疊加融合，再次內(nèi)插生成高分辨率的網(wǎng)格。具體方法：①分別利用源數(shù)據(jù)生成不同分辨率的模型，在稀疏數(shù)據(jù)源處，高分辨率網(wǎng)格中會存在大量的空值，隨著網(wǎng)格的增大，空值會逐漸減少，直到填滿源數(shù)據(jù)包圍的所有空隙；②把這些網(wǎng)格疊加成多分辨率的網(wǎng)格，最后插值生成高分辨率的網(wǎng)格模型。在疊加的過程中，要求保留高分辨率的數(shù)據(jù)，在沒有高分辨率值的空值節(jié)點，采用低一級分辨率的網(wǎng)格節(jié)點值，具體流程見圖1。

圖1 疊加融合多分辨率網(wǎng)格流程

這種方法可以針對多源不均勻數(shù)據(jù)獲得更高分辨率網(wǎng)格模型，網(wǎng)格分辨率可以由具有最高源數(shù)據(jù)密度的區(qū)域確定，最大程度的保留了高密度數(shù)據(jù)的細(xì)節(jié)，又盡可能減少了稀疏數(shù)據(jù)源區(qū)域內(nèi)插值偽影的產(chǎn)生。因為當(dāng)稀疏的數(shù)據(jù)直接生成高分辨率的模型時，一些標(biāo)準(zhǔn)的插值算法無法產(chǎn)生可靠的結(jié)果。如果先根據(jù)數(shù)據(jù)稀疏程度生成與之匹配的網(wǎng)格大小的模型，再采樣，出現(xiàn)假數(shù)據(jù)的機率會降低。

渤海是我國的內(nèi)海，位于117°35'～121°10'E、37°07'～41°00'N，南北長約480km，東西寬約300km。渤海被遼寧、河北、天津和山東包圍，只有東南向通過渤海海峽與黃海相接。渤海水淺，平均水深18m，地形平緩。研究區(qū)位于渤海中西部，具體位置見圖2。研究區(qū)從渤海灣到渤中淺灘，水深0～39m，既有復(fù)雜的沙脊溝槽地形、人工地形，也有平緩的平原地形，其典型地形特征為研究海底水深模型構(gòu)建提供了依據(jù)。

圖2 研究區(qū)位置和數(shù)據(jù)源

本研究的數(shù)據(jù)來源主要為2008—2010年我國近海海洋綜合調(diào)查與評價專項實測的部分多波束和單波束數(shù)據(jù)，數(shù)據(jù)分布情況見圖2。專項調(diào)查數(shù)據(jù)精度一致，滿足相關(guān)規(guī)范要求，深度基準(zhǔn)為理論深度基準(zhǔn)面。

研究區(qū)的水深數(shù)據(jù)密度極不均勻，多波束數(shù)據(jù)分辨率10～50m，單波束數(shù)據(jù)分辨率500～5000m，使用這些數(shù)據(jù)構(gòu)建成一個整體數(shù)字水深模型，需要適合的模型插值方法和融合方法。

采用不同的插值方法分別對單波束、多波束以及單波束和多波束混合數(shù)據(jù)進行插值并融合模擬計算，計算不同情況下的插值精度，分析數(shù)據(jù)密度對插值精度的影響，選擇典型地形區(qū)域計算插值誤差，根據(jù)插值誤差的分布特征，分析地形對DBM的影響，從而進一步分析插值方法的優(yōu)劣。

⒈結(jié)果

采用克里金法建立的研究區(qū)DBM形成的海底地勢見圖3，圖3展示了研究區(qū)的地形起伏形態(tài)。研究區(qū)西部為渤海灣曹妃甸近岸復(fù)雜地形區(qū)，水深變化劇烈，海底沖刷溝槽和潮流沙壩相間分布，曹妃甸外的海底深槽水深達39m，槽底呈平緩的波狀起伏形態(tài)。研究區(qū)中東部為渤海中央盆地，地形平坦開闊，坡度平緩，其東部為渤中淺灘，灘頂平緩。

圖3 研究區(qū)海底地勢

⒉精度評價

⑴評價方法

DBM的誤差主要有測量誤差和插值誤差，本文的精度評價主要對插值誤差結(jié)果進行評價，不涉及測量誤差。

精度評價方法采用交叉驗證法。交叉驗證法是抽取一定數(shù)量的原始深度值作為驗證數(shù)據(jù)，采用保留的數(shù)據(jù)作為訓(xùn)練數(shù)據(jù)來進行插值計算。將計算的深度和驗證的原始深度數(shù)據(jù)進行比對計算深度差異。分別從單波束、多波束源數(shù)據(jù)以及單波束和多波束混合數(shù)據(jù)隨機分離20%的數(shù)據(jù)作為驗證數(shù)據(jù)，采用剩下的80%的數(shù)據(jù)進行插值建模，然后計算訓(xùn)練數(shù)據(jù)和驗證數(shù)據(jù)之間的誤差。

⑵評價指標(biāo)

本研究的評價指標(biāo)采用均方根誤差(RMSE)，誤差計算公式如下：

式中：σ為研究區(qū)模型計算的RMSE；d_i為模型計算值和原始值的深度差；Z_(i)m為同一位置模型內(nèi)插計算的水深值；Z_(i)o為抽取的原始水深數(shù)據(jù)值。表1為不同數(shù)據(jù)源3種插值方法計算的插值誤差。

表1 不同插值方法和數(shù)據(jù)源的均方根誤差(m)

⒊討論分析

⑴數(shù)據(jù)密度對DBM的影響

由表1可知，多波束數(shù)據(jù)的計算精度最高，混合數(shù)據(jù)次之，單波束數(shù)據(jù)的插值計算精度最低。多波束數(shù)據(jù)區(qū)域雖然地形復(fù)雜，但是由于數(shù)據(jù)密度遠高于單波束數(shù)據(jù)，插值精度明顯高于單波束數(shù)據(jù)區(qū)域。由此可見，插值精度與原始數(shù)據(jù)密度有關(guān)，數(shù)據(jù)密度對DBM的精度起著重要的作用，數(shù)據(jù)密度密，插值精度高。

⑵地形對DBM的影響

有關(guān)文獻研究了曲率大小對插值精度的影響，但是只采用了部分剖面的曲率來分析，本研究采用地形變異性來研究地形對插值精度的影響，采用區(qū)域數(shù)據(jù)來分析，讓更多的數(shù)據(jù)參與分析，更具客觀性和代表性。

地形變異性被認(rèn)為是DBM不確定性的主要來源，為了探討其對DBM精度的影響，本文選取研究區(qū)4個地形特征不同的區(qū)域(圖3)，并計算這些區(qū)域的水深變異系數(shù)(CV)來衡量地形變異性。水深變異系數(shù)對海底地形起伏變化情況有較好的反映，系數(shù)越大，表明該區(qū)域的地形起伏越大，地形越復(fù)雜。水深變異系數(shù)(CV)公式：                         ⑸

 式中：z_i為區(qū)域內(nèi)點i的水深值；z為區(qū)域水深平均值；n為區(qū)域內(nèi)總的水深點數(shù)。

計算4個試驗區(qū)域DBM建立的均方根誤差，探討變異系數(shù)和均方根誤差的關(guān)系即可知道地形變異性對DBM建立精度的影響。表2為不同區(qū)域(圖3)的CV值和插值均方根誤差結(jié)果，區(qū)域1和區(qū)域2為多波束區(qū)域，區(qū)域3和區(qū)域4為單波束區(qū)域。表2的均方根誤差值包含了地形和數(shù)據(jù)密度的影響，故把數(shù)據(jù)密度不同的區(qū)域分開比較，區(qū)域1和區(qū)域2進行比較，區(qū)域3和區(qū)域4比較，以盡量消除數(shù)據(jù)密度帶來的誤差影響。區(qū)域1的CV值大于區(qū)域2，區(qū)域1的地形復(fù)雜程度要高于區(qū)域2，相應(yīng)的，區(qū)域1插值誤差大于區(qū)域2。

表2 不同區(qū)域的CV值和均方根誤差比較

區(qū)域3的CV值大于區(qū)域4，表明區(qū)域3的地形起伏較大，區(qū)域4的插值中誤差遠遠小于區(qū)域3，說明在數(shù)據(jù)稀疏的情況下，地形對插值精度的影響相當(dāng)大。從以上分析得知，CV對DBM的誤差有正向影響，不論采用何種插值方法，地形變異越大，其總體插值誤差越大，水深插值的不確定性增加。

⑶插值方法對DBM的影響

插值方法是影響DBM精度的重要因素。由表1可知，不論哪種數(shù)據(jù)源密度，這3種方法中，克里金法精度最高。數(shù)據(jù)源為混合的不均勻數(shù)據(jù)和稀疏數(shù)據(jù)時，連續(xù)曲率張力樣條法比三角網(wǎng)法精度稍高，數(shù)據(jù)源為密集的多波束數(shù)據(jù)，三角網(wǎng)法比連續(xù)曲率張力樣條法精度稍高。由表2可知，在平坦地形區(qū)，3種插值方法的精度相當(dāng)，克里金法無論是在地形復(fù)雜區(qū)還是平坦區(qū)，它的精度都是最高的，而在高密度的多波束復(fù)雜地形區(qū)，三角網(wǎng)法比張力樣條法更具優(yōu)勢，在稀疏的單波束復(fù)雜地形區(qū)，張力樣條法表現(xiàn)優(yōu)于三角網(wǎng)法。表1和表2的統(tǒng)計結(jié)果一致表明：克里金法精度最高，尤其在數(shù)據(jù)稀疏且地形復(fù)雜區(qū)，克里金法在插值精度方面具有明顯的優(yōu)勢，這個結(jié)果也體現(xiàn)了克里金算法相對于其他算法的優(yōu)勢，它充分考慮空間變量相關(guān)性，更能客觀地反映自然地形規(guī)律。連續(xù)曲率張力樣條法相對于三角網(wǎng)法更適合稀疏的數(shù)據(jù)，體現(xiàn)了連續(xù)曲率張力樣條法在尊重源數(shù)據(jù)基礎(chǔ)上，能生成圓滑表面的特點，三角網(wǎng)法適合稠密數(shù)據(jù)，這也體現(xiàn)了三角網(wǎng)法尊重源數(shù)據(jù)的特點。

對于數(shù)據(jù)量較大的大范圍的多源數(shù)據(jù)的插值建模，在選擇插值方法時，計算時間有時也在考慮之列。比較相同條件下的相對計算時間可以為選擇更有效的插值方法提供參考，表3是不同數(shù)據(jù)密度不同插值方法計算100m分辨率的DBM的時間，為了更有效地統(tǒng)計計算時間，這個DBM計算沒有采用模型融合方法，直接采用數(shù)據(jù)進行插值計算。克里金法在精度上具有優(yōu)勢，但是在計算效率上比連續(xù)曲率張力樣條法和三角網(wǎng)法低。

表3 不同插值方法計算時間(s)

本文采用克里金法、連續(xù)曲率張力樣條法和狄洛尼三角網(wǎng)法三種插值方法建立水深數(shù)字模型，并對模型進行融合，通過分析DBM的精度以及影響DBM建立的因素，認(rèn)為源數(shù)據(jù)密度、地形特征和插值方法都對DBM精度有顯著影響，主要結(jié)論如下：

①DBM的建立精度和數(shù)據(jù)密度直接相關(guān)，數(shù)據(jù)密度密，插值精度高。地形對DBM建立精度有影響，地形復(fù)雜區(qū)，插值精度低。這是多源數(shù)據(jù)建立DBM的規(guī)律性結(jié)論。建立地形復(fù)雜度、密度和插值誤差的關(guān)系是需要進一步研究的方向，可以幫助更好地選擇插值方法。

②本文選用的3種插值方法各有優(yōu)勢，適用不同的數(shù)據(jù)情況，所得結(jié)果體現(xiàn)了各插值方法的特點。由于多源水深數(shù)據(jù)的復(fù)雜性，建立DBM的插值方法的選擇要考慮諸多因素，如果DBM精度是考慮的首要因素，克里金法插值精度最高，該方法是較好的選擇。雖然計算機能力強大，如果計算范圍過大，效率也可在考慮范圍之內(nèi)，例如全球范圍內(nèi)的多源數(shù)據(jù)，數(shù)據(jù)極其不均勻，綜合考慮效率和精度，數(shù)據(jù)稀疏不均勻時，連續(xù)曲率張力樣條法是較好的選擇，這也是全球數(shù)字水深模型建立最終選擇連續(xù)曲率張力樣條法的原因；如果數(shù)據(jù)稠密均勻且地形復(fù)雜，狄洛尼三角網(wǎng)法是合適的選擇，當(dāng)數(shù)據(jù)密度滿足分辨率要求時，例如全覆蓋的多波束數(shù)據(jù)，3種插值方法的計算精度差異不大，建議選擇簡單高效的插值方法。

③本文針對多源水深數(shù)據(jù)特點，采用疊加多分辨率網(wǎng)格的融合方法解決了模型恒定網(wǎng)格大小和多種數(shù)據(jù)密度的矛盾，對多分辨率模型融合具有借鑒意義。