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

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

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

0 引言

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

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

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

1 研究區(qū)概況及數(shù)據(jù)源

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

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

2 方法原理

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

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

計算步驟如下：

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

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

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

(1)

計算單位高程增加對應(yīng)的體積增量ΔVi，i=1,2,…，n，建立面積Si與體積增量ΔVi的函數(shù)關(guān)系。

(2)

3 水儲量估算模型構(gòu)建

構(gòu)建過程包括湖泊水面提取，數(shù)據(jù)統(tǒng)計，回歸分析，高程—面積函數(shù)關(guān)系、面積—體積增量函數(shù)關(guān)系和面積—體積模型構(gòu)建等步驟(圖2)。

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

3.1 水面以上單位高程變化對應(yīng)體積計算

以納木錯湖盆范圍SRTM DEM作為基礎(chǔ)研究數(shù)據(jù)，對其進(jìn)行等高線分析，發(fā)現(xiàn)納木錯湖面高程為4 724 m。提取湖盆地中4 724 m等高線構(gòu)建納木錯湖泊水面，面積為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

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

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

3.2 周圍地形高程—面積—體積增量間的函數(shù)關(guān)系

統(tǒng)計分析表1中不同湖面高程對應(yīng)水平面的面積，結(jié)果表明，兩者呈強線性關(guān)系(R2=0.998 7)(圖4)，其表達(dá)式為

S=19.181H-88 672。

(3)

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

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

ΔV=0.010 1S-0.278 7。

(4)

圖5 不同高程對應(yīng)的平面面積與體積增量之間的函數(shù)關(guān)系Fig.5 Functional relationship between surface area and volume increment in different elevation

3.3 面積—體積模型的構(gòu)建

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

(5)

圖6 湖面面積與湖水體積的函數(shù)關(guān)系Fig.6 Relationship between area and volume

4 精度評價與誤差分析

4.1 模型精度評價

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

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

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

4.2 模型誤差分析

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

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

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

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

5 結(jié)論

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

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

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

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

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