基于DEM提取地性線研究

栗 業

（河北省地質礦產勘查開發局第二地質大隊，河北 唐山 063000）

1 研究背景概述

1.1 研究背景

數字高程模型（Digital Elevation Model， DEM）是依據高程數據對地形曲面的數字化模擬，是地表形態數字化的量化模型。DEM是地理空間地形分析的核心數據，是新一代的地形圖，同時也是GIS數據庫中重要的組成部分。DEM與傳統地形圖相比有著巨大優勢。首先，它可以在計算機中直接應用、計算，為各種輔助設計系統利用；其次，DEM是以數字形式存儲，便于數據的錄入、修改、更新、復制和管理，同時，可以方便地轉換成其他形式的材料文件和產品；最后，DEM中含有傳統地形圖中所無法容納和表達的垂直分布地物信息，即高程信息，所以DEM數據可以滿足眾多行業的發展需求。

1.2 研究意義

描述地貌形態的骨架線、控制線稱為地性線，包括山脊線、山谷線、鞍部、凹地等線性物體。山谷線和山脊線構成了地勢起伏變化的特征線，對地形地貌研究有重要意義。DEM數字高程模型中表達地形地貌的特征線，即地性線，其中包含了大量的地性線信息?；贒EM提取地性線的研究在區域水文分析、土木工程建設、地質勘探選定資源靶區、遙感影像自動配準等諸多領域具有重要的實踐應用價值和科學意義。

2 國內外研究現狀

多年來，國內外專家學者在數字高程模型中提取隱含著的山脊線和山谷線信息，及其如何將其快速自動提取，并應用于水文分析、地形分析、3D視域等多個研究領域方面，取得了一定的研究成果。

1984年，Yeoli提出了一個兩步算法。首先，使用斷面識別方法確定特征明顯的點，存儲于XMIN、YMIN、HMIN這3個矩陣中；然后從上述3個矩陣中去除已有山谷線上的點，用剩余的點中最大的點開始一條新山谷線的搜索。1991年，John Fairfield等在地表流水模型的基礎上，推出了隨機水流模型。該模型基本解決了DEM網格問題，但數據量成幾十倍的增加，計算耗時多，且只適用于區域較小且邊界較規則的范圍；1995年，Amnon Meisels等提出了骨架法。該算法是提取等值線上曲率值較大的點，然后把它們連接起來作為特征線。這種算法所需的參數少，缺點是不能在較大凹陷地區使用。

近年來，國內學者也在此領域取得了重要突破。張維軍應用道格拉斯-帕克提取山脊和山谷的特征點，并結合提取區域地形大致劃分和匯水線的幾何分析方法， 用概略的分水線和匯水線對特征點進行識別、歸類、順序提取各條山脊線和山谷線。姚慧敏等闡述了利用數學形態學的流域分割算法提取地性線的原理， 采用矢量柵格混合算法提取地性線的具體步驟， 最后實驗結果表明， 其提出的算法是可行的， 提取的地性線有較好的連續性。黃培之研究了特定實驗區域的山脊線和山谷線的幾何特性或物理特性的單一方面設計的提取山脊線和山谷線的算法，隨后又提出了一種基于地形表面流水分析與等高線幾何分析結合的提取地性線的方法。該方法把等高線幾何分析的方法與地形表面流水模擬分析的方法有機地結合起來， 能夠克服各自所具有的弊端。

3 DEM與地性線的數學表達

3.1 DEM的數學表達及數據結構

DEM數字高程模型是地表形態的數字化表達。從狹義角度說DEM是區域地表正高的數字化表達模型；廣義上來講，它是時空中的地理對象表面正高高度的數字化表達模型。

數學意義上的DEM是表示區域A上的三維向量有限序列，函數表達為：

式中，（Xi，Yi）是平面坐標，Zi是（Xi，Yi）對應的高程。

目前，DEM主要有4種數據結構：①不規則離散點數據結構。不規則離散點是最簡單的組織形式，它由不規則的離散點表示地表形態。這些點是通過測量獲得的獨立原始特征點數據。②規則格網（Grid）數據結構。在垂直、水平兩個方向上以相等間隔逐行逐列記錄每個網格的高程值稱為規則格網DEM。為方便行列號和平面坐標的轉換，規則網格中還標識了西南角坐標值、網格間距等。為壓縮存儲空間，常將規則網格每個單元值減去區域內最小高程值，再乘以一個常數是規則網格值全變為整數。規則格網DEM通常包含記錄西南角起點坐標、坐標類型、行列數的頭文件和記錄行列高程矩陣的數據體。規則格網的編碼結構有行程編碼、塊狀編碼、四叉樹編碼。③不規則三角網（TIN）數據結構。不規則三角網就是一系列互不重疊互不相交的三角形連接在一起表示地表形態。不規則三角網一般采用鏈表存儲，由結點和三角形頂點兩組記錄組成。TIN是典型的矢量拓撲結構。TIN克服了規則格網的缺點，數據量精簡，具有可變的分辨率，能表達復雜的地形。但該數據結構、算法都較復雜，管理困難，計算耗時長，適用于高精度小范圍的地形建模。④等高線數據結構。等高線是最常用的地形表達形式，在地形圖上高程相等的點所連成的閉合曲線即為等高線。它直觀地反映了高程、坡度、坡向、山脈走向等基本地貌。等高線以坐標點對序列形式存儲，在利用等高線重建地表模型時容易丟失細節信息，而且單條等高線無法反映地貌。

3.2 地性線的數學表達

3.2.1 山脊線

在實際地形地貌中條帶狀隆起的頂部線即山脊線。曲面上二階方向導數取得極小值，且一階導數為0的點是山脊點。山脊線則是相鄰山脊點的連線。

假設Ta和Tb正負型反斜坡，Lab=Ta∩Tb是兩斜坡的交線。在交線上存在一點P，其高程為Hp，則對任意iLab，存在CMi=Bi∩Ta∩Tb，使 得Hp-Hi≥0，得 出T=Ta∪Tb為山脊，Lab為山脊線或稱分水線，Bi為過點i的平面，CMi為山脊線等高線。如果地形曲面T的表達式為Z = A00+A10X+ A01Y+A11XY+A20X2+A02Y2+……（A20<0；A02<0且A11≠0） 時，判斷點P是否為山脊上一點，需滿足任意下面條件：坡度α=0；最大曲率Cmax>0；且最小曲率Cmin=0或坡度α>0；且斷面曲率Cs>0。

3.2.2 山谷線

山脊線和山谷線是兩個完全對應相反的概念。曲面上二階方向導數取得極大值，且一階導數為0的點是山谷點。山谷線則是相鄰山脊點的連線。山谷線收集兩側的水，屬于匯水區域，易形成河流。山谷線是純粹的局部地貌概念，因為它在端點處的一小段范圍內可能不會形成水系。由于它的局部性質，它比河流更適合，更便于計算機提取。所有我們近似地把它當成河流的一部分。

假設Ta和Tb正負型反斜坡，Lab=Ta∩Tb是兩斜坡的交線。在交線上存在一點P，其高程為Hp，則對任意iLab，存在CMi=Bi∩Ta∩Tb，使 得Hp-Hi≤0，得 出T=TaTb為 谷底，Lab為山谷線或稱匯水線，Bi為過點i的平面，CMi為山谷線等高線。如果地形曲面T的表達式為Z=A00+A10X+ A01Y+A11XY+A20X2+A02Y2+……（A20>0； A02>0且A11≠0） 時，判斷點P是否為山谷上一點，需滿足任意下面條件：坡度α=0；最大曲率Cmax=0；且最小曲率Cmin<0或坡度α>0；且斷面曲率Cs<0。

4 基于DEM提取地性線的數學算法

4.1 基于圖像處理的算法

DEM實質是一種柵格數據形式，網格大小即為像素大小，高程值為灰度值。DEM數據可以按數字圖像處理的方法來設計算子。在DEM中高程值越大，灰度越大，山脊線上的點比其他地方的點都要亮。按照各種濾波算子就能進行邊緣提取，將地形特征線提取出來。Chakrevavanieh提出用高斯算子從DEM上自動提取山脊線和山谷線上的點，并將其用于規則格網的數據壓縮。唐心紅等在提取地性線前利用小波分析進行多尺度表示，在多尺度DEM中提取。余生晨利用傅里葉變換在頻域中計算山脊線和山谷線，但該方法計算量大?；趫D像處理的技術還有以下兩個方面的難點：排除DEM中的噪聲影響；所提取特征點連接成線的算法問題。

4.2 基于地表幾何形態的算法

斷面極值法是地表幾何形態分析的基礎算法。該算法的基本思路是：在DEM垂直和水平橫截面上的極值，極大值就是分水點，極小值是匯水點；然后根據相應的準則劃分極值點歸屬形成地性線。地表幾何形態分析算法由于只判斷兩個斷面，提取過程中會出現遺漏現象。而且閾值過大會丟失特征點，形成的地性線也就斷裂較多，長度較短，閾值過小，特征點會出現過度提取現象，地性線冗余。

Fowler改進該算法，提出用4柵格掃面標記山脊點山谷點，這4柵格單元分別為（i， j）、（i， j+1）、（i+1， j）、（i+1， j+1），相比于8鄰域算法效率高而且不會遺漏山脊點和山谷點，缺陷是如果DEM采樣點密集會提取多余的特征點。Jenson認為任意一柵格與相鄰8單元格形成4對位置，至少一對高程都比該點高程值大，則該點可能是山谷點，如果至少有一對高程點比該點小，該點就可能為山脊點。這種改進算法易于編程實現，但產生地性線的斷裂較多。針對特征點丟失現象，黃培之提出再增加一組對角線上斷面的方法來克服這一缺陷，研究實驗證明該法有效，但同時增加了特征點噪聲。

4.3 基于地形表面流水物理模型分析算法

根據水往低處流，最終順谷而流，匯集成水系網的自然規律，Chakrevavanieh率先提出模擬地表流水物理模型。D8法是該模型中的最常用算法：假設只有一個單一的網格的水流進八個相鄰的網格單元。水流方向以最陡坡度法來確定，即在3×3的DEM柵格上，計算中心柵格與各相鄰柵格間的距離權落差，取距離權落差最大的柵格為中心柵格的流出柵格。以數值表示每個單元的流向，變化范圍是1~255。然后匯流柵格圖由流向柵格圖生成。匯流柵格上每個單元的值代表上游匯流區匯入該單元的流入路徑數較大者。從流域匯流柵格圖中可以輕易地提取流域各種參數特征值。在實際應用時需要解決以下問題：①水道起始點位置的確定；②偽負地形的識別及填充；③凹陷與平坦處水流方向的設定。而且基于地表水流模擬方法存在兩個缺陷：①洼地的填充和平坦地區水流方向的設定十分麻煩；②該方法所計算的匯水量與高程有關，計算的結果為高程值大的區域匯水量小，高程值小的區域匯水量大，由于水流遵從高處連續性的向下流動的自然特征，這就使得在一個河道水系或山谷線上的相對位置較高的某一點因其匯水量較小而被遺漏，而再非山谷處處于低位置的某點因其匯水量較大而被誤判為山谷線上的點。在后來的研究中，許多學者提出利用平滑DEM，墊高填平凹陷區，標識洼地集水區的方法解決上述問題，并取得了良好效果。

4.4 基于地形表面幾何形態分析和地表流水物理模型分析的算法

由于地表幾何形態分析和地表流水分析模型的算法有一些不足之處，所以專家們提出了將兩種模式結合提取地形線。這種算法的基本思路是首先按照地表流水模型在稀疏的DEM中概略提取地性線；然后對周圍地域進行幾何地理形態分析，在精確定位地性線。該方法的關鍵是求出概略地性線與DEM網格的交點，在該點周圍的小范圍內進行幾何分析，找出正交方向地形斷面上的高程極值點，正是這一點為精確點。這種綜合算法在特征點提取閾值的選取、DEM網格大小的確定方面還有待進一步研究。

4.5 平面曲率與坡形組合法

求取DEM平面曲率和地表的正負地形，正地形上曲率的大值是山脊，負地形上曲率的大值是山谷。這種方法中平面曲率的大小可以用來調教提取寬度，方法簡便效果好，但提取的圖形是柵格形式，在轉為矢量線性方面還需進一步研究改進。

5 結語

綜上所述，基于DEM提取地性線是項復雜的工作，目前提取算法相對成熟，但不同算法各有利弊，針對實際工作的不同地域地形和收集到的DEM數據質量，應選擇最佳的算法進行提取。