等高線曲率計(jì)算方法研究

張 亮，楊勤科，2，蘭 敏

（1.中國(guó)科學(xué)院 水利部 水土保持研究所，陜西 楊凌712100；2.西北大學(xué) 城市與環(huán)境學(xué)院，西安710069；3.中煤西安設(shè)計(jì)工程有限責(zé)任公司，西安710054）

數(shù)字高程模型（Digital Elevation Models，DEMs）是一組表示景觀內(nèi)基于某一基準(zhǔn)面的高程點(diǎn)在空間上的分布狀況的數(shù)字陣列［1］。DEM包括規(guī)則格網(wǎng)、不規(guī)則三角網(wǎng)和等高線三種模型，規(guī)則格網(wǎng)DEM數(shù)據(jù)結(jié)構(gòu)簡(jiǎn)單的優(yōu)點(diǎn)使其在數(shù)字地形分析中的應(yīng)用最為廣泛［2－3］?；诘匦螠y(cè)量成果數(shù)據(jù)（等高線、高程點(diǎn)和河流等）通過(guò)插值建立規(guī)則格網(wǎng)DEM，是目前最為成熟和經(jīng)濟(jì)實(shí)用的DEM建立方法［4－6］。在基于等高線插值建立格網(wǎng)DEM過(guò)程中，只有設(shè)置合理的分辨率才能保證地形信息不損失或少損失，又不引入過(guò)多的數(shù)據(jù)冗余［7－9］。曲率是表示曲線彎曲程度的參數(shù)，它可以直接反映等高線的復(fù)雜程度，是等高線信息量的重要指標(biāo)［10］。因而計(jì)算等高線曲率，可以通過(guò)衡量等高線建立DEM過(guò)程中信息的轉(zhuǎn)換率，為準(zhǔn)確設(shè)置DEM分辨率提供參考。為滿足實(shí)際應(yīng)用的需求，前人研究者們提出了多種地形曲率指標(biāo)及其算法［11］。然而這些算法均基于格網(wǎng)DEM，通過(guò)局部窗口的曲面擬合函數(shù)，來(lái)模擬地形表面的方式實(shí)現(xiàn)；由于等高線節(jié)點(diǎn)數(shù)據(jù)分布不規(guī)則等原因，如何基于等高線直接計(jì)算地形曲率則缺乏必要的研究。另外，格網(wǎng)DEM本身是等高線的派生數(shù)據(jù)，基于其計(jì)算的地形曲率則是二次派生數(shù)據(jù)，精度比直接基于等高線計(jì)算的結(jié)果低，所以研究用等高線計(jì)算地形曲率有其必要性。本文將討論矢量格式等高線曲率的計(jì)算方法，以期為DEM分辨率的確定提供科學(xué)依據(jù)，對(duì)數(shù)字地形分析中地形曲率的研究也具有理論意義。

1 研究方法

1.1 等高線曲率算法

等高線曲率，是指通過(guò)地表上某一點(diǎn)的等高面（水平面）與地表的交線（即等高線）在該點(diǎn)的曲率，曲率是衡量曲線彎曲程度的指標(biāo)；等高線曲率表達(dá)了地表物質(zhì)運(yùn)動(dòng)的匯合和發(fā)散模式。

為了量化表示等高線的彎曲程度，本文用等高線上相鄰三個(gè)節(jié)點(diǎn)擬合曲率圓，來(lái)計(jì)算等高線曲率。不在一條直線上的三點(diǎn)可以唯一確定一個(gè)圓，中間一點(diǎn)處的曲率用半徑的倒數(shù)表示。如三點(diǎn)共線（即等高線的局部為直線），則直接賦曲率值為零。

如圖1所示，設(shè)r為三點(diǎn)（Pi－1，Pi，Pi＋1）所組成三角形的外接圓半徑，那么1／r即可認(rèn)為是點(diǎn)Pi處的曲率。根據(jù)海倫公式、正弦定理、三角形面積公式可以推導(dǎo)出，圓擬合法計(jì)算等高線曲率的公式為：

式中：Cc——節(jié)點(diǎn) Pi處的等高線曲率；a，b，c——線段的Pi－1Pi，PiPi＋1，Pi＋1Pi－1長(zhǎng)度，等高線曲率單位為弧度每米（rad／m）。由于該方法以單條等高線上的節(jié)點(diǎn)為基礎(chǔ)來(lái)計(jì)算曲率，因而暫不能表示地形的凹凸變化。本文主要討論等高線曲率，為簡(jiǎn)便起見(jiàn)，下文中未指明的曲率一律特指等高線曲率。

圖1 等高線曲率計(jì)算示意圖

1.2 等高線曲率的計(jì)算

本研究基于Microsoft.NET平臺(tái)用ArcEngine組件庫(kù)開(kāi)發(fā)定制程序，用于計(jì)算等高線曲率。圖2為等高線曲率計(jì)算的數(shù)據(jù)流程圖，表1為數(shù)據(jù)點(diǎn)的存儲(chǔ)結(jié)構(gòu)，主要步驟包括：（1）對(duì)等高線數(shù)據(jù)（CNT）做拓?fù)浜瓦B接等方面的錯(cuò)誤修復(fù)，拓?fù)溴e(cuò)誤包括等高線疊合、交叉、自相交等，修復(fù)后的等高線圖層標(biāo)記為CNT＿conn；（2）使用ArcGIS的Feature Vertices To Points工具將CNT＿conn轉(zhuǎn)換成點(diǎn)圖層，并用Add XY coordinates工具提取等高線節(jié)點(diǎn)坐標(biāo)，結(jié)果標(biāo)記為CNT＿points；（3）編寫(xiě)程序，以CNT＿points中的“ORIG＿ID”字段（記錄了所在等高線的索引號(hào)）為依據(jù)讀取同一條等高線上相鄰三節(jié)點(diǎn)的坐標(biāo)數(shù)據(jù)，用公式（1）計(jì)算等高線曲率并存儲(chǔ)于CNT＿points屬性表的curvature字段中。

圖2 等高線曲率計(jì)算數(shù)據(jù)流程

表1 等高線節(jié)點(diǎn)數(shù)據(jù)結(jié)構(gòu)

1.3 試驗(yàn)數(shù)據(jù)

本文選取黃土丘陵區(qū)、黃土塬區(qū)和東北漫崗丘陵區(qū)三個(gè)典型試驗(yàn)樣區(qū)。黃土丘陵樣區(qū)（圖3a，樣區(qū)1）位于陜北安塞縣縣南溝小流域，平均坡度29.9°，梁地形比較簡(jiǎn)單平緩，溝坡則較復(fù)雜，海拔1 070.3～1 338.5m，面積2.41km2。黃土塬樣區(qū)（圖3b，樣區(qū)2）位于陜西長(zhǎng)武縣，平均坡度19.2°，塬面地形比較簡(jiǎn)單，溝坡地形比較復(fù)雜，海拔范圍977.9～1 197.0m，面積2.44km2。東北漫崗丘陵區(qū)（圖3c，樣區(qū)3）位于黑龍江省拜泉縣羅家屯溝小流域，平均坡度3.9°，整體呈波狀緩坡丘陵，少部分較陡峭，海拔范圍254.8～305.8m，面積2.40km2。研究選用1∶1萬(wàn)等高線作為試驗(yàn)數(shù)據(jù)，其中黃土丘陵區(qū)等高距5m（CNT1），黃土塬區(qū)等高距10m（CNT2），東北漫崗丘陵區(qū)等高距2.5m（CNT3）。

圖3 試驗(yàn)樣區(qū)等高線

1.4 地形類(lèi)型劃分及相關(guān)地形指標(biāo)提取

地形類(lèi)型劃分：用手工的方式分別勾繪黃土高原兩個(gè)試驗(yàn)樣區(qū)的溝沿線和塬邊線，將黃土丘陵區(qū)地形簡(jiǎn)單區(qū)分為梁和溝壑；黃土塬區(qū)分為現(xiàn)代溝坡、古代溝坡及塬面三種；根據(jù)坡度特征將漫崗丘陵樣區(qū)劃分為較緩坡和較陡坡丘陵（圖4）。這種劃分并非嚴(yán)格的地貌類(lèi)型劃分，僅用于分析比較不同地形中等高線曲率的統(tǒng)計(jì)分析。

圖4 地形分區(qū)圖

相關(guān)地形指標(biāo)提?。河脴訁^(qū)1等高線數(shù)據(jù)（5m等高距）建立2.5m分辨率格網(wǎng)DEM（DEM1），樣區(qū)2等高線（10m等高距）建立5m分辨率DEM（DEM2），樣區(qū)3等高線（2.5m 等高距）建立1m分辨率DEM（DEM3）。然后分別提取其坡度（SLOPE1，SLOPE2，SLOPE3）、坡 長(zhǎng) （SLPLEN1，SLPLEN2，SLPLEN3）、剖 面 曲 率 （PROCURV1，PROCURV2，PROCURV3）和 平 面 曲 率 （PLANCURV1，PLANCURV2，PLANCURV3），用于比較分析等高線曲率與相關(guān)地形指標(biāo)之間的空間格局分布關(guān)系。

2 結(jié)果與分析

2.1 空間格局分析

等高線曲率的計(jì)算結(jié)果（圖5）可見(jiàn)：（1）等高線比較彎曲的地方曲率值大，等高線比較平直的地方曲率值小；（2）在山谷線和山脊線附近等高線曲率計(jì)算值比較大，而一般斜坡上曲率值小。可見(jiàn)本文建議的圓擬合曲率算法計(jì)算的等高線曲率，結(jié)果是正確的。但仔細(xì)觀察可見(jiàn)，在比較平值的等高線上，曲率計(jì)算值略有波動(dòng)，這可能是由于等高線上的線段節(jié)點(diǎn)的局地誤差所致，對(duì)此有待進(jìn)一步研究。

圖5 試驗(yàn)樣區(qū)等高線曲率圖（局部）

按照前述地形類(lèi)型劃分（圖4）對(duì)等高線曲率進(jìn)行地形分區(qū)統(tǒng)計(jì)（表2），結(jié)果表明：（1）黃土丘陵區(qū)溝壑地形等高線曲率要普遍高于梁地形；（2）黃土塬區(qū)現(xiàn)代溝坡等高線曲率比古代溝坡和塬面均大；（3）漫崗丘陵區(qū)較陡坡地形也是比較緩坡地形等高線曲率大。而這種大小變化，均由等高線彎曲程度不同所致，這進(jìn)一步說(shuō)明本文提出的算法，計(jì)算結(jié)果符合實(shí)際情況。

2.2 統(tǒng)計(jì)分布分析

黃土丘陵區(qū)、黃土塬區(qū)和漫崗丘陵區(qū)等高線曲率相關(guān)統(tǒng)計(jì)數(shù)據(jù)（圖6）表明：（1）等高線曲率值集中在0.2以下（92%～98%），而且曲率集中于小值的趨勢(shì)依黃土丘陵區(qū)、黃土塬區(qū)、漫崗丘陵區(qū)的次序增強(qiáng)；（2）等高線曲率的平均值和范圍方面，黃土丘陵區(qū)均大于黃土塬區(qū)和漫崗丘陵區(qū)。這是由于黃土丘陵區(qū)等高線彎曲程度普遍較大，而黃土塬區(qū)在塬面和古代溝坡上，特別是漫崗丘陵區(qū)等高線相對(duì)比較平直。

2.3 與常見(jiàn)地形指標(biāo)的比較

提取常見(jiàn)地形指標(biāo)，包括坡度、坡長(zhǎng)、剖面曲率和平面曲率如圖7—9所示（等高線曲率圖空白為背景色），其與等高線曲率的對(duì)比可以看出：（1）等高線曲率與平面曲率的空間格局最為相似，特征線條沿地形特征線方向，都是在山脊線和山谷線處曲率大；但數(shù)值上并不完全相同。（2）在丘陵區(qū)的溝壑、塬區(qū)的現(xiàn)代溝和漫崗丘陵區(qū)的較陡坡區(qū)，等高線曲率、剖面曲率和坡度等數(shù)值均較大，但從空間上并不具有十分精確的對(duì)應(yīng)關(guān)系。（3）剖面曲率反映的是地形在垂直方向上二階導(dǎo)數(shù)，與坡度的變率相關(guān)，而坡度主要由等高線之間的距離決定，因此剖面曲率圖和坡度圖在空間格局和紋理特征上都最為相似。（4）等高線曲率與坡長(zhǎng)的空間分布關(guān)系似乎不太大。

表2 試驗(yàn)樣區(qū)等高線曲率分區(qū)統(tǒng)計(jì) rad／m

圖6 等高線曲率頻率曲線

圖7 黃土丘陵區(qū)各地形因子

圖8 黃土塬區(qū)各地形因子

圖9 漫崗丘陵區(qū)各地形因子

3 結(jié)論與討論

3.1 結(jié) 論

（1）本文設(shè)計(jì)的圓擬合算法，從曲率的數(shù)學(xué)含義入手，用曲線擬合的方法，可實(shí)現(xiàn)對(duì)等高線曲率的計(jì)算。結(jié)果表明在等高線比較彎曲的地方，等高線曲率比較大；從地形來(lái)看，溝壑地帶（溝沿線到坡腳線之間）等地形比較復(fù)雜的地方等高線曲率比較大。等高線曲率計(jì)算值的空間分布與地形特征具有良好的對(duì)應(yīng)關(guān)系，結(jié)果合理可信。

（2）不同復(fù)雜程度的樣區(qū)，曲率計(jì)算結(jié)果的空間和統(tǒng)計(jì)分布有所不同，地形越復(fù)雜等高線曲率值越大；黃土丘陵區(qū)的溝壑處等高線曲率最大，黃土塬區(qū)的現(xiàn)代溝坡次之，古代溝坡區(qū)較小，塬面和緩坡丘陵區(qū)則最小。等高線曲率可直接、客觀地反映等高線的復(fù)雜程度，在一定程度上也可以反映地形的復(fù)雜程度。

3.2 討 論

本文初步實(shí)現(xiàn)了對(duì)等高線曲率的計(jì)算，但是還有一些問(wèn)題有待深入研究，包括：（1）已經(jīng)實(shí)現(xiàn)的等高線曲率計(jì)算方法是基于單條等高線來(lái)計(jì)算曲率，沒(méi)有使用等高線拓?fù)潢P(guān)系信息，因而尚不能區(qū)分曲率的正負(fù)，也就不能區(qū)分地形的的凸凹變化。（2）等高線曲率與柵格數(shù)據(jù)計(jì)算的坡度、坡長(zhǎng)、平面曲率和剖面曲率等地形指標(biāo)之間，從概念、空間格局、統(tǒng)計(jì)分布等方面，具體的相同處和不同處，有什么量化關(guān)系。（3）比較平滑的等高線上，本文計(jì)算的曲率值存在的局部變異，有待通過(guò)對(duì)算法容錯(cuò)能力的提高，或者通過(guò)對(duì)等高線的修復(fù)處理來(lái)解決。

［1］ Moore I D，Grayson R B，Ladson A R.Digital terrain modelling：a review of hydrological，geomorphological，and biological applications［J］.Hydrological Processes，1991，5（1）：3－30.

［2］ 周啟鳴，劉學(xué)軍.數(shù)字地形分析［M］.北京：科學(xué)出版社，2006：57－63.

［3］ 程崢，姚志宏，張婧，等.南方紅壤農(nóng)作區(qū)數(shù)字流域建立和應(yīng)用［J］.水土保持研究，2012，19（3）：12－15.

［4］ Hutchinson M F.A new procedure for gridding elevation and stream line data with automatic removal of spurious pits［J］.Journal of Hydrology，1989，106（3）：211－232.

［5］ 楊勤科，Mcvicar T R，李領(lǐng)濤，等.ANUDEM：專(zhuān)業(yè)化數(shù)字高程模型插值算法及其特點(diǎn)［J］.干旱地區(qū)農(nóng)業(yè)研究，2006，24（3）：36－41.

［6］ Wilson J P.Digital terrain modeling［J］.Geomorphology，2012，137（1）：107－121.

［7］ 沈晶玉，史明昌.DEM網(wǎng)格尺寸對(duì)地形因子影響研究：以北京市延慶縣八達(dá)嶺小流域?yàn)槔跩］.水土保持研究，2006，13（5）：66－69.

［8］ 楊勤科，張彩霞，李領(lǐng)濤，等.基于信息含量分析法確定DEM分辨率的方法研究［J］.長(zhǎng)江科學(xué)院院報(bào)，2006，23（5）：21－23.

［9］ 王峰，王春梅.地形因子與DEM分辨率關(guān)系的初步研究：以蒙陰縣為例［J］.水土保持研究，2009，16（4）：225－229.

［10］ 陳其明.在自動(dòng)制圖中對(duì)曲線信息量的數(shù)學(xué)描述與壓縮［J］.測(cè)繪學(xué)報(bào)，1981，10（3）：239－254.

［11］ Evans I S.General geomorphometry，derivatives of altitude，and descriptive statistics［M］∥Chorley R J.Spatial Analysis in Geomorphology.London：Methuen＆Co.，1972.