地表面積計算方法研究

韓買俠，程傳錄，王　維

(國家測繪地理信息局大地測量數據處理中心, 陜西 西安 710054)

HAN Maixia，CHENG Chuanlu，WANG Wei

地表面積計算方法研究

韓買俠，程傳錄，王維

(國家測繪地理信息局大地測量數據處理中心, 陜西 西安 710054)

Research on the Calculation Method of Surface Area

HAN Maixia，CHENG Chuanlu，WANG Wei

摘要：通常面積是指地圖上的面積，而地理國情監測中面積是指地表面積，即實際面積。地表面積計算涉及坐標換算、投影、積分問題，計算復雜，本文研究了地表面積計算的理論，分析了影響地表面積計算精度的因素，確定了地表面積計算的方法步驟，并利用實際數據對該計算方法進行了驗證，取得了良好的效果。

引文格式： 韓買俠，程傳錄，王維. 地表面積計算方法研究[J].測繪通報，2015(9)：72-74.DOI:10.13474/j.cnki.11-2246.2015.0283

關鍵詞：地表面積；地理國情監測；圖斑面積；計算方法

中圖分類號：P237

文獻標識碼：B

文章編號：0494-0911(2015)09-0072-03

收稿日期：2014-08-01

作者簡介：韓買俠(1963—)，女，高級工程師，主要從事數據庫設計及4D產品坐標轉換方法研究。E-mail：hmxmaill@163.com

一、引言

目前我國正在實施的地理國情監測項目就是綜合利用全球導航衛星系統、航空航天遙感技術、地理信息系統技術等現代測繪技術，綜合各時期測繪成果檔案資料，對自然地理要素或地表人工設施要素等進行動態和定量化、空間化的監測，并統計分析其變化量、變化頻率、分布特征、地域差異、變化趨勢等，形成反映各類資源、環境、生態、經濟要素的空間分布及其發展變化規律的監測數據、地圖和研究報告等，從地理空間的角度客觀、綜合展示國情國力。地理國情監測中的一項重要內容是統計計算面狀要素的地表面積。由于三維地形是隨機起伏的，因此在這個面域內三維地形的表面也是隨機高低起伏的[1]，在計算面狀要素地表面積時，需要結合地形起伏的高程數據才能得到實際地表面積。

本文使用數字高程模型(DEM)和似大地水準面模型實現地表三維空間幾何恢復，將地表覆蓋的任意區域按DEM分辨率分為整格網與非整格網，對整格網和非整格網的地表面積計算方法分別進行了研究，分析了影響地表面積計算精度的因素。

二、基本原理

1. 地表面積定義

地理國情監測項目地形表面區域面積是指面要素沿地球表面累計的面積。在實際計算過程中以2000國家大地坐標系參考橢球為基礎，并考慮實際地形起伏下沿地球表面積累的植被覆蓋(面)、水系(面)、行政區劃單元(面)等面狀要素的面積。

2. 地表三維空間幾何模型恢復

DEM是通過有限的地形高程數據實現對地形曲面的數字化模擬(即地形表面形態的數字化表示)，反映了地表形態空間地理信息，DEM數據是采用二維平面坐標附加正常高高程實現的一種三維空間地理信息。在計算地形表面區域面積時，可以使用DEM正常高和似大地水準面模型中恢復地形表面的三維空間幾何模型，如圖1所示。

1) 將高斯平面坐標(x,y)反算為大地坐標(B,L)[2]。

2) 大地高(H)=正常高(h)+高程異常(ζ)[2]。

3) 將B、L、H換算為空間直角坐標X、Y、Z[2]。

圖1　地形在不同參考基準下的表達

3. 地形表面理論面積

地形表面區域面積是指面要素沿地形表面積累的面積。在空間直角坐標系下，沿地表區域對面積微元進行累加計算(如圖2所示)，即可獲得地形表面區域面積，即

式中,ds是沿地表區域的面積微元。

圖2　地表面積微元計算

三、地表區域面積計算方法

DEM數據有多種分辨率，不同分辨率DEM數據計算同一區域地形表面區域面積不同。理論上，DEM分辨率越高，計算的地形表面區域面積越接近真實值，但DEM分辨率越高計算地形表面區域面積效率越低。以1∶10 000比例尺DEM(分辨率為5 m)為基礎，計算任意區域地表面積方法如下：

1) 對于任意封閉區域(如圖3所示)，讀取范圍，再按DEM分辨率對封閉區域范圍外擴一個像素。

2) 根據封閉區域外擴范圍，計算覆蓋的DEM最大行列數(Row,Col)。

3) 根據范圍最小值及行列數，生成Row×Col個封閉矩形區域(如圖4所示)。

4) 根據生成封閉矩形4個角點坐標，讀取所在DEM正常高，根據地表三維空間幾何模型恢復方法，換算為空間直角坐標(X，Y，Z)。

圖3　任意封閉區域

圖4　封閉矩形區域

5) 計算圖3、圖4交集[3-5]得圖5，在圖5中包含大量的整格網(陰影)和非整格網(非陰影)。分別計算整格網和非整格網的面積，累加得到要素地表面積。

6) 整格網地表面積計算：基于整格網構建4個三角形(如圖6所示)，由海倫定理(即式(1))計算三角形地表面積[6]，三角形如圖7所示。

(1)

圖5　封閉區域與封閉矩形交集

圖6　整格網構建三角形

圖7　三角形

(2)

式中，(xi,yi,zi)、(xj,yj,zj)為兩點的空間直角坐標。

7) 計算非整格網地表面積：非整格網區域如圖8所示，計算地表面積方法就是將非整格網區域所有界址點構三角網[7-9]，計算三角網中三角形地表面積累計之和，即得到非整格網地表面積，如圖9所示。

圖8　非整格網區域

圖9　三角網

四、試驗及結果分析

筆者的研究以陜西省吳堡縣地理國情監測成果的地表覆蓋和DEM數據為基礎，進行了地表面積計算與分析。在本文中，以其中兩個比較復雜且具有代表性的地表覆蓋面域為例，計算其地表面積，并進行詳細的分析，面域如圖10所示。

1) 將1∶10 000分辨率為5 m的DEM數據內插為4 m、3 m、2 m分辨率數據。使用DEM正常高和似大地水準面模型恢復地形表面的三維空間幾何模型，計算圖10地表面積，結果見表1。

圖10　地表覆蓋面域

DEM分辨率/m面積(S1)/m2耗時/s面積(S2)/m2耗時/s51183651.6826173.281472121.3442103.15641168083.0308281.344466585.9800150.82831170208.2786410.938467486.9408217.37521184241.5754711.031472176.1035361.406

從表1可以看到，DEM分辨率不同，計算的地表面積值及效率也不同，分辨率越高，計算效率越低。分辨率2 m數據計算的地表面積與分辨率5 m數據計算的地表面積值相差比較小，前者計算效率比后者慢3～4倍。綜合計算精度與效率，建議使用5 m分辨率DEM。

2) 使用1∶10 000分辨率為5 m的DEM數據正常高和似大地水準面模型恢復地形表面的三維空間幾何模型，計算961個地表覆蓋面域的地表面積，見表2中采用大地高計算地表面積。

3) 使用1∶10 000分辨率為5 m的DEM數據正常高，忽略似大地水準面模型，計算961個地表覆蓋面域的地表面積，見表2中采用正常高計算地表面積。

表2　地表覆蓋面域地表面積　m2

從表2可以看到，似大地水準面模型對計算地表面積影響非常小，在實際計算地表面積時，對似大地水準面比較平緩的地區，可以不考慮似大地水準面模型對地表面積計算的影響。

4) DEM正常高高程對地表面積計算影響。一個整格網，4個點正常高相同，假設將其中一個格網點的正常高逐步抬高0.5 m，用似大地水準面模型恢復地形表面的三維空間幾何模型，計算其整格網地表面積，計算結果見表3。

表3　正常高變化對計算地表面積影響

從表3可以看到，正常高對計算地表面積影響比較大，也進一步證明表1 DEM分辨率不同，計算的地表面積值也不同。建議在實際計算時， 用1∶10 000分

辨率為5 m的DEM成果計算地表面積。

五、結束語

本文從地表面積的定義、三維空間恢復方法、面狀要素地表面積計算思路到計算的具體實現、計算結果分析等進行了較為全面的闡述，從實際計算結果分析了計算的效率和精度，并對計算中的DEM分辨率、高程、高程異常等因素的影響進行了分析和比較。經過實際驗證，文中提出的計算方案高效可行，可對地理國情監測中的地表面積計算與統計提供有價值的參考。

參考文獻：

[1]魏東,張秀程.基礎遞歸算法的三維地形面積計算方法研究[J].工程地質計算機應用，2007(3):5-7.

[2]孔祥元,郭際明,劉宗泉.大地測量學基礎[M].武漢：武漢大學出版社,2009.

[3]吳秀琴,張洪巖,李瑞改,等.ArcGIS 9地理信息系統應用與實踐[M].北京:清華大學出版,2008:231-233.

[4]樊建華,黃有群,劉嘉敏.帶空洞的多邊形求交集算法[J].沈陽工業大學學報,2001,23(5):429-431.

[5]楊維芳.兩個復雜多邊形求交的矢量算法[J].蘭州鐵道學院學報：自然科學版,2002,21(1):108-110.

[6]四川礦業學院數學教研組.數學手冊[M].北京：科學出版社,1978.

[7]武曉波,王世新,肖生春.Delaunay三角網的生成算法研究[J].測繪學報,1999,28(1):28-34.

[8]張渭軍.改進Delaunay生長算法研究[J].測繪科學,2011,36(4):58-59.

[9]余代俊，蒲朝旭，朱逍賢.一種Delaunay三角剖分的改進算法[J].測繪通報,2014(6):51-54.

[10]薛樹強, 黨亞民, 秘金鐘,等. 顧及非線性地形因子的地表面積計算[J]. 測繪學報， 2015, 44(3): 330-337.