李 闖,趙建輝,2,3,戚雯雯,常繼科
1(河南大學 計算機與信息工程學院,開封 475004)
2(河南大學 河南省大數(shù)據(jù)分析與處理重點實驗室,開封 475004)
3(河南大學 河南省智能技術與應用工程技術研究中心,開封 475004)
近年來,由于氣候變化,極端天氣增多,洪災仍是給人類社會帶來損失最大的自然災害[1].水庫作為一種重要的水利樞紐,對攔截洪水、保障供水、調配水資源等工作具有積極作用,也為社會經濟發(fā)展和生態(tài)文明建設提供著不可忽略的助力.水庫蓄水量是水庫健康運行的重要指標之一,與區(qū)域的氣候、濕度、降水等密切相關.因此,研究水庫蓄水量對合理保護和利用水資源以及制定相關政策十分必要[2].
傳統(tǒng)的水庫蓄水量測量方法主要依靠人力到現(xiàn)場布設各種監(jiān)測儀器來獲取數(shù)據(jù)[3],通過繪制庫區(qū)地圖并將其分割成等高線、三角網等,進行水域面積的提取,然后結合等高距,計算水庫蓄水量.傳統(tǒng)方法不僅耗費大量的人力、物力和財力,而且水域面積的測量精度受制圖精度的影響較大,所利用的儀器也常常會因為天氣惡劣和日常消耗而損壞,從而導致數(shù)據(jù)的缺失.相對于傳統(tǒng)方法,遙感方法具有監(jiān)測范圍廣、收集信息快和工作平臺多樣化等優(yōu)勢.單一遙感方法主要利用遙感技術對水庫的一種參數(shù)進行提取[4],其他參數(shù)則采用現(xiàn)場調查的方式獲取,這種方法研究手段單一,在同一時期獲取的可用遙感數(shù)據(jù)少,精度相對較低.而對于我國國土面積遼闊、水庫分布零散的現(xiàn)狀而言,并不是所有的水庫都適合到現(xiàn)場開展調查.例如青藏高原,由于海拔高、氣候惡劣、地形地貌復雜等原因,會在一定時期造成人類難以到達的無人區(qū),使之不宜開展現(xiàn)場調查.采用多源遙感能夠豐富數(shù)據(jù)源,使不同遙感技術優(yōu)勢互補,獲取比單一遙感更好的監(jiān)測效果.因此,基于多源遙感信息進行水庫蓄水量估算研究具有十分重要的現(xiàn)實意義.
目前已有許多基于多源遙感信息的水資源監(jiān)測研究.2016年,文獻[5]基于Envisat、Jason-2和Terra/Aqua MODIS 等多源衛(wèi)星遙感數(shù)據(jù),獲取了河流水位和河寬數(shù)據(jù),利用改進的曼寧公式估算了亞馬遜和尼羅河的徑流量.2017年,文獻[6]利用合成孔徑雷達(Synthetic Aperture Radar,SAR)影像、光學遙感影像和機載激光雷達(Light Detection And Ranging,LiDAR)數(shù)據(jù),建立了3 個數(shù)據(jù)集的融合分類模型,采用多級決策樹,通過最小化單個數(shù)據(jù)集的模型之間的差異實現(xiàn)水體的高精度提取.2019年,文獻[7]綜合運用衛(wèi)星測高技術和光學遙感岸線監(jiān)測技術,研制了高時間分辨率的2000年-2017年青藏高原湖泊水位、水量變化數(shù)據(jù)集.多源遙感衛(wèi)星技術已在地表水資源監(jiān)測研究中得到廣泛應用,為地表水資源監(jiān)測提供了可靠的理論和技術支撐,但目前基于多源遙感信息的水庫蓄水量估算研究非常少見,基于多源遙感信息的水庫蓄水量估算系統(tǒng)更是鮮有報道.
近年來,已有研究人員開發(fā)了水庫防洪調度[8,9]、庫容計算[10]、水庫地質災害視頻監(jiān)測預警[11,12]等多種水庫信息管理系統(tǒng),但是這些系統(tǒng)所利用的數(shù)據(jù)多是地面實測數(shù)據(jù)或是第三方提供的數(shù)據(jù),很少利用多源遙感技術直接在軟件系統(tǒng)中提取水庫參數(shù).2018年,文獻[13]設計和實現(xiàn)了一款基于SAR 遙感監(jiān)測的水庫參數(shù)計算系統(tǒng),利用從SAR 影像提取的水域面積、從國家水利部“全國水雨情信息網”獲取的實測水位和標準蓄水量等數(shù)據(jù),通過數(shù)據(jù)同化,對丹江口水庫的面積和蓄水量等參數(shù)進行了估算研究,但該文的蓄水量估算模型僅利用SAR 遙感影像自動提取水域面積參數(shù),模型所需的實測水位等關鍵參數(shù)仍需從第三方網站獲取,對實測數(shù)據(jù)的依賴較大,而且受限于網絡傳輸速率等方面的限制,一旦缺乏數(shù)據(jù)支持(如重大災害導致數(shù)據(jù)網絡中斷等),該模型和系統(tǒng)將無法獨立工作.因此,開發(fā)出一款具有獨立工作能力的多源遙感水庫蓄水量估算系統(tǒng)具有十分重要的研究意義和應用價值.
針對上述問題,本文利用星載雷達高度計(Rader Altimetry,RA)[14]提取水庫水位,利用優(yōu)勢互補的星載多光譜成像儀(Multi-Spectral Instrument,MSI)影像和星載SAR 影像提取水庫的水域面積,采用數(shù)據(jù)同化法建立了多源遙感水域面積差異的定量函數(shù)關系,構建了多源遙感水庫蓄水量估算模型,并基于該模型,設計和實現(xiàn)了一款多源遙感水庫蓄水量估算系統(tǒng).
本系統(tǒng)在水庫參數(shù)提取中,基于多源遙感技術獲取水位和水域面積等參數(shù),利用數(shù)據(jù)同化法建立了多源遙感水庫蓄水量估算模型,其技術路線如圖1所示.在系統(tǒng)實現(xiàn)上,采用Matlab和Java 混合編程技術來完成數(shù)據(jù)處理功能;利用虛擬水文站點技術使水位以點標形式呈現(xiàn)在百度地圖API 上;使用ECharts 技術將系統(tǒng)中數(shù)字化信息轉化成可視化圖表信息.
1.1.1 提取遙測水位
系統(tǒng)利用RA 數(shù)據(jù)提取水庫的遙測水位.提取過程中,要進行一系列的誤差校正,包括干對流層Rangedt、濕對流層Rangewt、電離層Rangei、固體潮Rangeset和極潮Rangept等測距誤差校正,之后需要進行水面橢球高GeoEGM96轉化,使遙測水位與實測水位一樣,都以中國黃海水準面為基準,如式(1)、式(2)所示,其中H為遙測水位,Altitude為衛(wèi)星到地球橢球面距離,Range為衛(wèi)星到水面距離,RangeCorrected為誤差校正距離.
圖1 多源遙感水庫蓄水量估算模型的技術路線
1.1.2 提取水域面積
系統(tǒng)使用MSI 影像和SAR 影像作為獲取水域面積的數(shù)據(jù)源,對多源遙感數(shù)據(jù)進行預處理,利用水陸分割算法提取水庫的水體輪廓,根據(jù)水體輪廓計算面積[15].水域面積提取原理為:先統(tǒng)計整幅影像中總像元個數(shù)M,再統(tǒng)計水庫水域面積的像元個數(shù)N,設整幅影像對應到現(xiàn)實地面中的實際大小為S,則實際水域面積P的計算公式如式(3).
1.1.3 估算蓄水量
多源遙感水域面積提取后,對其精度進行定量分析,基于數(shù)據(jù)同化技術,構建兩者水域面積差異的定量函數(shù)關系,以SAR 影像為基準,修正水域面積差異,構建“遙測水位-多源遙感水域面積”模型和“多源遙感水域面積-標準蓄水量”模型,然后將上述兩種模型進行數(shù)學轉化,建立本系統(tǒng)所使用的多源遙感水庫蓄水量估算模型,如式(4)所示,其中Impoundage為估算蓄水量,H為遙測水位.
在上述過程中,標準蓄水量用于輔助建模,模型建成后,在模型使用過程中無需再使用標準蓄水量,克服了當前蓄水量估算模型建立之后仍需輸入實測數(shù)據(jù)的不足,在一定程度上完善了實測數(shù)據(jù)缺失情況下的研究工作,具有一定現(xiàn)實意義.
系統(tǒng)采用C/S (Client/Server,客戶端-服務器)架構,使用Matlab和Java 混合編程技術來完成數(shù)據(jù)處理功能.利用Matlab 強大的影像處理功能,進行遙感影像水陸分割算法以及RA 提取遙測水位算法的編程實現(xiàn).利用Java 統(tǒng)計水陸分割之后水庫水體輪廓的像素信息,進而獲得具體水域面積的值.在Java 與Matlab 混合編程過程中,先將各種處理算法的Matlab 程序打包,通過設置調用方法的Class 類名導入Java 工程中,然后配置MCR (Matlab Compiler Runtime,Matlab 編譯運行)環(huán)境,并安裝MCR 插件,即可完成相關的調用.系統(tǒng)調用框架圖如圖2所示.
圖2 系統(tǒng)調用框架圖
處理RA 原始數(shù)據(jù)時,僅根據(jù)處理之后遙測數(shù)據(jù)的經緯度數(shù)值,難以判斷所提取位置是否在水庫水體上.系統(tǒng)利用百度地圖API 抽取RA 數(shù)據(jù)中每一個數(shù)值數(shù)據(jù)的經緯度,以點標形式將其顯示在百度地圖上,從而達到虛擬水文站的效果.百度地圖API 使用獨立的坐標體系,因此在地圖上使用經緯度坐標前還需要調用百度地圖API 接口進行坐標轉換,之后就可以以點標形式顯示在地圖上,如圖3所示,用戶能夠直觀看到RA 數(shù)值數(shù)據(jù)的位置信息,并可以通過點擊水庫水體上的各個點標,來查看相應的水位值.
圖3 在百度地圖中顯示的虛擬水文站點
ECharts是一個使用JavaScript 實現(xiàn)的開源可視化庫,不僅提供了折線圖、雷達圖、和弦圖等各種常用的圖表類型來滿足不同用戶的需求,還提供了一些相對專業(yè)的圖表來滿足專業(yè)性較高的要求.本系統(tǒng)中ECharts 可視化處理流程如圖4所示.
圖4 ECharts 可視化處理流程
本文從系統(tǒng)框架、流程以及數(shù)據(jù)庫等3 個方面對系統(tǒng)的總體設計進行分析與介紹.
系統(tǒng)分為數(shù)據(jù)管理、數(shù)據(jù)處理和結果展示等3 個功能模塊,模塊結構如圖5所示.
圖5 系統(tǒng)功能模塊結構圖
(1)數(shù)據(jù)管理模塊
本模塊主要包括水文數(shù)據(jù)管理、遙感數(shù)據(jù)管理和模型管理等功能.水文數(shù)據(jù)管理主要是記錄水文站的實測數(shù)據(jù);遙感數(shù)據(jù)管理主要進行RA 數(shù)據(jù)、MSI 影像數(shù)據(jù)和SAR 影像數(shù)據(jù)的上傳、檢索以及處理后的下載與刪除等操作;模型管理主要實現(xiàn)系統(tǒng)內各算法和模型的增加、刪除和修改等功能.
(2)數(shù)據(jù)處理模塊
本模塊主要實現(xiàn)各種數(shù)據(jù)處理功能,包括利用RA 提取算法提取遙測水位;利用水陸分割算法對預處理后的SAR 影像和MSI 影像進行水陸分割,從而得到水體輪廓;調用面積統(tǒng)計算法計算具體的水域面積值;結合遙測水位、水域面積、標準蓄水量,建立蓄水量估算模型,并估算水庫蓄水量.
(3)結果展示模塊
本模塊主要實現(xiàn)水庫水位、水域面積和蓄水量等結果的展示功能,利用建好的相關數(shù)據(jù)庫,實現(xiàn)水位、水域面積和蓄水量之間的相互匹配,并以曲線圖和表格等方式予以展示.
系統(tǒng)業(yè)務流程如圖6所示,通過數(shù)據(jù)輸入、數(shù)據(jù)處理和結果展示等3 大步驟,實現(xiàn)基于多源遙感數(shù)據(jù)的水庫蓄水量估算功能.
圖6 系統(tǒng)業(yè)務流程圖
(1)數(shù)據(jù)輸入
該步驟用于輸入模型處理所需的各種數(shù)據(jù),主要包括實測數(shù)據(jù)、RA 數(shù)據(jù)和遙感影像數(shù)據(jù)等3 種類型.
實測數(shù)據(jù):用戶上傳每日記錄的Excel 文件,或者直接在系統(tǒng)里利用爬蟲技術爬取互聯(lián)網上的相關數(shù)據(jù),包括實測水位和蓄水量.
RA 數(shù)據(jù):從網上下載原始的RA 數(shù)據(jù),經過解壓之后存儲到數(shù)據(jù)庫里.
遙感影像數(shù)據(jù):包括預處理后的MSI 影像和SAR影像,上傳時需要填寫時間、衛(wèi)星名稱、所屬水庫等相關信息.
(2)數(shù)據(jù)處理
該步驟通過調用各種算法來計算水位、水域面積、蓄水量等水庫參數(shù),是水庫蓄水量估算模型的主要實現(xiàn)步驟.
水位計算:遙測水位是調用RA 提取算法處理RA原始數(shù)據(jù)而得到的水位,在百度地圖中顯示了遙測水位對應在RA 衛(wèi)星星下點的位置,方便用戶選取合適的虛擬水文站點.
水域面積計算:輸入系統(tǒng)的SAR 影像已進行斑點濾波和輻射定標等預處理操作,在系統(tǒng)中只需調用水陸分割算法區(qū)分水體和非水體就能完成水庫輪廓提取;輸入系統(tǒng)的MSI 影像已進行水體信息增強、水體與非水體區(qū)分等預處理步驟,在系統(tǒng)中僅需進行連通域識別,以去除非水體小連通域.通過以上步驟獲得MSI 影像和SAR 影像中的水體輪廓后,再調用面積統(tǒng)計算法,計算水庫水域面積值.
蓄水量估算:基于數(shù)據(jù)同化技術,將遙測水位、水域面積和標準蓄水量輸入蓄水量估算模型,估算水庫蓄水量,并繪制多種過程曲線圖和表格來展示結果.
(3)結果展示
該步驟主要以曲線圖和表格等方式展示各種處理結果.系統(tǒng)里的過程曲線主要有“水位-日期”、“水位-水域面積”和“水位-蓄水量”等過程曲線,其中,當用戶點擊“水位-日期”過程曲線圖時,在“水位-水域面積”和“水位-蓄水量”等過程曲線圖里相應的標記會隨之改變,在結果表格中會顯示具體的數(shù)值,包括日期、實測水位、遙測水位、水域面積和蓄水量等.
系統(tǒng)數(shù)據(jù)庫用于存放系統(tǒng)處理過程中所涉及到的各種數(shù)據(jù).該系統(tǒng)的實體-聯(lián)系圖(Entity Relationship diagram,E-R 圖)如圖7所示.根據(jù)系統(tǒng)E-R 圖,共設計了水庫信息表、實測水位結果表、RA 原始數(shù)據(jù)表、遙測水位結果表、MSI 影像信息表、SAR 影像信息表、影像處理算法表、水域面積表、蓄水量模型表、水庫蓄水量表等共計10 個數(shù)據(jù)表.
圖7 系統(tǒng)E-R 圖
系統(tǒng)采用C/S 架構,客戶端提供數(shù)據(jù)管理、功能選擇、結果顯示等用戶交互界面,服務器端提供數(shù)據(jù)處理、參數(shù)計算等后臺響應功能.
數(shù)據(jù)管理模塊主要對系統(tǒng)中各種數(shù)據(jù)進行管理,實現(xiàn)實測數(shù)據(jù)、遙測水位、水域面積、影像管理、模型管理、水庫管理等數(shù)據(jù)的上傳、查詢、修改、下載和刪除等功能.在圖8所示的系統(tǒng)首頁中,點擊左側的“數(shù)據(jù)管理”按鈕,即可進入圖9所示的數(shù)據(jù)管理頁面,對系統(tǒng)數(shù)據(jù)進行各種管理操作.其中,SAR 影像的上傳功能如圖10所示,用戶填好相關信息之后,選擇SAR影像路徑或直接把SAR 影像文件拖拽到上傳控件,點擊“下一步”按鈕,客戶端即可調用文件傳輸接口,以流式文件方式向服務器端上傳SAR 影像文件.影像下載功能實現(xiàn)方式類似,在圖9數(shù)據(jù)管理頁面中點擊右側“下載”按鈕后,客戶端會向服務器端發(fā)送SAR 影像記錄ID,服務器端收到請求后會將SAR 影像流式文件傳輸至客戶端,并彈出如圖11所示的下載任務確認窗口,用戶可以選擇存放SAR 影像的路徑,點擊“下載”按鈕后,完成SAR 影像的下載.
圖8 多源遙感水庫蓄水量估算系統(tǒng)首頁
圖9 數(shù)據(jù)管理頁面
數(shù)據(jù)處理模塊主要包含水位提取、水域面積提取、蓄水量估算等3 個功能.系統(tǒng)通過調用RA 提取算法的Matlab 包來處理RA的原始數(shù)據(jù)進而提取水庫的遙測水位;調用水陸分割算法的Matlab 包實現(xiàn)對遙感影像中水體輪廓的提取,然后調用面積統(tǒng)計算法的Java 程序對水體輪廓所代表的水域面積進行計算;最后結合數(shù)據(jù)同化算法估算出水庫蓄水量.
圖10 上傳SAR 影像
圖11 下載SAR 影像
(1)水位提取
系統(tǒng)中使用的水位數(shù)據(jù)主要有兩種,一種是水文站的實測水位,另一種是從RA 數(shù)據(jù)中提取的遙測水位,其中實測水位的作用是對遙測水位進行驗證.對于水文站實測水位,系統(tǒng)提供了兩種獲取方式,一種是人工導入方式,由用戶每天登錄國家水利部“全國水雨情信息網”進行記錄,并以Excel 文件形式人工導入系統(tǒng);另一種是自動爬取方式,系統(tǒng)利用爬蟲技術自動從“全國水雨情信息網”爬取水文站實測水位及官方公布的標準蓄水量,并自動存入系統(tǒng)數(shù)據(jù)庫.圖12是利用爬蟲技術爬取實測水位的開始頁面,用戶點擊圖中右上角“獲取今日實測數(shù)據(jù)”按鈕后,系統(tǒng)將自動爬取相關數(shù)據(jù)并直接存入系統(tǒng)數(shù)據(jù)庫.
圖12 利用爬蟲技術爬取實測水位的開始頁面
(2)水域面積提取
SAR 影像和MSI 影像的水域面積提取過程類似,這里以SAR 影像提取水域面積過程為例予以說明.首先,設置調用方法的class 類名與Program()方法名,將.m 程序打包成.jar 文件,與JavaBuilder.jar 文件一起導入Java 工程.然后,調用接口訪問上傳至系統(tǒng)的SAR影像,若影像為空,則用戶需再次提交SAR 影像;非空時,先隨機初始化模糊矩陣并根據(jù)SAR 影像建立目標函數(shù),然后根據(jù)目標函數(shù)建立聚類中心,之后更新模糊矩陣并分割出水體輪廓,再統(tǒng)計水體輪廓和整個SAR 影像像素點個數(shù).最后,利用經緯度信息計算出SAR 影像水域面積.水陸分割后的SAR 影像數(shù)據(jù)示例如圖13所示,SAR 影像水域面積提取結果示例如圖14所示.
圖13 水陸分割后的SAR 影像數(shù)據(jù)示例
圖14 SAR 影像水域面積提取結果示例
(3)蓄水量估算
基于不同數(shù)據(jù)源數(shù)據(jù)的不同組合方式,系統(tǒng)提供了多種面積提取模型和蓄水量估算模型,方便用戶針對目標水庫實際情況對比分析和擇優(yōu)選出精度更高的蓄水量估算模型.對于系統(tǒng)提供的各種模型,用戶可以進行添加、修改和刪除等編輯操作,模型編輯界面如圖15所示.用戶對編輯界面中的各模型參數(shù)進行調整并保存后,系統(tǒng)會使用保存的各模型參數(shù)進行后續(xù)計算.各模型計算結果對比示例如圖16所示.
圖15 模型編輯界面
圖16 各模型計算結果對比示例
結果展示模塊調用數(shù)據(jù)庫中的數(shù)據(jù)并以圖表形式將其繪制在系統(tǒng)頁面上,實現(xiàn)“水位-日期”、“水位-水域面積”和“水位-蓄水量”等數(shù)據(jù)的可視化展示.系統(tǒng)讀取數(shù)據(jù)庫里的數(shù)據(jù)后,利用ECharts 可視化技術將數(shù)字化信息轉化成可視化的圖表信息,最后把轉化后的圖表信息呈現(xiàn)在系統(tǒng)的前端頁面上,如圖17所示.
圖17 結果展示頁面
本文以丹江口水庫為例進行系統(tǒng)實例測試,以從星載RA 數(shù)據(jù)提取的遙測水位、從星載MSI和SAR數(shù)據(jù)提取的水域面積,作為多源遙感水庫蓄水量估算模型的輸入?yún)?shù),將輸出的蓄水量與“全國水雨情信息網”上公布的標準蓄水量進行對比,結果如圖18所示,日期范圍是2019年1月1日至2019年12月31日.由圖18可知,本系統(tǒng)所使用的多源遙感水庫蓄水量估算模型得到的蓄水量和標準蓄水量的波動趨勢具有較高吻合度.在7~10月的汛期,丹江口水庫地區(qū)降雨較多,蓄水量增長迅速,在11月17日左右達到最高峰,然后蓄水量開始慢慢減少,逐漸進入枯水期,在5月初達到最低值.整體上,隨著旱季和汛季的轉換,丹江口水庫的蓄水量呈現(xiàn)規(guī)律性變化.
圖18 蓄水量估算結果對比圖
采用擬合優(yōu)度(R-square,R2)和均方根誤差(Root Mean Squared Error,RMSE)對多源遙感水庫蓄水量估算模型得到的蓄水量和標準蓄水量進行定量分析.分析結果顯示R2=0.99,說明該模型所得蓄水量與標準蓄水量具有很高的吻合度,能夠真實地反映水庫蓄水量的變化情況;此外RMSE=1.055,表明該模型估算的水庫蓄水量標準差在1.055 億m3以內,相對誤差控制在1%以內.
本文開發(fā)了一款基于多源遙感技術的水庫蓄水量估算系統(tǒng),在數(shù)據(jù)同化法的基礎上,通過建立多源遙感水庫蓄水量估算模型,實現(xiàn)水庫蓄水量估算功能.以丹江口水庫為例進行實驗,結果表明本文建立的估算模型所估算的蓄水量與官方公布的標準蓄水量具有很高的吻合度,能夠真實地反映水庫蓄水量的變化情況.該系統(tǒng)能夠以多源遙感方式獲取水位、水域面積、蓄水量等水庫參數(shù),可以很好地彌補目前單一遙感水庫蓄水量估算系統(tǒng)的不足,有望應用于防汛抗旱、水庫調度等領域,為管理部門科學化、高效化管理水庫提供技術支撐.
本文所提出的水庫蓄水量估算模型在建模階段僅使用了數(shù)據(jù)同化法,為了使蓄水量模型的估算結果更加精確,在未來的研究工作中,還可以嘗試使用水文公式法、積分法、均值法等其他估算方法,通過對比分析與實驗驗證,構建和選用更合適的蓄水量估算模型,以進一步提高蓄水量估算精度.同時,在系統(tǒng)的功能實現(xiàn)方面,還可以添加出/入庫徑流量監(jiān)測、水庫水質監(jiān)測等功能,為多角度、全方位、立體化的水庫監(jiān)測研究提供更多的技術支撐.