基于地形參數的湖泊水儲量估算方法
——以納木錯為例

馬小奇，盧善龍，馬 津，朱立平

(1.中國科學院遙感與數字地球研究所遙感科學國家重點實驗室，數字地球院重點實驗室，北京 100101；2.中國地質大學(北京)地球科學與資源學院，北京 100083；3.山東農業大學信息科學與工程學院，泰安 271018；4.中國科學院青藏高原研究所，北京 100101)

0 引言

湖泊的水儲量及其隨時間的變化影響著湖泊生態系統的物理、化學及生物過程，可作為降水和蒸發之間的平衡及地表水和地下水系統之間相互作用的指示性指標[1-2]。準確監測和估算湖泊水儲量對開展湖泊生態學、水文學研究和水資源調度管理具有重要意義。湖泊水儲量的計算通常利用水下地形數據，基于“3S”(遙感(remote sensing，RS)、全球定位系統(global positioning system，GPS)和地理信息系統(geographic information system，GIS))技術進行估算。Gleason等[3]利用現場測量GPS數據構建了面積—體積預測模型；Lane等[4]利用激光雷達數據和ArcGIS軟件中的不規則三角網(triangulated irregular network,TIN)多邊形體積模型計算了美國佛羅里達州中北部獨立湖泊的水量；Zhang等[5]利用多年遙感影像、納木錯水下地形和GIS技術，給出了1976—2009年間納木錯水儲量年度變化。上述方法可以準確地計算湖泊水儲量及其變化，但對于無資料地區而言，受自然環境、資金和人力等因素的影響，難以獲取水下地形和水位等實測數據。截止目前，全球范圍內還有很大一部分湖泊水體沒有實測數據用于其水儲量的計算。

湖泊水面以下的湖盆和水面以上的湖盆環境是由相同的地質演化過程(如冰川作用、河流沖刷沉積等)形成的，湖泊水下地形的起伏應該是水面以上周圍地形的一個延續，具有相似特征。因此，湖泊周圍地形高度的變化特征可能在其表面以下具有類似的變化特征[6]。Sobek等[7]將瑞典湖泊數據集與GIS導出的地理和地形參數相結合，基于地形圖上獲取的地形參數預測湖泊的體積和深度，發現湖泊面積和湖泊邊界以外50 m寬區域的最大坡度是湖泊體積和深度最重要的預測因子；Messager等[8]開發了地理統計模型并基于湖泊周圍的地形信息計算了全球142萬個面積小于500 km2的天然湖泊的平均深度，通過平均深度與湖面面積的乘積來計算湖泊體積。

青藏高原分布著海拔高、數量多、面積大的高原湖泊群。至今這些湖泊大多處于自然狀態，受人類活動影響較小；同時由于高寒缺氧，自然條件差，人員設備難以到達，許多湖泊水位與水儲量之間的關系仍難以確定。近年來雖然有很多關于高原湖泊水儲量變化的研究，但鮮有關于湖泊凈水儲量估算的研究。Liao等[9]利用1970年的地形圖和Landsat MSS，TM和ETM+衛星影像，提取高原主要湖泊的表面積，分析青藏高原40 a間湖面面積變化與氣候變化之間的關系；Yang等[10]和朱立平等[11]基于航天飛機雷達地形測繪任務(shuttle Radar topography mission，STRM)數字高程模型(digital elevation model，DEM)和Landsat影像，估算了1976—2013年間青藏高原大于50 km2湖泊的水量變化；Qiao等[12]和Song等[13]將光學遙感影像、衛星測高數據和重力測量數據集成在一起，研究分析發現青藏高原的大部分湖泊水儲量呈現出明顯的上升趨勢。本文以納木錯為例，介紹一種利用湖泊水面以上地形和湖泊水面面積數據構建面積—體積模型的湖泊水儲量估算方法。

1 研究區概況及數據源

納木錯是青藏高原第三大的湖泊，地理范圍在E90°16′～91°03′，E30°30′～30°55′之間[14]，位于西藏自治區那曲地區當雄縣，2000年經STRM DEM測得湖面海拔為4 724 m，水域面積為1 944.30 km2，最大深度為98 m[15]。2005年以前，納木錯流域內沒有氣象觀測站或水文站；2005年后中國科學院青藏高原研究所設立了納木錯觀測研究站，逐步開展了一些常規參數的監測。研究人員在2005—2007年間進行了3次綜合調查，得出了納木錯的水深資料[15]。

本研究使用的DEM數據為90 m空間分辨率的STRM DEM數據,由美國國家航空航天局噴氣推進實驗室與美國國家影像與制圖局聯合測量獲得[16-18]，來源于中國科學院計算機網絡信息中心國際科學數據鏡像網站(http://www.gscloud.cn)；納木錯湖泊水面數據來源于盧善龍等[19]以MODIS MOD09Q1為數據源提取的2000—2012年間青藏高原每8 d的湖泊水面數據集；驗證數據參考Wang等[15]發表的納木錯水下DEM,該數據是基于2005—2007年間湖泊測深數據(共305 721個測量點)，通過等深線圖繪制及水下TIN模型構建獲得。

2 方法原理

假定湖泊水下為規則湖盆，如錐形湖盆(水面與水深呈線性變化)或杯型湖盆(水面與水深呈非線性變化)，則湖泊水下地形可以認為是水上地形的延續，具有相似性特征。湖泊水面以上一定范圍內的地形特征參量(如地形高差、平均高程和坡度等)可以作為水下地形的表征指標。因此，通過分析和構建湖泊水面以上高程、平面面積、體積增量之間的函數關系能推導出水面以下水量與對應水面面積之間的數值關系模型(圖1)。

圖1 基于湖盆地形相似特征的湖泊水儲量估算方法示意圖Fig.1 Lake water volume estimation model based on lake basin similarity

計算步驟如下：

1)以單個湖盆DEM數據為基礎，計算湖面H0以上等間距的不同高程平面的面積Si(i=1,2,…，n)，建立面積Si與高程Hi的函數關系。

2)計算湖面以上等間距的不同高程平面至湖面間的體積Vi(i=1,2,…，n)，利用

ΔVi=Vi-Vi-1，i=1,2,…，n，

(1)

計算單位高程增加對應的體積增量ΔVi，i=1,2,…，n，建立面積Si與體積增量ΔVi的函數關系。

(2)

3 水儲量估算模型構建

構建過程包括湖泊水面提取，數據統計，回歸分析，高程—面積函數關系、面積—體積增量函數關系和面積—體積模型構建等步驟(圖2)。

圖2 本文方法流程Fig.2 Process of the proposed method

3.1 水面以上單位高程變化對應體積計算

以納木錯湖盆范圍SRTM DEM作為基礎研究數據，對其進行等高線分析，發現納木錯湖面高程為4 724 m。提取湖盆地中4 724 m等高線構建納木錯湖泊水面，面積為1 944.30 km2(圖3)。

圖3 納木錯湖面及周圍等高線圖Fig.3 Contour of 4 724 m and 4 748 m and of Nam Co Lake and the maximum water surface in 2008

青藏高原的現代湖泊絕大多數為封閉的咸水湖和鹽湖,在其周圍廣泛分布著湖泊收縮時殘留的湖濱相砂礫石和湖蝕崖坎構成的古湖岸線或粉細砂黏土組成的湖積平原和階地 ,依此可以重建古湖面的高度。納木錯湖泊擴張時形成的最高古湖岸線重建湖面的高度為4 748 m[20-21]，以高程4 724 m對應的水平面為起始參考平面，利用表面體積功能求取高程4 724～4 748 m間每上升1 m所在平面的面積及該平面至起始參考面的體積Vi，這部分體積即是水位上升至該高程時水儲量的變化量。利用式(1)計算相鄰2個參考面的體積增量ΔVi，具體如表1所示，其中將4 724 m處體積視為0.01億m3。

表1 納木錯湖水面以上單位高程變化對應的體積及體積增量Tab.1 Calculated volumes in different elevation related to 4 724 m of Nam Co Lake

3.2 周圍地形高程—面積—體積增量間的函數關系

統計分析表1中不同湖面高程對應水平面的面積，結果表明，兩者呈強線性關系(R2=0.998 7)(圖4)，其表達式為

S=19.181H-88 672。

(3)

圖4 不同湖面高程與平面面積之間的函數關系Fig.4 Functional relationship between elevation and surface area in different location above 4 724 m

采用同樣的方法計算得到不同湖面高程對應的面積與體積增量之間的函數關系，兩者之間同樣呈強線性關系(R2=0.999 9)(圖5),其表達式為

ΔV=0.010 1S-0.278 7。

(4)

圖5 不同高程對應的平面面積與體積增量之間的函數關系Fig.5 Functional relationship between surface area and volume increment in different elevation

3.3 面積—體積模型的構建

V=0.000 3S2-0.009 5S+0.083 1。

(5)

圖6 湖面面積與湖水體積的函數關系Fig.6 Relationship between area and volume

4 精度評價與誤差分析

4.1 模型精度評價

以4 724 m等高線提取的納木錯水面面積(1 944.30 km2)為輸入，利用式(5)估算求得納木錯湖水儲量為1 115.70億m3。

參考2005—2007年間納木錯湖測深數據建立的水下地形圖(圖7)，通過ArcGIS10.1軟件功能性表面工具集的表面體積功能計算得到納木錯的水儲量為1 019.50億m3。以水下地形圖求得水儲量為實測值，可知本文模型估算的誤差(實測值與模擬值的差值)為-96.20億m3，相對誤差為9.40%。

圖7 納木錯水下地形圖Fig.7 Underwater topography of Nam Co Lake

4.2 模型誤差分析

為了進一步分析模型模擬結果隨水深變化的情況，計算并分析了水下95 m至水面范圍內不同水面面積及對應體積增量的模型估算結果與實測結果之間的差異，如圖8—9所示。水深在0～20 m范圍內，模型對水面面積與體積增量的估計值比實測值大，且隨著深度的增加差異逐漸減少；水深在20～60 m范圍內，模型對水面面積和體積增量的估計均比實測小,且在40 m水深左右隨著深度的增加差異逐漸減少，之后呈緩慢增加趨勢；水深在60～80 m范圍內，不同深度對應的面積誤差與體積增量誤差出現突增，且隨著深度的增加誤差逐漸增大；水深在80～95 m范圍內，不同深度的面積誤差與體積增量誤差出現回落，且隨著深度的增加誤差逐漸減少。

圖8 不同水深與對應水平面面積誤差Fig.8 Errors between model calculated and TIN derived surface areas in different water depths

圖9 不同水深與對應單位深度范圍內的體積增量誤差Fig.9 Errors between the model calculated and TIN derived volumes in different water depths and related changed volumes

整體而言，當水深在60 m以內時，湖泊地形呈單一變化趨勢，模型能較好地反映地形特征，模型估算的體積誤差較小；當水深超過60 m后，地形變化逐漸復雜，特別是在湖泊西北和東北3個水下地形凸起區域(圖7黑框區域)，地形變化顯著，單一趨勢變化的模型無法正確刻畫這種起伏變化的水下地形特征，模型估算結果誤差急劇增加。而當深度超過83 m時，因實際地形又趨于單一，模型結果誤差再次減小。

5 結論

1)本研究基于湖盆地形的相似特性，構建了一種基于湖泊水面以上地形高程、面積和體積增量數值關系的湖泊水儲量估算方法。該方法具有不依賴于實測水深測量數據的特點，為無資料地區水上和水下的地形特征參量基本一致湖泊水儲量計算提供了一種簡單而高效實用的技術途徑。

2)利用該方法基于納木錯湖盆DEM求得水面面積為1 944.30 km2時水儲量為1 115.70億m3，該結果與基于實測水深數據建立的納木錯水下地形DEM計算的水儲量1 019.50億m3相比，絕對誤差為96.20億m3，相對誤差為9.40%。

3)本文構建的納木錯水儲量估算模型，可用于納木錯多年水儲量監測，在水下地形和水位數據匱乏的無資料地區，利用該方法構建的水儲量估算模型只要獲得湖泊水面面積即可估算出水儲量，對于水量變化監測具有實用意義。

4)對于水下和水下地形相似性較弱的湖泊，因模型模擬方式難以準確刻畫其復雜的水下地形特征，本文方法的應用仍具有一定的局限性。但對于這類湖泊利用本文構建的水儲量模型只需有可靠的水面面積變化數據，便可求得精度較高的水儲量變化數據。

5)青藏高原面積大于1 km2的湖泊有1 055個[22]，但目前只有30多個湖泊有實測水下地形數據。未來的研究中，將以這些有實測數據的湖泊為參考，針對不同面積大小、不同類型的湖泊研究和構建相應的湖泊水儲量估算模型，并最終形成可用于整個高原區湖泊水儲量估算的技術方法體系。