王夢(mèng)琦 張 文 孟令奎

1武漢大學(xué)遙感信息工程學(xué)院，湖北 武漢，430079

隨著遙感技術(shù)不斷發(fā)展，遙感影像的數(shù)據(jù)質(zhì)量和分辨率不斷提升，遙感技術(shù)也逐漸成為水庫(kù)時(shí)空變化監(jiān)測(cè)的重要手段。孫建蕓等［1］利用GF-1影像提取丹江口水庫(kù)面積進(jìn)行分析；馬艷敏等［2］基于MODIS數(shù)據(jù)對(duì)2000—2017年月亮泡水庫(kù)、向海水庫(kù)和查干湖面積進(jìn)行時(shí)空變化分析；龔文峰等［3］利用Landsat8 OLI數(shù)據(jù)提取磨盤(pán)山水庫(kù)面積進(jìn)行分析；郝小翠等［4］利用MODIS數(shù)據(jù)分析甘肅省主要水庫(kù)的面積變化特征。但目前對(duì)北京密云水庫(kù)和官?gòu)d水庫(kù)時(shí)間序列的面積變化分析較少，本文以北京地區(qū)密云水庫(kù)和官?gòu)d水庫(kù)為研究對(duì)象，采用GF-1數(shù)據(jù)對(duì)水庫(kù)面積進(jìn)行提取和變化分析，結(jié)合氣溫降雨等因素分析南水北調(diào)工程對(duì)密云水庫(kù)和官?gòu)d水庫(kù)的變化影響。

1 研究區(qū)概況及數(shù)據(jù)預(yù)處理

1.1 研究區(qū)概況

密云水庫(kù)位于北京密云區(qū)，是北京最大的唯一水源供應(yīng)地，位于潮河、白河中游偏下，主要水流來(lái)源為潮白河［5］，負(fù)責(zé)供水、防洪和灌溉。官?gòu)d水庫(kù)橫跨北京和河北，位于河北省張家口市懷來(lái)縣和北京市延慶縣界內(nèi)，官?gòu)d水庫(kù)是我國(guó)建國(guó)后所建設(shè)的第一座大型水庫(kù)，主要水流來(lái)源為河北懷來(lái)永定河，在防洪、灌溉、發(fā)電等方面發(fā)揮了巨大的作用。密云水庫(kù)和官?gòu)d水庫(kù)在北京地區(qū)的分布如圖1所示。

圖1 密云、官?gòu)d水庫(kù)分布圖Fig.1 The Distribution of Miyun Reservoir and Guanting Reservoir

1.2 數(shù)據(jù)選取及預(yù)處理

本文所采用的數(shù)據(jù)源是GF-1衛(wèi)星的寬幅相機(jī)WFV所獲取的1A級(jí)多光譜遙感影像（影像來(lái)源：中國(guó)資源衛(wèi)星應(yīng)用中心http：//www.cresda.com/cn/）。選擇2014年1月至2019年12月每月15日左右水庫(kù)處無(wú)云或云量較少的數(shù)據(jù)，最終選取了130景GF-1衛(wèi)星的WFV數(shù)據(jù)。其中密云水庫(kù)2018年6月影像缺失，官?gòu)d水庫(kù)2014年9月影像缺失。對(duì)獲取的遙感影像進(jìn)行波段合成、輻射定標(biāo)、大氣校正和正射校正，本文以2017年1月影像為基準(zhǔn)對(duì)其余時(shí)間序列的影像進(jìn)行配準(zhǔn)。

2 研究方法

2.1 水體提取方法

由于密云水庫(kù)是丘陵型水庫(kù)，官?gòu)d水庫(kù)為山谷型水庫(kù)，一年中3月至10月水庫(kù)周?chē)恢脖桓采w。為了更好的區(qū)分水庫(kù)和植被，基于歸一化水體指數(shù)（normalized difference water index，NDWI）進(jìn)行水體提取，NDWI可以有效消除植被信息，突出水體，NDWI計(jì)算公式為［6］：

式中，Green表示綠光波段；NIR表示近紅外波段。

一般的，對(duì)于NDWI影像進(jìn)行水體和非水體分割的閾值選取為0，但由于本文所選的水庫(kù)分布位置不同且時(shí)間跨度較大，選擇NDWI=0作為水體和非水體的分割閾值得到的水庫(kù)提取效果較差，因此對(duì)不同影像選取不同閾值進(jìn)行水庫(kù)提取，本文選取的水庫(kù)分布于山脈間，水庫(kù)周?chē)脖惠^多，在NDWI影像上，水體的值較大主要集中在NDWI≥0處，而背景值較小主要集中在NDWI＜0處，對(duì)NDWI影像進(jìn)行像頻值統(tǒng)計(jì)常常會(huì)出現(xiàn)雙峰現(xiàn)象，在NDWI值較大處集中的峰可以認(rèn)為是水體，而在NDWI值較小處集中的峰可以認(rèn)為是背景，一般兩個(gè)峰值間的像頻值最小的位置可以認(rèn)為是最佳分割閾值，用來(lái)區(qū)分水體和非水體，為了實(shí)現(xiàn)自動(dòng)化閾值選取，需要找到影像的兩個(gè)像頻值最大的峰值，兩個(gè)峰值之間的像頻值最小處的NDWI值即是影像最佳分割閾值［7］，實(shí)現(xiàn)水體自動(dòng)化提取。

2.2 結(jié)冰水面提取方法

京津冀地區(qū)冬季平均氣溫僅為0～8℃，在1月、2月和12月會(huì)出現(xiàn)部分水庫(kù)結(jié)冰現(xiàn)象，結(jié)冰水面和水體的反射特性差異較大，使用NDWI無(wú)法提取出結(jié)冰水面。結(jié)冰水面和積雪具有相似的反射特性，但傳統(tǒng)的積雪指數(shù)需要短波紅外波段，GF-1影像不具備短波紅外波段，無(wú)法使用積雪指數(shù)進(jìn)行結(jié)冰水面的提取。在遙感影像中結(jié)冰水面在各個(gè)波段的反射率均較高，本文采用四波段相加之和自適應(yīng)閾值法進(jìn)行影像分割，利用結(jié)冰水面在4個(gè)波段均有較高反射率的特性區(qū)分冰面和背景，實(shí)現(xiàn)對(duì)水庫(kù)中結(jié)冰處反射率較高、在影像上顏色較淺的冰面進(jìn)行提?。?，9］。

2.3 水庫(kù)面積變化指標(biāo)

對(duì)水庫(kù)面積變化采用定量分析和定性分析，定量分析采用水庫(kù)面積變化幅度指標(biāo)和水庫(kù)面積變化強(qiáng)度指標(biāo)。水庫(kù)變化幅度指標(biāo)可以反映水庫(kù)在一定時(shí)間范圍內(nèi)的面積變化幅度，有正變化和負(fù)變化之分，水庫(kù)面積變化幅度指標(biāo)的表達(dá)式為［10］：

式中，L1表示在研究時(shí)間范圍內(nèi)水庫(kù)面積變化總幅度；Ua表示研究期初水庫(kù)面積；Ub表示研究期末水庫(kù)面積；T表示研究時(shí)間間隔年數(shù)；L2表示在研究時(shí)間范圍內(nèi)水庫(kù)面積年變化幅度。當(dāng)L1＞0時(shí)，水庫(kù)面積增加，水庫(kù)面積為正變化，L1＜0時(shí)，水庫(kù)面積減小，水庫(kù)面積為負(fù)變化。

水庫(kù)面積變化強(qiáng)度指標(biāo)可以描述水庫(kù)在不同時(shí)期的相對(duì)面積變化，使用該指標(biāo)可以定量描述水庫(kù)的面積變化強(qiáng)弱［11，12］，其表達(dá)式為：

式中，C表示水庫(kù)面積變化強(qiáng)度；ΔA ba表示研究期a～b時(shí)間內(nèi)水庫(kù)的面積變化；A a表示a年水庫(kù)的面積；Δt表示研究時(shí)間間隔年數(shù)。

3 結(jié)果與分析

3.1 密云水庫(kù)面積變化分析

密云水庫(kù)在2014年2月、2017年1月和2019年2月影像出現(xiàn)結(jié)冰現(xiàn)象，利用NDWI自適應(yīng)閾值方法無(wú)法提取出結(jié)冰水面，以2014年2月影像為例，如圖2所示，而采用四波段相加之和的自適應(yīng)閾值法可以很好地提取出結(jié)冰水面范圍，結(jié)合NDWI的提取結(jié)果實(shí)現(xiàn)對(duì)結(jié)冰水庫(kù)的完整提取。

圖2 密云水庫(kù)2014年2月結(jié)冰水面提取結(jié)果Fig.2 Extraction Results of Miyun Reservoir Frozen Water in Feb.2014

對(duì)2014—2019年各年內(nèi)的水庫(kù)面積變化進(jìn)行統(tǒng)計(jì)得到的結(jié)果如圖3所示，2014—2019年密云水庫(kù)面積總體變化趨勢(shì)如圖4所示。

圖3 密云水庫(kù)2014—2019年各年面積變化Fig.3 Area Change Each Year of Miyun Reservoir in 2014 to 2019

圖4 密云水庫(kù)2014—2019年總體面積變化Fig.4 Area Change of Miyun Reservoir in 2014 to 2019

1）密云水庫(kù)2014年年內(nèi)面積變化總體呈現(xiàn)逐月減小的趨勢(shì)；2015年年內(nèi)面積變化總體呈現(xiàn)先減小后增大的趨勢(shì)，變化波動(dòng)范圍較?。?016年、2017年和2018年年內(nèi)面積變化總體呈現(xiàn)逐月增加的趨勢(shì)，其中2016年水庫(kù)面積變化波動(dòng)范圍較大；2019年年內(nèi)面積變化呈現(xiàn)波動(dòng)增加的趨勢(shì)，變化波動(dòng)范圍很小。

2）2014—2019年間密云水庫(kù)面積總體呈現(xiàn)先減小后增大的趨勢(shì)，如表1所示，2019年面積較2014年極大增加。在2014—2019年間，2019年6月密云水庫(kù)面積達(dá)到最大值為145.65 km2，2015年7月面積達(dá)到最小值為64.39 km2，2015年平均面積最小，2015年之后面積逐漸增加。2014—2015年水庫(kù)面積呈現(xiàn)負(fù)變化，其余年份間均呈現(xiàn)正變化，2015—2016年面積變化強(qiáng)度最大。

表1 2014—2019年密云水庫(kù)年變化幅度、年變化強(qiáng)度Tab.1 Annual Change Range and Intensity of Miyun Reservoir in 2014 to 2019

3）密云水庫(kù)在夏季出現(xiàn)暴雨、大暴雨時(shí)，水庫(kù)面積會(huì)出現(xiàn)一個(gè)增加的峰值，如2014年7月和2015年8月，其他時(shí)間密云水庫(kù)面積變化無(wú)明顯規(guī)律，受人為因素影響較大。

4）密云水庫(kù)是北京重要的供水水庫(kù)，承擔(dān)北京的生活、灌溉用水，自南水北調(diào)工程通水后，極大緩解了密云水庫(kù)的供水壓力，由于近年來(lái)南水北調(diào)工程輸水、周?chē)畮?kù)向密云水庫(kù)輸水、密云水庫(kù)流域降雨增加、供水減少的“三增一減”，有效解決了密云水庫(kù)干涸問(wèn)題［13］，近年來(lái)密云水庫(kù)面積逐年增加，2019年面積達(dá)到近年來(lái)最大值。

3.2 官?gòu)d水庫(kù)面積變化分析

對(duì)2014—2019年各年年內(nèi)的官?gòu)d水庫(kù)面積變化進(jìn)行統(tǒng)計(jì)得到的結(jié)果如圖5所示，2014—2019年官?gòu)d水庫(kù)面積總體變化趨勢(shì)如圖6所示。

圖5 官?gòu)d水庫(kù)2014—2019年各年面積變化Fig.5 Area change each year of Guanting Reservoir in 2014 to 2019

圖6 官?gòu)d水庫(kù)2014—2019年總體面積變化Fig.6 Area Change of Guanting Resevoir in 2014 to 2019

1）2014—2018年年內(nèi)官?gòu)d水庫(kù)面積變化呈現(xiàn)先增大后減小再增大的趨勢(shì)，2019年年內(nèi)面積變化總體呈現(xiàn)先減小后增大的趨勢(shì)，除2019年外，官?gòu)d水庫(kù)各年年內(nèi)面積變化波動(dòng)范圍均不大，如表2所示。官?gòu)d水庫(kù)面積一般1月至4月開(kāi)始增加，隨著溫度升高，冰雪消融，河流結(jié)冰處融化，水量增加，水庫(kù)得到補(bǔ)充，加上水庫(kù)水量蒸發(fā)小，灌溉用水需求量不大，水庫(kù)面積增加；水庫(kù)面積一般5月至8月逐漸減小，由于溫度繼續(xù)升高，蒸發(fā)量增加，且夏季灌溉用水需求量增加，水庫(kù)面積減小。由于京津冀地區(qū)降雨主要集中在夏季，會(huì)發(fā)生多次中到大雨，甚至暴雨和大暴雨，會(huì)在降雨時(shí)間結(jié)束后出現(xiàn)一個(gè)短暫的面積增加，但由于夏季蒸發(fā)量大于降雨量，降雨結(jié)束后一段時(shí)間水庫(kù)面積會(huì)繼續(xù)縮減；水庫(kù)面積9月至12月逐漸增加，由于氣溫逐漸降低，蒸發(fā)量減小，官?gòu)d水庫(kù)面積有所增加。

表2 2014—2019年官?gòu)d水庫(kù)年變化幅度、年變化強(qiáng)度Tab.2 Annual Change Range and Intensity of Guanting Reservoir in 2014 to 2019

2）在2014—2019年間官?gòu)d水庫(kù)面積總體呈現(xiàn)波動(dòng)增加趨勢(shì)，最大面積出現(xiàn)在2019年1月，為101.26 km2，最小面積出現(xiàn)在2015年9月，面積為55.93 km2，在2014—2019年間，官?gòu)d水庫(kù)面積最大值一般出現(xiàn)在1月至4月或9月至12月，面積最小值一般出現(xiàn)在6月至9月。2014—2018年官?gòu)d水庫(kù)平均面積逐年增加，2019年平均面積較2018年略微減小，2018—2019年官?gòu)d水庫(kù)面積為負(fù)變化，其余均為正變化，2015—2016年年變化強(qiáng)度最大。

3）官?gòu)d水庫(kù)面積變化季節(jié)性規(guī)律較強(qiáng)，受氣溫和降水影響較大，官?gòu)d水庫(kù)自1997年不再向城市供應(yīng)生活用水，主要負(fù)責(zé)灌溉和防洪，因此面積變化受人為因素影響較小。

4）隨著近年來(lái)官?gòu)d水庫(kù)流域降雨量增加和南水北調(diào)工程的供水，緩解了官?gòu)d水庫(kù)對(duì)工農(nóng)業(yè)用水的供水壓力，使得官?gòu)d水庫(kù)面積增加，2019年面積達(dá)到近年來(lái)最大值。

4 結(jié)束語(yǔ)

本文利用NDWI自適應(yīng)閾值法實(shí)現(xiàn)對(duì)水體的有效提取，對(duì)于結(jié)冰水面，采用4波段相加之和的自適應(yīng)閾值法進(jìn)行提取，結(jié)合提取結(jié)果對(duì)水庫(kù)面積進(jìn)行統(tǒng)計(jì)及分析，得到的結(jié)果如下：

2014—2019年間密云水庫(kù)流域冬季氣溫較官?gòu)d水庫(kù)流域氣溫更低，密云水庫(kù)2014年2月、2017年1月和2019年2月均出現(xiàn)了結(jié)冰現(xiàn)象，官?gòu)d水庫(kù)冬季沒(méi)有出現(xiàn)結(jié)冰現(xiàn)象；由于密云水庫(kù)具有供應(yīng)城市生活用水、防洪灌溉等作用，面積變化受人為因素影響較大，因此密云水庫(kù)各年內(nèi)面積變化季節(jié)性規(guī)律不強(qiáng)；由于周?chē)畮?kù)輸水、南水北調(diào)輸水，密云水庫(kù)流域降雨量增加和密云水庫(kù)本身供水減少，密云水庫(kù)近年來(lái)面積逐漸增加。由于1997年之后官?gòu)d水庫(kù)不再供應(yīng)生活用水，官?gòu)d水庫(kù)面積變化受氣候和降水影響較大，受人為因素影響較小，因此官?gòu)d水庫(kù)各年內(nèi)面積變化季節(jié)性規(guī)律較強(qiáng)，近年來(lái)，官?gòu)d流域降雨量增加，南水北調(diào)工程供水，緩解了官?gòu)d水庫(kù)的工農(nóng)業(yè)用水壓力，官?gòu)d水庫(kù)面積增加。

本文所選取的遙感影像為月度數(shù)據(jù)，為了實(shí)現(xiàn)精確的水庫(kù)變化監(jiān)測(cè)，后續(xù)可以選取時(shí)間間隔更小的遙感影像進(jìn)行水庫(kù)變化分析。