基于星載雷達(dá)測高資料估計(jì)博斯騰湖水位-水量變化研究

吳紅波

(1.陜西理工大學(xué) 地理科學(xué)系， 陜西 漢中 723000； 2.西北大學(xué)陜西省地表系統(tǒng)與環(huán)境承載力重點(diǎn)實(shí)驗(yàn)室，陜西 西安 710127； 3.中國科學(xué)院青藏高原地球科學(xué)卓越創(chuàng)新中心， 北京 100101)

1 研究背景

湖泊作為區(qū)域陸地水循環(huán)和氣候變化的重要載體之一[1]，是大氣圈、冰凍圈、生物圈、巖石圈和陸地水圈相互作用的連接點(diǎn)[2]，是物質(zhì)循環(huán)、能量流動(dòng)和信息傳遞的重要場所。內(nèi)陸湖泊水量主要依靠大氣降水、地表徑流、冰雪融水和地下水補(bǔ)給。湖泊水量受氣候、地貌、土壤、植被等自然因素和人為活動(dòng)的共同作用，其變化過程存在一些確定性規(guī)律，同時(shí)也表現(xiàn)出強(qiáng)烈的隨機(jī)性。在氣候變化和人類活動(dòng)的共同干預(yù)下，湖泊水位和面積變化相繼發(fā)生了一系列波動(dòng)，人口增長和經(jīng)濟(jì)社會發(fā)展對內(nèi)陸河流徑流和湖泊水量的依賴程度愈加明顯[3]。因此，內(nèi)陸湖泊的變化不僅能夠揭示區(qū)域氣候變化特征，也有助于理解內(nèi)陸湖泊水循環(huán)過程、湖泊水量驅(qū)動(dòng)機(jī)理及水量平衡，應(yīng)對我國西北地區(qū)淡水資源短缺具有重要意義。

目前，湖泊水位、流量、面積等實(shí)時(shí)監(jiān)測資料的獲取，主要依靠地面定點(diǎn)觀測和水文調(diào)查方法，而水文調(diào)查方法時(shí)效性差、成本較高，地點(diǎn)觀測法較難反映湖泊水位、面積變化的時(shí)空特征。隨著遙感技術(shù)、空間定位、衛(wèi)星通信、地理信息系統(tǒng)等對地觀測技術(shù)的發(fā)展，星載遙感技術(shù)能為全球氣候變化提供全實(shí)時(shí)監(jiān)測的數(shù)據(jù)源。利用多光譜遙感技術(shù)在湖泊面積信息動(dòng)態(tài)監(jiān)測方面，已有相關(guān)報(bào)道[4]。劉瑞霞等[5]根據(jù)水體與地物光譜特性的反射率差異，利用NOAA/AVHRR資料定量反演出了博斯騰湖面積信息。玉素甫江·如素力等[6]采用多波段比值法對Landsat ETM影像中水體進(jìn)行分類識別與信息提取，并在博斯騰湖流域水體和非水體信息提取的總體精度達(dá)到99.29%。孫愛民等[7]利用長序列的Landsat影像提取1988-2014年博斯騰湖面積變化，并分析了湖水面積年際變化及空間變化趨勢。伊麗努爾·阿力甫江等[8]對博斯騰湖水位變化的自然要素與人為要素進(jìn)行驅(qū)動(dòng)力影響分析，認(rèn)為不同時(shí)期的氣候要素變化對水位波動(dòng)影響是有差異的。白淑英等[9]利用MODIS積雪產(chǎn)品和SSM/I雪深數(shù)據(jù)，系統(tǒng)分析了博斯騰湖流域雪蓋、雪深的時(shí)空變化規(guī)律及其對水量波動(dòng)影響。王濤等[10]利用博斯騰湖流域2001-2013年16天合成MODIS NDVI影像，分析了植被變化對水文響應(yīng)及時(shí)空差異。而且，對于湖泊水位變化估計(jì)及不確定性，星載雷達(dá)和激光雷達(dá)測高技術(shù)在垂直結(jié)構(gòu)和地表高度測量上的精度可達(dá)到分米級，可以滿足湖泊水位信息提取的精度要求，尤其是較大面積的湖泊水位估計(jì)具有較大潛力[11]。

近年來，國內(nèi)外研究人員利用TOPEX/Poseidon、ICESat、CryoSat-2、ENVISat、Jason-1&2等在湖泊的水位變化中進(jìn)行了一些相關(guān)研究，張國慶等[12]利用ICESat和Landsat數(shù)據(jù)對中國十大湖泊水量平衡進(jìn)行了估算，并結(jié)合湖水表面高程數(shù)據(jù)探討了湖泊水位變化趨勢。Hwang等[13]利用TOPEX/Poseidon測高資料嘗試構(gòu)建中國6個(gè)內(nèi)陸湖泊的水位時(shí)變序列，發(fā)現(xiàn)水位估計(jì)值與測量結(jié)果具有很好的一致性。但是，從已有的星載遙感技術(shù)發(fā)展、應(yīng)用和湖泊信息獲取手段來看，主要是針對湖泊面積、水位的遙感監(jiān)測與估計(jì)，而對湖泊水量變化參數(shù)估計(jì)未引起足夠重視，也鮮有相關(guān)報(bào)道。

為了檢驗(yàn)星載雷達(dá)測高數(shù)據(jù)在湖泊變化監(jiān)測的可行性，通過湖泊面積與水位、水位與水量變化的相關(guān)關(guān)系，聯(lián)合Landsat多光譜遙感數(shù)據(jù)，反演出博斯騰湖湖泊水位、面積、水量變化信息以及趨勢，為湖泊水文過程模型、水資源管理提供理論基礎(chǔ)和數(shù)據(jù)支持。文中基于ICESat、ENVISat、Jason-1&2、ERS-1&2測高數(shù)據(jù)和LandsatTM/ETM/OLI資料，首先，借助水體歸一化指數(shù)(Normalized Difference Water Index，NDWI)提取1990-2015年博斯騰湖湖泊水域面積；其次，利用雷達(dá)測高數(shù)據(jù)提取湖泊水面的高程信息和水位瞬時(shí)值，并結(jié)合大湖區(qū)附近水文站的水位、水量數(shù)據(jù)對ICESat、ENVISat、Jason-1&2、ERS-1&2的水位估計(jì)值進(jìn)行對比與驗(yàn)證；再次，根據(jù)湖泊面積-水位-水量的相關(guān)關(guān)系和遙感數(shù)據(jù)源不確定性誤差，構(gòu)建1990-2015年博斯騰湖水域面積、水位和水量波動(dòng)曲線；最后，分析博斯騰湖湖泊面積、水位和水量年際變化特征，旨在為博斯騰湖流域水資源合理配置、生態(tài)環(huán)境可持續(xù)發(fā)展和氣候變化研究提供技術(shù)支持。

2 數(shù)據(jù)與方法

2.1 研究區(qū)概況

博斯騰湖位于天山南麓焉耆盆地東南部，新疆維吾爾族自治區(qū)博湖縣境內(nèi)，屬于山間陷落湖，地理范圍45°56′～42°14′ N、86°40′～87°56′E。博斯騰湖包括大湖區(qū)、小湖群和湖濱濕地3部分。其中，博斯騰湖湖水補(bǔ)給主要來自開都河、黃水溝、清水河流域山區(qū)冰雪融水、大氣降水和地表徑流，湖泊出水補(bǔ)給孔雀河見圖1。流域?qū)儆谂瘻貛Ц苫哪畾夂?，受海拔高差影響，降水量、氣溫等差異較大。根據(jù)中國第二次冰川編目資料顯示[14]，流域內(nèi)冰川數(shù)目為700條，冰川面積達(dá)332.89 km2。

圖1 博斯騰湖流域及星載測高數(shù)據(jù)的分布

2.2 數(shù)據(jù)來源

(1)多光譜遙感數(shù)據(jù)。選取1990-2015年間湖泊水面無冰封、無云覆蓋(≤10%)的Landsat TM/ETM/OLI遙感影像，共297景(見表1)，其中，Path/Row為142/31的影像190景，Path/Row為143/31的影像107景。Landsat TM/ETM/OLI資料通過美國地質(zhì)勘探局(United States Geological Survey，USGS)Global Visualization Viewer (GloVis)查詢下載[15]。

(2)多源星載雷達(dá)測高數(shù)據(jù)。利用2003年2月至2009年10月ICESat-GLAS測高資料、2002年1月至2010年12月期間的Jason-1&2衛(wèi)星雷達(dá)高度計(jì)測高數(shù)據(jù)、2002年3月至2010年4月期間的ENVISat衛(wèi)星雷達(dá)測高資料，提取湖泊水面高度變化序列。其中，ICESat-GLAS測高資料中使用GLA01和GLA14產(chǎn)品[16]，數(shù)據(jù)版本分別為33和34，研究時(shí)段為2003年10月22日至2009年10月11日，湖泊大湖區(qū)內(nèi)GLAS光斑共計(jì)2 743個(gè)，見表2；GLAS數(shù)據(jù)可從美國冰雪數(shù)據(jù)中心(National Snow and Ice Data Center，NSIDC)網(wǎng)站http://nsidc.org/data/icesat/獲取。

選用1990-2015年ENVISat、ERS-1&2、Jason-1&2衛(wèi)星的GDR數(shù)據(jù)(Geophysical Data Record)，用于博斯騰湖瞬時(shí)水位估計(jì)，GDR數(shù)據(jù)通過法國圖盧茲大學(xué)LEGOS(Laboratoire d'Etudes en Géophy-sique etceanographie Spatiales)實(shí)驗(yàn)室(URL：http://www.legos.obs-mip.fr)/獲取，觀測時(shí)段為2002年10月1日- 2015年1月29日。

(3)水文觀測資料。博斯騰湖流域內(nèi)大山口、黃水溝、克爾古提、焉耆水文站的入湖年徑流量、年引水量數(shù)據(jù)，通過中國水文信息網(wǎng)(URL：http://www.hydroinfo.gov.cn/)獲取。寶浪蘇木、揚(yáng)水站位于湖泊水域周邊，其當(dāng)日水位觀數(shù)據(jù)用于與ICESat、Jason-1&2、ENVISat、ERS-1&2測高資料的水位估計(jì)值進(jìn)行驗(yàn)證與對比。塔什店水文站位于孔雀河上游，用于監(jiān)測湖泊流入孔雀河水量。

2.3 研究方法

2.3.1 湖泊水域面積信息提取

(1)選取1990年3月15日的Landsat TM遙感影像的湖泊界限和水域范圍，作為湖泊面積變化的參考。

(2)多源多時(shí)相Landsat遙感數(shù)據(jù)經(jīng)過幾何校正、地理配準(zhǔn)、輻射校正和大氣校正等處理后，可用于湖泊水域范圍識別和面積提取。

(3)利用歸一化水體指數(shù)RNDWI可突出水體信息，抑制地表土壤、植被信息的特性進(jìn)行閾值分割[17]。

(1)

式中:rgreen為Landsat TM/ETM/OLI影像的綠色波段反射率；rnir為近紅外波段的反射率。當(dāng)RNDWI值大于0.3時(shí)，該像元的地物類型為水體。

(4)當(dāng)確定博斯騰湖水域范圍后，統(tǒng)計(jì)出任一時(shí)期湖泊水域面積Alake,ti和面積誤差Se，湖泊水域面積誤差用公式(2)計(jì)算[18]：

(2)

式中：Se為任一時(shí)期面積誤差，km2；λ為波段空間分辨率，像元大小近似采用30 m × 30 m；εgeo為配準(zhǔn)誤差，m。通過增加地面控制點(diǎn)的數(shù)量，均勻分布控制點(diǎn)，使湖泊水域邊緣區(qū)的配準(zhǔn)誤差絕對值小于5 m。

(5)統(tǒng)計(jì)出任一ti時(shí)期Landsat影像中博斯騰湖水域范圍內(nèi)的像元數(shù)量，采用公式(3)估算湖泊水域面積Alake,ti：

(3)

式中:p為任一時(shí)期水體像元數(shù)量，p=1，2… ,n；alake,ti為水體像元面積，m2。

2.3.2 湖泊水位與水量變化估算

(1)湖泊水位與平均水位估計(jì)。星載雷達(dá)測高資料對湖泊瞬時(shí)水位Hlake的估計(jì)見公式(4)：

Hlake=Hsat-Crange+Cdelay+Cpresure+Cwet+

Cst+Cpt+e

(4)

式中:Hsat為代表的衛(wèi)星飛行高度，km；Crange為衛(wèi)星到地表距離，km；Cdelay為電離層傳播延遲校正，m；Cpresure為大氣氣壓變化所引起的信號延遲校正，m；Cwet為大氣濕度變化引起的信號延遲校正，m；Cst為地殼運(yùn)動(dòng)所引起的垂直高度修正值，m；Cpt為潮汐變化所引起的高程修正值，m；e為在計(jì)算過程中未考慮的不確定性誤差，m。

表1 1990-2015年博斯騰湖湖泊邊界提取所用的Landsat TM/ETM/OLI遙感數(shù)據(jù)

表2 2003-2009年博斯騰湖大湖區(qū)GLAS光斑數(shù)據(jù)匯總

(5)

(6)

(7)

(2)湖泊水量變化估計(jì)。根據(jù)內(nèi)流湖泊水量平衡方程，任一時(shí)段Tti+1-ti=ti+1-ti內(nèi)，湖泊水量變化值ΔVlake由湖泊水位和湖泊水域面積決定，估算關(guān)系式見公式(8)：

(8)

式中:Vti+1為ti+1時(shí)刻的湖泊容積，km3；Vti為ti時(shí)刻的湖泊容積，km3。

考慮湖泊面積和水位的遙感和估計(jì)偏差的共同影響，湖泊水量變化值ΔVlake簡化為：

(9)

(10)

2.3.3 湖泊水量變化的估計(jì)誤差 若忽略人工調(diào)水量和地下水交換等影響，內(nèi)陸湖泊水量變化的偏差主要由湖泊面積估計(jì)誤差和水位誤差共同決定[22]，湖泊水量變化的估計(jì)偏差δvlake計(jì)算式為：

(11)

鑒于湖盆的形態(tài)及其變化對湖泊容積和水量變化估計(jì)產(chǎn)生不確定誤差，對于同一湖泊的水量誤差，湖盆地形誤差ebed可描述為

ebed=e1+e2+e3+e4+e5

(12)

式中：e1、e2、e3、e4、e5分別為湖岸、沿岸帶、湖岸邊淺灘、水下斜坡、湖盆底地形變化所引起的水量估計(jì)偏差。

3 結(jié)果與分析

3.1 水位估計(jì)值與當(dāng)日水位觀測值對比

為了檢驗(yàn)ENVISat&ERS、Jason-1&2和ICESat-GLAS的瞬時(shí)水位估計(jì)值的準(zhǔn)確性和可比性，既要對星下點(diǎn)做同一投影坐標(biāo)系和高程參考，也需要對系統(tǒng)誤差做一元線性回歸擬合和地球物理修正，使湖泊水面高程值具有可比性。由圖2可知，ENVISat&ERS、Jason-1&2和ICESat-GLAS的當(dāng)日水位估計(jì)值與附近揚(yáng)水站、寶浪蘇木水文站的湖面水位觀測值的簡單相關(guān)系數(shù)分別為0.95、0.97和0.98，絕對誤差分別小于0.21、0.18和0.15 m。隨著ICESat衛(wèi)星運(yùn)行時(shí)間的推移和儀器特性下降，星下點(diǎn)的高程誤差有節(jié)律性增加，此類系統(tǒng)性誤差可通過線性回歸平移模型進(jìn)行修正。由于ENVISat和ERS衛(wèi)星的測高雷達(dá)的測量精度限制，二者的當(dāng)日水位估計(jì)值與當(dāng)日觀測值存在一定的高度偏差，可通過增加地表控制點(diǎn)數(shù)量和布點(diǎn)區(qū)域，降低水位估計(jì)誤差。雖然ENVISat&ERS、Jason-1&2和ICESat-GLAS衛(wèi)星航跡對應(yīng)的星下點(diǎn)和GLAS光斑與揚(yáng)水站、寶浪蘇木站存在一定距離，但是當(dāng)日水位估計(jì)值與觀測值存在較強(qiáng)的相關(guān)性和一致性，能滿足博斯騰湖湖泊水位估計(jì)的要求。

圖2 湖泊水位觀測值與GLAS數(shù)據(jù)、Jason-1&2和ENVISat資料的水位估計(jì)值

3.2 湖泊水位與面積關(guān)系

湖泊水位是連接湖泊面積與水量的重要參數(shù)之一，同時(shí)，湖泊面積與水位、水位與水量存在較好的相關(guān)關(guān)系。為此，文中根據(jù)測高衛(wèi)星過境時(shí)間和當(dāng)日水位觀測值，隨機(jī)選取2002-2015年60景無云覆蓋、無冰封的Landsat TM/ETM影像提取湖泊水域面積，根據(jù)當(dāng)日湖泊水位觀測值和湖泊水域面積，利用冪函數(shù)對湖泊水域面積與揚(yáng)水站當(dāng)日水位觀測值進(jìn)行擬合(見圖3(a))，復(fù)相關(guān)系數(shù)R2為0.9201，通過0.01的顯著水平檢驗(yàn)。利用一元線性回歸關(guān)系式將當(dāng)日水位觀測值與湖泊水量變化進(jìn)行擬合(見圖3(b))，復(fù)相關(guān)系數(shù)R2為0.9058，通過0.01的顯著水平檢驗(yàn)。根據(jù)1990-2015年Landsat衛(wèi)星的297景中提取的湖面水域面積時(shí)變序列，將面積-水位擬合關(guān)系式、水位-湖泊水量變化擬合關(guān)系式，用于構(gòu)建1990-2015年湖泊水量變化波動(dòng)曲線。

圖3 湖泊當(dāng)日水位值與湖泊面積、水量變化的擬合曲線

3.3 湖泊水位-面積-水量波動(dòng)變化

1990-2015年博斯騰湖湖泊面積、水位、水量變化曲線見圖4。由圖4(a)可知，1990-2015年博斯騰湖水域面積波動(dòng)較大，1990-2000年期間，博斯騰湖水域面積呈現(xiàn)增長趨勢，1990年湖泊面積為(937.8±3.5)km2，到2000年增長到1 072.9 km2，以年均(13.5±3.0)km2的速度增長。2002年年均水位達(dá)到1956年以來湖泊年均水位的最大值(1048.6±0.25)m。與此同時(shí)，湖泊最大水域面積為(1112.4±3.5)km2。2002-2015年期間，博斯騰湖水域面積波變化呈現(xiàn)縮減趨勢。到2015年，水域面積為(905.60±3.5)km2，比2005年面積減少了(63.9±3.5)km2，比1990年減少了(32.20±3.5)km2。

圖4 1990-2015年博斯騰湖湖泊面積、水位、水量變化曲線

由圖4(b)、4(c)進(jìn)一步對博斯騰湖年內(nèi)逐月水位、水量變化情況分析如下：

湖泊水位年內(nèi)變化曲線呈現(xiàn)出雙峰特征，3-4月份期間水位逐漸升高，出現(xiàn)第一個(gè)峰值；5-6月份緩慢下降，7-8月份形成第二個(gè)峰值；隨后湖泊水位開始回落，11月份或者12月份達(dá)到最低值，湖面會出現(xiàn)封凍現(xiàn)象。由于冬季降雪、湖面封凍，雷達(dá)測高數(shù)據(jù)對湖面冰雪的回波信號異常，提取的湖面泊水高程值存在較大偏差，一般在湖泊水位年際變化計(jì)算中予以剔除。

博斯騰湖年內(nèi)的水量波動(dòng)較大，湖泊水量波動(dòng)受開都河、清水河上游山區(qū)降水變化、冰雪融水的作用，6-8月份氣溫上升，上游山區(qū)冰雪消融；6-8月份也是降水集中的時(shí)段，地表徑流增加了入湖水量，因此，7-8月份湖泊水容量達(dá)到當(dāng)年最大值。

3.4 年際變化分析

1990-2015年博斯騰湖年均水位、面積、水量變化值的年際變化見圖5。由圖5可知，1990-2015年博斯騰湖湖泊面積、年均水位和水量變化整體上波動(dòng)較大，1990-2002年呈現(xiàn)出上升過程，2002-2015年呈現(xiàn)出下降過程；其中，1990-2002年博斯騰湖年均水位上升了(3.24±0.20)m，平均每年上升0.27 m。2015年博斯騰湖湖泊水位較2002年水位下降了(4.04±0.23)m，平均每年下降0.31 m，而且湖泊面積減少了約(32.20±3.5)km2。

圖5 1990-2015年博斯騰湖年均水位、面積、水量變化值的年際變化

與1990年湖泊年均水位相比，2015年博斯騰湖的水位和、水量呈現(xiàn)下降趨勢(見圖5(c))，2015年博斯騰湖湖泊年均水位下降了(0.81±0.19)m；2015年湖泊水量比1990年湖泊水量減少了(9.49±0.22)×108m3。2002-2015年博斯騰湖水位和水量均呈減少趨勢，雖然在2010年湖泊水量略有增加，但2005年之后，湖泊水量出現(xiàn)急速減少。博斯騰湖作為中國西北干旱區(qū)內(nèi)流湖，其水量波動(dòng)不僅受氣候變化影響，也在一定程度上收到人類活動(dòng)的干擾。尤其是在氣候變化背景下，圍湖造田、不合理開發(fā)淡水資源、水利工程調(diào)蓄等人類活動(dòng)不僅在一定程度上造成入湖徑流減少、湖泊水位下降，而且也會加劇湖泊面積萎縮和湖區(qū)生態(tài)環(huán)境破壞[23]。

4 結(jié) 論

(1)在考慮衛(wèi)星航跡的星下點(diǎn)投影坐標(biāo)系、地理配準(zhǔn)、地球物理修正后，ICESat-GLAS、ENVISat&ERS、Jason-1&2測高數(shù)據(jù)在博斯騰湖湖泊水位估計(jì)中具有一定可比性，ICESat-GLAS、ENVISat&ERS、Jason-1&2的當(dāng)日水位估計(jì)值與附近揚(yáng)水站、寶浪蘇木水文站的湖面水位觀測值絕對誤差分別小于0.21、0.18、0.15 m，而且具有較強(qiáng)的相關(guān)性和一致性，能滿足博斯騰湖湖泊水位估計(jì)的要求。

(2)湖泊水位作為湖泊面積、水量之間關(guān)聯(lián)的重要參數(shù)，同時(shí)，博斯騰湖湖泊面積與水位、水位與水量存在顯著的相關(guān)關(guān)系，復(fù)相關(guān)系數(shù)R2分別為0.9201、0.9058。

(3)1990-2015期間年博斯騰湖湖泊面積、年均水位和水量變化整體上波動(dòng)較大，1990-2002年水位、水量處于增加過程，2002-2015年水位、水量處于下降趨勢。與1990年湖泊水位相比，2015年湖泊的年均水位下降了(0.81±0.19)m，2015年湖泊水量比1990年減少了(9.49±0.22)×108m3，湖泊面積縮減了(32.20±3.5)km2。

此外，湖泊水量波動(dòng)直接反映出人類活動(dòng)和氣候變化共同作用的結(jié)果，湖泊水位是地表水體水量波動(dòng)的遙感反演的重要參數(shù)，而且星載雷達(dá)測高技術(shù)的工作模式和數(shù)據(jù)采集方式正在打破地表遙感反演的時(shí)空限制，使測量精度和數(shù)據(jù)覆蓋能夠滿足面積較小的湖泊、河流、水庫等地表水體的水位、水量監(jiān)測需求。