張 珂，吳 南，徐國鑫，劉林鑫1，，范亞洲1，，張企諾，周佳奇1，，劉 挺

(1.河海大學水文水資源與水利工程科學國家重點實驗室，江蘇 南京 210098；2.河海大學水文水資源學院，江蘇 南京 210098； 3.長江保護與綠色發(fā)展研究院，江蘇 南京 210098； 4.中國氣象局-河海大學水文氣象研究聯(lián)合實驗室，江蘇 南京 210098； 5.陜西省引漢濟渭工程建設(shè)有限公司，陜西 西安 710010)

隨著空間技術(shù)的發(fā)展，通過遙感影像提取水體在洪澇災(zāi)害、水資源監(jiān)測以及水環(huán)境保護等領(lǐng)域得到廣泛應(yīng)用[1-3]，水體提取已經(jīng)成為遙感技術(shù)研究領(lǐng)域的重要分支之一[4]。目前常見的水體提取方法有單波段閾值法、譜間關(guān)系法、水體指數(shù)法、決策樹法、密度分割法以及面向?qū)ο蠓ǖ萚5-9]，其中水體指數(shù)法是一種方便高效的水體識別法，應(yīng)用也比較成熟[10]。由于衛(wèi)星的分辨率不同，水體提取結(jié)果也有所差異，鄧勁松等[11]基于SPOT高分辨率影像提出一種波段間的提取方法，有效地將居民地和水體分開；周成虎等[12]針對AVHRR影像提出了自動提取識別水體的描述模型，并在太湖、淮河等多個地區(qū)得到應(yīng)用；王培培[13]利用ETM遙感影像數(shù)據(jù)，通過光譜特征分析和歸一化植被指數(shù)識別水體，實現(xiàn)了水體的自動提取。以上研究所使用的數(shù)據(jù)源均來自國外遙感衛(wèi)星，利用我國高分衛(wèi)星進行水體提取的研究和算法還不多見。

高分一號衛(wèi)星(GF-1)于2013年4月26日發(fā)射成功，是我國第一顆考核壽命大于5 a的衛(wèi)星，衛(wèi)星除裝載有兩臺2 m分辨率全色和8 m分辨率多光譜相機外，還搭載了4臺16 m分辨率多光譜相機。6臺相機的配合使得我國首次實現(xiàn)了單顆衛(wèi)星高分辨率與大幅寬技術(shù)的有效結(jié)合[14]，大大提高了我國高分辨率影像數(shù)據(jù)的自給率[15]，將為我國的科技發(fā)展帶來不可估量的影響[16]。

受天氣影響，遙感影像在地表水體提取過程中易受到云層的影響。遙感影像上受云層影響的區(qū)域覆蓋了大面積的地表信息，使得圖像信息不完善，給影像數(shù)據(jù)的處理和應(yīng)用造成了極大的障礙[17-18]。如何消除遙感影像上云層的影響，至今仍是研究的熱門話題。目前已有不少國內(nèi)外相關(guān)專家和學者對此展開研究，研究結(jié)果表明：同態(tài)濾波法適用于大面積有薄云存在的遙感影像[19-21]，對厚云區(qū)域基本不產(chǎn)生效果；多源數(shù)據(jù)融合法使用不同傳感器的遙感影像，將無云部分替代其他影像的云干擾部分[22]；時間序列傅立葉分析法簡化了壓縮時間序列遙感數(shù)據(jù)，對遙感影像中的云進行檢測、去除和替代處理[23]。除此之外，將遙感影像與地形信息有效結(jié)合的研究也有很多[24-25]，王鑫蕊等[26]利用湖泊邊界位于同一高度的特點，提出了一種利用湖泊的已知部分邊界結(jié)合DEM數(shù)據(jù)，還原完整邊界的自動方法，有效克服了云層的干擾。

屯溪流域為中小河流重點研發(fā)項目的典型示范流域，東方紅水庫是屯溪流域最大的一座中型水庫，其庫容、面積變化影響著該流域的應(yīng)急防控技術(shù)。由于遙感數(shù)據(jù)周期間隔較長，無法實時獲取水庫的面積信息，且研究區(qū)內(nèi)影像常年或多或少地受到云層的影響，降低了遙感數(shù)據(jù)的實用性。為了利用我國現(xiàn)有的光學遙感影像提高遙感數(shù)據(jù)的質(zhì)量，本文以東方紅水庫為研究對象，提出一種GF-1遙感影像結(jié)合DEM消除云層干擾的連續(xù)水體重建法，對于連續(xù)水體局部受云層影響的影像，結(jié)合GF-1遙感影像和高精度等高線數(shù)據(jù)，通過地理配準還原遙感影像水體面積，從而實現(xiàn)云層干擾下水體的自動提取。

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

1.1 研究區(qū)概況

東方紅水庫位于新安江水系二級支流東亭河(虞山溪)上(圖1)，是以防洪、灌溉為主，結(jié)合發(fā)電、旅游供水等綜合利用的中型水利樞紐工程。水庫集水區(qū)域處在117°E～119°E、29°N～31°N之間，集水面積為60.0 km2，主河道長11.8 km。東方紅水庫總庫容0.244 5億m3，設(shè)計灌溉面積0.3萬hm2，電站總裝機容量400 kW，多年平均發(fā)電量131萬kW·h。工程規(guī)模為安徽省重點中型水庫，工程等別為Ⅲ等，攔河壩、溢洪道、引水涵洞，以及由導流隧洞改建的放水底孔等均為3級建筑物。

圖1 研究區(qū)位置和地形Fig.1 Location and topography of the study area

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

遙感數(shù)據(jù)來自中國資源衛(wèi)星應(yīng)用中心(http：//www.cresda.com/CN/index.shtml)的高分一號寬幅傳感器(GF1_WFV)數(shù)據(jù)，影像空間分辨率為16 m。通過實地查勘的方式，收集了東方紅水庫高精度等高線數(shù)據(jù)及2015年水庫詳細運行資料。2015年為豐水年，東方紅水庫集水區(qū)域平均降水量為2 411.6 mm。通過2020年2月陸面過程資料分布式數(shù)據(jù)存檔中心LP DAAC(Land Processes Distributed Active Archive Center)公布的1 rad/s分辨率NASADEM數(shù)據(jù)集產(chǎn)品(https：//lpdaac.usgs.gov/news/release-nasadem-data-products/)提取研究流域包含的水體區(qū)域。為體現(xiàn)水庫在不同時期的蓄水面積信息，選取2015年非雨期(4、10月)、雨期(5、8、9月)共7景衛(wèi)星影像。

2 研究方法

2.1 數(shù)據(jù)處理

2.1.1 遙感影像數(shù)據(jù)處理

選擇中國資源衛(wèi)星應(yīng)用中心公布的1A級影像數(shù)據(jù)產(chǎn)品。為降低遙感影像在成像過程中大氣傳輸、傳感器本身、地球曲率等多種因素的影響，需要對原始遙感影像進行預處理。結(jié)合IDL(interactive data language)與ENVI(the environment for visualizing images)對數(shù)據(jù)進行批量預處理，預處理的方式主要包括：輻射定標、大氣校正[27]、正射校正等。預處理結(jié)束后，將研究區(qū)域?qū)?yīng)的矢量文件與遙感影像配準，完成研究區(qū)域的裁剪。2015年7景水體局部受云層影響的遙感影像數(shù)據(jù)見表1。

表1 高分一號衛(wèi)星遙感影像數(shù)據(jù)集

2.1.2 高精度等高線數(shù)據(jù)處理

利用線性均分方式將東方紅水庫對應(yīng)的高精度等高線數(shù)據(jù)進行插值，并選取多幅高保證率歷史遙感影像，通過最大類間方差膨脹迭代法提取水體，從而得到水庫不同水位值對應(yīng)的水體邊緣輪廓線，以此為參考，利用GIS軟件對插值得到的高精度等高線數(shù)據(jù)進行人工解譯，精準修正等高線數(shù)據(jù)，并確保高程數(shù)據(jù)完全覆蓋水體范圍。

2.2 水體提取與云層消除

2.2.1 最大類間方差膨脹迭代法提取水體

歸一化差分水體指數(shù)(normalized difference water index，NDWI)法[9]利用遙感影像的綠波與近紅外波兩個特定波段進行歸一化差值處理，從而突顯影像中的水體。將經(jīng)過預處理的遙感影像利用NDWI法初步提取水體。

最大類間方差法由日本學者Ostu等[28]提出，適用于雙峰情況的閾值分割，可以實現(xiàn)自動求取閾值的目的。該方法根據(jù)圖像的灰度值，將圖像分為背景和目標兩部分。分割閾值從最小灰度值到最大灰度值遍歷每個值，不同的閾值分割使得背景和目標之間的方差不同，類間方差越大說明被閾值分割的兩部分差別越大。若出現(xiàn)目標與背景的錯分、漏分等現(xiàn)象，均會導致兩部分差別變小。因此，利用類間方差最大時對應(yīng)的閾值進行分割，被錯分的概率應(yīng)該是最小的。

利用NDWI法提取未受云層影響的局部水體，將初始閾值T0設(shè)為-0.02，若歸一化差分水體指數(shù)INDWIT0，得到的結(jié)果為初步提取的水體。以此基于形態(tài)學膨脹算法建立水體的緩沖區(qū)。當膨脹算法中的結(jié)構(gòu)元素大小發(fā)生變化時，可得到不同范圍的緩沖區(qū)，結(jié)構(gòu)元素過大、過小都會造成水體的誤提。根據(jù)歷史經(jīng)驗，該研究區(qū)選取大小為15×15的結(jié)構(gòu)元素膨脹提取效果最佳，故以該結(jié)構(gòu)元素大小對初步提取的部分水體進行膨脹，建立緩沖區(qū)，對緩沖區(qū)內(nèi)包含的區(qū)域應(yīng)用最大類間方差法確定水體分割閾值，從而實現(xiàn)水體的提取。以前一次提取的水體區(qū)域作為后一次建立緩沖區(qū)的初始范圍，重復利用最大類間方差法進行多次迭代，自適應(yīng)地確定水體與周邊地物的分割閾值，由此得到每次迭代確定的水體面積，當相鄰兩次水體面積差值絕對值滿足要求精度時，迭代結(jié)束。

|Sk-Sk-1|<0.01 (k=1，2,3,…)

(1)

S0=Sk

(2)

式中：Sk——第k次迭代得到的水體面積；Sk-1——第k-1次迭代得到的水體面積；S0——最大類間方差膨脹迭代法提取到的最優(yōu)水體面積。

2.2.2 基于高精度等高線數(shù)據(jù)的云層影響消除法

通過最大類間方差膨脹迭代得到東方紅水庫的局部水體范圍，與處理后的高精度等高線數(shù)據(jù)進行地理坐標配準，可得到每個水體柵格對應(yīng)的水位高程值。某一時刻連續(xù)水體水位應(yīng)為唯一值，而由于云層的干擾以及誤差影響導致遙感影像提取的邊界水位并不唯一，故對云層干擾下由遙感數(shù)據(jù)得到的水體不同水位進行數(shù)學統(tǒng)計計算。

定義受云層影響遙感影像得到的局部水體為W類水體，局部水體邊界網(wǎng)格以內(nèi)的內(nèi)部水體網(wǎng)格為W1，水庫原始外輪廓邊緣線對應(yīng)的網(wǎng)格為W2，受云層影響而產(chǎn)生的邊界線對應(yīng)的水體網(wǎng)格為W3，其他非水體網(wǎng)格為W0。

若W2類水體網(wǎng)格對應(yīng)的不同水位值共n個，對此進行均值計算：

(3)

SM⊕SV=SV

(4)

式中：SM——受云層影響遙感影像得到的局部水體面積；SV——還原水體水位值所包圍的區(qū)域面積。兩類水體經(jīng)地理配準后可實現(xiàn)云層干擾下的水體還原提取。

3 結(jié)果與分析

在IDL + ENVI編程環(huán)境中實現(xiàn)遙感影像數(shù)據(jù)的預處理。將經(jīng)過預處理的圖像在MATLAB中進行最優(yōu)閾值的確定，再通過GIS軟件確定最優(yōu)閾值對應(yīng)的水體范圍，與水體范圍對應(yīng)的高精度等高線數(shù)據(jù)進行地理配準，從而實現(xiàn)云層干擾下水體的還原處理。選取GF-1_2015-09-08景遙感影像(云量為11%)進行有云水體提取，各步驟處理結(jié)果如圖2所示。

圖2 水體GF-1遙感影像云消除過程Fig.2 Cloud removal process in water body extraction for the GF-1 satellite imageries

圖2(a)為東方紅水庫集水區(qū)域2015-09-08GF-1_2015-09-08景遙感影像經(jīng)過輻射定標、大氣校正、正射校正、配準及研究區(qū)裁剪等預處理后的圖像。由圖2(a)可以看出，紅色區(qū)域為水庫水體區(qū)域，部分水體受到云層影響。圖2(b)為云層干擾下通過最大類間方差膨脹迭代提取的東方紅水庫水體面積圖，由于受云層影響較為嚴重，水庫外輪廓線提取并不完整，且水庫內(nèi)部也缺失部分水體。圖2(c)為東方紅水庫云干擾水體提取圖與其高精度等高線數(shù)據(jù)的地理配準圖，圖中水體與等高線數(shù)據(jù)產(chǎn)生輕微偏移，水位值并不唯一。圖2(d)為通過地理配準還原云層干擾下的東方紅水庫水體面積圖，較為清晰地還原了水庫的基本輪廓。

根據(jù)東方紅水庫289～297 m水位區(qū)間對應(yīng)的水位-面積曲線(圖3)以及觀測資料記錄的每日水庫實測水位值(表2)，可以推導出對應(yīng)時間水庫的水體面積，以此作為SV的參考值。

圖3 東方紅水庫水位-面積關(guān)系曲線Fig.3 Water level-area relationship curve of Dongfanghong Reservoir

從表2可以看出，東方紅水庫GF-1_2015-09-08景遙感影像實測水位值h為291.35 m，根據(jù)水位-面積關(guān)系曲線得到的參考面積S為1.995 6 km2。然而由于云層的干擾，根據(jù)最大類間方差膨脹迭代得到的水體面積為1.060 7 km2，面積差值為0.934 9 km2，差距較大。按照本文提出的通過地理配準還原云層干擾水庫水體面積，可得到水庫水位還原值為291. 95 m，還原后庫區(qū)水體面積為2.042 km2，面積差值為0.046 4 km2，可以看出，該方法簡單有效地還原了水庫水體面積。

表2 消除云層干擾后的庫區(qū)水體遙感提取面積

另選取東方紅水庫2015年6景影像并使用上述方法進行云消除處理。圖4為最大類間方差膨脹迭代提取的水體范圍與其等高線數(shù)據(jù)的地理配準圖，藍色區(qū)域為S0。從圖4可以看出云層干擾對東方紅水庫水體提取范圍影響較大，嚴重時水庫外輪廓基本無法確定。除此之外還可以看出，庫區(qū)水體外輪廓線對應(yīng)的水位值并不唯一，因此需要通過數(shù)學統(tǒng)計計算方式確定需要還原的水位值。

圖4 地理配準后的衛(wèi)星影像Fig.4 Satellite imageries after geographic referencing

圖5為經(jīng)地理配準后還原的水體，與等高線數(shù)據(jù)相結(jié)合經(jīng)均值計算后，水庫對應(yīng)的水位值可唯一確定，最終結(jié)果表明該方法簡單有效地消除了云層的干擾，較準確地還原了東方紅水庫在非雨期(4、10月)、雨期(5、8、9月)云層影響下的水體。從圖5(a)(b)(f)可以看出，由于水庫在4、10月來水量較少，隨著下游水量的補給，庫區(qū)水體范圍逐漸減小，提取的水體范圍較??；圖5(c)(d)(e)表明，5、8、9月降水量較多，水庫來水量大于出水量，蓄水量逐漸增多，故水體面積增大，水位升高，庫區(qū)之間連接的細小河道蓄水位也逐漸升高。

圖5 水體還原結(jié)果Fig.5 Results of extracted water bodies

從表2可以看出，S0與S相比差距較大，其中GF-1_20150426景遙感影像水體面積相差最大，達1.368 4 km2，難以反映真實水體。選取的7景影像經(jīng)地理配準后得到的水體還原值與實測水位參考值接近，還原后的水體水位平均差值為0.438 6 m，有效水體面積平均差值為0.080 9 km2，平均相對誤差在5%以內(nèi)，較為精確地還原了庫區(qū)的水體。其中GF-1_20150512景遙感影像水體庫區(qū)面積差值僅為0.008 3 km2，較好地消除了云層干擾，還原了水體水位值及對應(yīng)的面積。

從表2還可以看出，GF-1_20150401、GF-1_20150426和GF-1_20151003景遙感影像還原的水位值與實測水位值差距略微偏高，這可能是由于所用影像分辨率較低，為16 m，在預處理過程中水體范圍整體產(chǎn)生輕微偏移，使得提取的水體存在部分誤差，且該方法對地形數(shù)據(jù)精度要求較高，本文將1 m等高線經(jīng)線性插值、人工修訂后得到的等高線與提取的遙感信息進行地理配準，不同數(shù)據(jù)源之間易出現(xiàn)誤差疊加現(xiàn)象。

4 結(jié) 語

高分一號衛(wèi)星遙感影像受云層影響限制了數(shù)據(jù)的實用性，且在云層干擾下通過最大類間方差膨脹迭代法得到的水體并不完整，面積平均差值為1.009 4 km2，無法準確反應(yīng)水體的真實信息，難以直接用于陸地表面水體的提取。通過本文提出的地理配準還原云層干擾下水體的方法，還原后的水體水位平均差值為0.438 6m，有效水體面積平均差值為0.080 9 km2，最小僅為0.008 3 km2，平均相對誤差在5%以內(nèi)，有效地克服了云層干擾，較為精確地還原了庫區(qū)水體。該方法簡單有效，與高精度等高線數(shù)據(jù)相結(jié)合，通過簡單的圖像計算便可克服云層的干擾，提高了衛(wèi)星遙感影像的實用性，可廣泛用于地表水資源的動態(tài)提取，為水資源調(diào)查及地表水體動態(tài)監(jiān)測等領(lǐng)域提供高質(zhì)量的遙感影像數(shù)據(jù)。

通過本文方法得到的水體面積雖然有效克服了云層的干擾，但該方法對地形數(shù)據(jù)精度要求較高，對于缺少水下高程數(shù)據(jù)的水體來說，水位動態(tài)變化值常在某范圍內(nèi)上下變動，故可根據(jù)多景無云歷史遙感影像(包含豐水期、枯水期較多高程值變化)得到的水位值對水下空缺的高程值進行填補，再進一步進行人工修正，但精度方面問題有待提高。除此之外，對于一些復雜的水體情況，如湖心島等，該方法需要結(jié)合人工解譯才能進行水體提取，對此仍需要進一步的深入研究。