多源水深數據融合的近海數字水深模型構建

陳義蘭,唐秋華,劉曉瑜,王燕紅

(1.自然資源部 第一海洋研究所,山東 青島 266061;2.自然資源部 海洋測繪重點實驗室,山東 青島 266590)

數字水深模型(Digital Bathymetric Model,DBM)是在一定范圍內描述水深信息的數據集,用于實現對水下地形地貌形態的數字化表達[1-4]。數字水深模型是海洋基礎地理信息數字成果的重要內容,是海洋經濟建設、海洋災害防治、海洋生態環境以及海洋科學研究等的基礎地理信息資源。我國近海的水深測量均為不同時間采用不同的單波束和多波束設備進行,測量比例尺不盡相同,導致數據源密度不同,所以給近海DBM的建立造成較大的困擾。為了最大程度地利用這些數據,有必要尋求合適的方法把這些多源多時相的水深數據進行集成和融合建立DBM,以滿足現階段對DBM的需求。

多源多時相水深數據融合建立DBM是目前國際上水深模型建立研究的熱點,國內外學者對多源多時相水深數據建立DBM有著諸多研究。目前國際上對多源數據建立DBM的插值方法主要采用連續曲率張力樣條法、三角網法、克里金法以及最近鄰法[5-13]。General Bathymetric Chart of the Oceans(GEBCO)組織全球相關專家研制了全球海陸一體數字高程和水深模型[7-8],目前的版本為2020年發布的GEBCO_2020網格,是全球海洋工作者廣泛采用的模型,它的建模方法也受到了較廣泛的關注。該模型中陸地采用SRTM的數據,海洋部分采用全球公開的多波束、單波束測深數據,以及衛星反演的水深數據,用連續曲率張力樣條法對數據進行插值,采用“移去-恢復”法進行模型融合,從而構建全球的數字高程-水深模型。B?CKSTR?M[9]基于單波束、多波束、海圖數據和衛星反演數據建立了瑞士韋特恩湖的數字水深模型,認為連續曲率張力樣條法是最合適的插值方法。Amante[10]評估了不同插值技術構建DBM的精度,認為連續曲率張力樣條法是最精確的插值技術。我國的研究者也采用張力樣條插值算法建立DBM[11-13]。可見,學者們對多源數據的DBM的建立采用最廣泛的插值方法是連續曲率張力樣條法,但是對建模的精度以及影響因素鮮少給出定量的客觀評價和系統分析。目前多源水深數據建立DBM的研究主要集中在深海大洋,沒有明確的精度要求,數據分辨率及精度和我國近海的情況不盡相同,它采用的方法是否一定適合我國近海需求,值得探討。作為國家地理信息數字成果的近海淺水區DBM建立是我們關注的重點,其對分辨率和精度要求有明確的規定,所以建模精度評價是不能缺少的,而對建模影響因素的分析對提高建模精度有很大的意義。

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

1 數字水深模型構建

1.1 插值方法

1.1.1 克里金(Kriging)法

Kriging數學模型是基于地統計學的網格插值方法,在許多領域有廣泛應用。這種方法與傳統插值方法的不同之處在于估計元觀測樣本數值時,不僅考慮待插值點與鄰近有觀測數據點的空間位置,還考慮了已知點之間的位置關系,而且利用已有觀測值空間分布的結構特點,使其估計比傳統方法更精確,更符合實際。Kriging方法用協方差函數和變異函數來確定高程變量隨空間距離變化的規律,在有限區域內對區域化變量進行無偏最優估計,使內插函數處于最優狀態[4,14]。

式中:z0為待插值點;zi為第i個采樣點的實測值;wi為第i個采樣點的權重系數,取決于測量點、預測位置的距離和預測位置周圍的測量值之間空間關系的擬合模型。本文采用線性模型為變異函數擬合模型。克里金法有較大的計算要求,需要復雜的編程,在處理較大數據量時效率相對較低。該方法被廣泛用于數字地形模型的建立,并被認為是精度較高的模型建立方法[15-18]。

1.1.2 連續曲率張力樣條插值法

Smith和Wesse[19]通過在最小曲率網格算法中引進張力參數,將最小曲率格網插值算法概括為一個更具普遍適用性的算法,即連續曲率張力樣條法,該方法可以消除最小曲率法存在較大的波動和無關變形點缺點[20]:

式中:TI為張力參數,下標I表示擬合區域內部,TI可在[0,1)區間取值。當TI=0時,式(2)為最小曲率格網化方程。邊界條件為:

式中:TB為邊界張力參數,在[0,1)區間取值。

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

1.1.3 狄洛尼三角網

狄洛尼三角網法將三角不規則網中的已知數據點連接起來。“空外接圓法則”是狄洛尼三角網生成原理,這個法則指如果三角網中每個三角形的外接圓只包含建立三角網的3個點,不包含點集中的任何其他點,則被認為是狄洛尼三角網[2]。狄洛尼三角網是目前受廣泛應用的三角網構網形式,它結構良好、數據結構簡單、存儲效率高,適應各種分布密度的數據。該方法是建立數字水深模型的常用插值方法,在許多研究中取得了較好的結果[17-18,21-22]。

1.2 模型融合方法

無論采用哪種插值方法構建模型,都會出現相同的問題,即由于數據密度的不一致,若只采用統一分辨率的網格大小會導致:網格過大,體現不出高密度數據表達細節的優勢;網格過小,低密度數據區內插會產生空白或者突變的偽值。為解決這個矛盾,最好的辦法是設計一個可變分辨率網格的模型,根據數據密度調整網格大小,使網格的大小和數據密度匹配,但是目前大多數的地形分析和可視化軟件不支持可變分辨率網格模型。所以,需要將不同分辨率的模型進行融合,形成一種恒定分辨率的模型。多數文獻選擇采用“移去-恢復”法對模型進行融合來解決上述問題,取得了較好的效果[11-13,23-24]。“移去-恢復”法首先根據數據密度分別建立低分辨率和高分辨率的兩種網格模型,低分辨率網格保證稀疏數據有正確的內插值且填充空白,但是移去了地形細節,高分辨率的網格保證高密度數據的細節,但是在稀疏區會出現空白;然后將低分辨率的網格重采樣內插成和高分辨率網格分辨率一致的網格作為基礎網格,將基礎網格和高分辨率網格比較得出二者差值網格,用基礎網格加上差值網格得到最終的恢復細節的網格。這種方法針對數據密度只分為2種類型,如果研究區的數據密度有多種,則可以采用疊加融合方法。本文采用一種疊加多分辨率網格的方法建立恒定網格大小的模型,該方法是在“移去-恢復”法的基礎上發展而來。

疊加多分辨率網格方法是根據研究區數據密度分成多個不同分辨率的網格,將不同分辨率網格模型進行疊加融合,再次內插生成高分辨率的網格。具體方法:①分別利用源數據生成不同分辨率的模型,在稀疏數據源處,高分辨率網格中會存在大量的空值,隨著網格的增大,空值會逐漸減少,直到填滿源數據包圍的所有空隙;②把這些網格疊加成多分辨率的網格,最后插值生成高分辨率的網格模型。在疊加的過程中,要求保留高分辨率的數據,在沒有高分辨率值的空值節點,采用低一級分辨率的網格節點值,具體流程見圖1。這種方法可以針對多源不均勻數據獲得更高分辨率網格模型,網格分辨率可以由具有最高源數據密度的區域確定,最大程度的保留了高密度數據的細節,又盡可能減少了稀疏數據源區域內插值偽影的產生。因為當稀疏的數據直接生成高分辨率的模型時,一些標準的插值算法無法產生可靠的結果[23]。如果先根據數據稀疏程度生成與之匹配的網格大小的模型,再采樣,出現假數據的機率會降低。

圖1 疊加融合多分辨率網格流程Fig.1 Process of stacking and fusing grids with multi-resolutions

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

圖2 研究區位置和數據源Fig.2 Study area location and data sources

本研究的數據來源主要為2008—2010年我國近海海洋綜合調查與評價專項實測的部分多波束和單波束數據,數據分布情況見圖2。專項調查數據精度一致,滿足相關規范要求,深度基準為理論深度基準面。

研究區的水深數據密度極不均勻,多波束數據分辨率10～50 m,單波束數據分辨率500～5 000 m,使用這些數據構建成一個整體數字水深模型,需要適合的模型插值方法和融合方法。

3 結果與分析

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

3.1 結果

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

圖3 研究區海底地勢Fig.3 Submarine topographic map of the study area

3.2 精度評價

(1)評價方法

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

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

(2)評價指標

本研究的評價指標采用均方根誤差(Root Mean Squared Errors,RMSE),誤差計算公式如下[26]:

式中:σ為研究區模型計算的RMSE;di為模型計算值和原始值的深度差;Z(i)m為同一位置模型內插計算的水深值;Z(i)o為抽取的原始水深數據值。表1為不同數據源3種插值方法計算的插值誤差。

表1 不同插值方法和數據源的均方根誤差(m)Table 1 RMSE statistics of different Interpolation methods and data sources(m)

3.3 討論分析

(1)數據密度對DBM的影響

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

(2)地形對DBM的影響

文獻[10]研究了曲率大小對插值精度的影響,但是只采用了部分剖面的曲率來分析,本研究采用地形變異性來研究地形對插值精度的影響,采用區域數據來分析,讓更多的數據參與分析,更具客觀性和代表性。地形變異性被認為是DBM不確定性的主要來源[27],為了探討其對DBM精度的影響,本文選取研究區4個地形特征不同的區域(圖3),并計算這些區域的水深變異系數(coefficient of variation,CV)來衡量地形變異性。水深變異系數對海底地形起伏變化情況有較好的反映,系數越大,表明該區域的地形起伏越大,地形越復雜。水深變異系數(CV)公式[27]:

式中:zi為區域內點i的水深值;為區域水深平均值;n為區域內總的水深點數。

計算4個試驗區域DBM建立的均方根誤差,探討變異系數和均方根誤差的關系即可知道地形變異性對DBM建立精度的影響。表2為不同區域(圖3)的CV值和插值均方根誤差結果,區域1和區域2為多波束區域,區域3和區域4為單波束區域。表2的均方根誤差值包含了地形和數據密度的影響,故把數據密度不同的區域分開比較,區域1和區域2進行比較,區域3和區域4比較,以盡量消除數據密度帶來的誤差影響。區域1的CV值大于區域2,區域1的地形復雜程度要高于區域2,相應的,區域1插值誤差大于區域2。區域3的CV值大于區域4,表明區域3的地形起伏較大,區域4的插值中誤差遠遠小于區域3,說明在數據稀疏的情況下,地形對插值精度的影響相當大。從以上分析得知,CV對DBM的誤差有正向影響,不論采用何種插值方法,地形變異越大,其總體插值誤差越大,水深插值的不確定性增加。

表2 不同區域的CV值和均方根誤差比較Table 2 Comparison of CV values and RMSEs in different regions regions

(3)插值方法對DBM的影響

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

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

表3 不同插值方法計算時間(s)Table 3 Computation time of different methods(s)

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

①DBM的建立精度和數據密度直接相關,數據密度密,插值精度高。地形對DBM建立精度有影響,地形復雜區,插值精度低。這是多源數據建立DBM的規律性結論。建立地形復雜度、密度和插值誤差的關系是需要進一步研究的方向,可以幫助更好地選擇插值方法。

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

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