陸 瑤 ,Eric GAUME ,Olivier PAYRASTRE ,顏合想 ,信昆侖
(1.同濟(jì)大學(xué) 環(huán)境科學(xué)與工程學(xué)院,上海 200092;
2.法國交通及國土整治與網(wǎng)絡(luò)研究所 水與環(huán)境實(shí)驗(yàn)室,法國 南特 44344)
法國南部地中海沿岸是山洪災(zāi)害易發(fā)區(qū),據(jù)統(tǒng)計(jì),1988 — 2015年,該地區(qū)因山洪死亡人數(shù)超過250人[1]。2017年初,法國建成Vigicrues - Flash山洪預(yù)警系統(tǒng)[2],該系統(tǒng)利用雷達(dá)定量降水估測進(jìn)行水文模擬,根據(jù)模擬流量的重現(xiàn)期劃分洪水預(yù)警等級,一旦有險(xiǎn)情,系統(tǒng)將會以郵件、短信等方式向居民發(fā)送預(yù)警信息。Vigicrues - Flash系統(tǒng)已對法國本土超過3萬km的河道進(jìn)行實(shí)時監(jiān)測,包括缺乏水文資料的小流域。
為進(jìn)一步提高防洪能力,法國交通及國土整治與網(wǎng)絡(luò)研究所IFSTTAR提出了PICS項(xiàng)目[3],旨在開發(fā)一個山洪及其社會經(jīng)濟(jì)影響的即時預(yù)報(bào)預(yù)警系統(tǒng)。即時預(yù)報(bào)的含義為在山洪發(fā)生前幾小時對災(zāi)害及其影響做出精確預(yù)報(bào),影響的表現(xiàn)形式為受災(zāi)房屋數(shù)、被洪水切斷的道路等。相比于Vigicrues - Flash系統(tǒng),該系統(tǒng)將能夠?qū)ι胶楫a(chǎn)生的影響做出預(yù)報(bào),并依此劃分預(yù)警等級,使預(yù)警信息更加直觀,以方便救援、交通、市政等部門的人員使用。PICS項(xiàng)目已于2018年初啟動,歷時4 a完成。本文綜合PICS項(xiàng)目前期研究成果,開發(fā)一個山洪實(shí)時預(yù)報(bào)預(yù)警系統(tǒng)。
系統(tǒng)開發(fā)的理論依據(jù)為Le Bihan提出的“降雨 — 徑流 — 影響”預(yù)報(bào)鏈[4],流程見圖1。
圖1 “降雨—徑流—影響”預(yù)報(bào)鏈流程圖
預(yù)報(bào)鏈?zhǔn)褂美走_(dá)定量降水估測作為輸入數(shù)據(jù),調(diào)用分布式水文模型CINECAR進(jìn)行降雨 — 徑流模擬得到模擬流量。SHYREG水文頻率數(shù)據(jù)庫覆蓋全法國,包括缺乏水文資料的小流域,在該數(shù)據(jù)庫中可以找到與模擬流量對應(yīng)的洪水重現(xiàn)期。根據(jù)CARTINO方法,可以制作洪水淹沒范圍地圖,與用地類型圖重疊,即可得到各個重現(xiàn)期洪水淹沒的房屋數(shù)。Le Bihan運(yùn)用該方法對法國Alès地區(qū)2002年9月和Draguignan地區(qū)2010年6月的山洪事件進(jìn)行分析,模擬得到的受影響房屋數(shù)與保險(xiǎn)公司保單理賠數(shù)進(jìn)行比較,結(jié)果表明該方法是有效和可行的。
2015年10月3日晚間,暴雨襲擊位于法國南部地中海沿岸的Alpes - Maritimes省并引發(fā)山洪。共有20人在此事件中遇難,經(jīng)濟(jì)損失為5.5億 ~ 6.5億歐元。洪水主要發(fā)生在3個流域:Argentiere(47 km2)、Brague(68 km2)和Frayere(22 km2)。利用數(shù)字地形模型DTM提取流域河網(wǎng)并進(jìn)行子流域劃分,DTM文件的分辨率為75 m。Argentiere劃分為39個子流域,Brague劃分為53個子流域,F(xiàn)rayere劃分為10個子流域,每個子流域包含1條河段,河段長度為0.2 ~ 3.0 km(見圖2a)。本文將結(jié)合該事件闡述“降雨—徑流—影響”預(yù)報(bào)鏈預(yù)報(bào)過程。
圖2 Argentiere、Brague和Frayere數(shù)字流域和河網(wǎng)圖
3.1.1 CINECAR分布式水文模型
CINECAR分布式水文模型[5]是專為缺乏水文資料的小流域開發(fā)的。由于使用徑流曲線數(shù)法(SCS - CN)計(jì)算產(chǎn)流過程,只有1個參數(shù)CN需要率定,CN可用降雨事件開始前5 d的累計(jì)降雨量進(jìn)行率定。
Argentiere、Brague和Frayere流域無水文監(jiān)測站,為了驗(yàn)證水文模擬的可靠性,采用法國地中海地區(qū)2個水文站(Orbiel和Herault)的洪水過程進(jìn)行模型檢驗(yàn)(見圖3)。Orbiel水文站的模型檢驗(yàn)結(jié)果表明,峰現(xiàn)時間差為0.75 h,洪峰流量相對誤差為-26.04%,確定性系數(shù)為0.89;Herault水文站的模型檢驗(yàn)結(jié)果表明,峰現(xiàn)時間差為0.25 h,洪峰流量相對誤差為5.00%,確定性系數(shù)為0.87。結(jié)果表明,模型具有較好的適用性,可用于法國地中海地區(qū)無資料流域的水文模擬。
3.1.2 水文模擬
Brague流域的子流域A位于上游,子流域B位于下游(見圖2b)。圖4展示了2個子流域的水文模擬結(jié)果,對于子流域A,降雨量峰值發(fā)生在20∶30,降雨強(qiáng)度為116.59 mm/h;模擬流量峰值發(fā)生在21∶15,為74.04 m3/s。對于子流域B,降雨量峰值發(fā)生在20∶45,降雨強(qiáng)度為223.40 mm/h;模擬流量峰值發(fā)生在21∶30,為338.00 m3/s。
由于該地區(qū)無水文監(jiān)測站,為了驗(yàn)證模擬結(jié)果,研究人員在洪水事件后進(jìn)行實(shí)地考察,測量斷面形狀、坡度和最高水位留下的痕跡等數(shù)據(jù),并據(jù)此估算峰值流量區(qū)間(見圖4中“工”字型標(biāo)注)。對于子流域A,模擬洪峰流量落在測算區(qū)間中;對于子流域B,模擬洪峰流量稍大于測算區(qū)間。
圖4 Brague流域的2個子流域的模擬水文過程線圖
SHYREG水文頻率數(shù)據(jù)庫[6]用1個分辨率為50 m的格柵覆蓋法國本土,格柵文件中存儲著各單元格所在地區(qū)重現(xiàn)期為2,5,10,20,50,100,500,1 000 a的洪水流量。將數(shù)字河網(wǎng)與水文頻率格柵重疊,可提取每條河段所覆蓋的格柵單元格,取其中最小流量作為該河段的洪水流量數(shù)據(jù)。例如,表1為子流域B的水文頻率數(shù)據(jù)。完成降雨 —徑流模擬后,即可找出各時刻模擬流量對應(yīng)的重現(xiàn)期。
表1 子流域B水文頻率數(shù)據(jù)表
CARTINO方法利用一維水力學(xué)模型,根據(jù)流量值計(jì)算河流斷面水位,其中斷面情況提取自高分辨率的DTM文件,各斷面間的水位數(shù)據(jù)利用插值法估算,并根據(jù)DTM中的高程數(shù)據(jù)進(jìn)行合理性分析。由CARTINO方法,可制作各洪水重現(xiàn)期對應(yīng)的洪水淹沒范圍地圖。圖5為Brague流域重現(xiàn)期為10 a和1 000 a的洪水淹沒范圍地圖。把洪水淹沒范圍地圖與用地類型圖重疊,可數(shù)出受洪水影響的房屋數(shù)。表2為子流域B受影響房屋數(shù)。
圖5 Brague流域洪水淹沒范圍地圖用地類型圖重疊圖
表2 子流域B受影響房屋數(shù)表
根據(jù)洪水重現(xiàn)期和受影響房屋數(shù),可繪制影響地圖。河段的顏色表示該時刻的洪水重現(xiàn)期,實(shí)心圓點(diǎn)的大小表示該時刻河段周圍受影響的房屋數(shù)(見圖6)。根據(jù)受影響房屋數(shù),可劃分4個洪水預(yù)警等級:綠色(1 ~ 10幢房屋受影響)、黃色(11 ~ 20幢)、橙色(21 ~ 50幢)、紅色(>50幢)。從圖6可以看出,如果只以洪水重現(xiàn)期來判斷山洪災(zāi)害的嚴(yán)重性是不夠準(zhǔn)確的,例如20∶30,下游河段的洪水重現(xiàn)期僅為5 ~ 20 a,而受影響房屋數(shù)超過了50幢;21∶30,中上游河段洪水重現(xiàn)期都超過了100 a,而一些河段周圍卻沒有房屋受到影響。
圖6 Brague流域影響地圖
本文以“降雨 — 徑流 — 影響”預(yù)報(bào)鏈為基礎(chǔ),以Argentiere、Brague和Frayere流域?yàn)槔_發(fā)山洪實(shí)時預(yù)報(bào)預(yù)警系統(tǒng)。系統(tǒng)采用B/S結(jié)構(gòu),服務(wù)器端運(yùn)行的主程序?qū)崟r下載雷達(dá)定量降水估測,完成“降雨 — 徑流 — 影響”模擬,得到影響數(shù)據(jù),并以Excel表格的形式存儲在服務(wù)器中;用戶從瀏覽器端訪問網(wǎng)頁,讀取Excel表格中的數(shù)據(jù),制作影響地圖并展示給用戶。此外,用戶還可通過網(wǎng)頁查看歷史洪水事件的影響地圖,系統(tǒng)流程見圖7。系統(tǒng)主程序用R語言編寫,交互式網(wǎng)頁用R語言Shiny包搭建,影響地圖用R語言Leaflet包制作。
圖7 系統(tǒng)流程示意圖
系統(tǒng)所用的實(shí)時降水?dāng)?shù)據(jù)來自于法國氣象局的雷達(dá)定量降水估測產(chǎn)品MOSAIQUE NATIONALE,數(shù)據(jù)格式為 BUFR(Binary Universal Form for the Representation of Meteorological Data),雷達(dá)圖像覆蓋全法國,空間分辨率為1 km,時間分辨率為5.0 min,延遲時間為2.5 min,例如,可在 18∶37∶30 拿到 18∶30∶00 — 18∶35∶00 的 5.0 min 累計(jì)降水量數(shù)據(jù)。
法國氣象局實(shí)時上傳降水?dāng)?shù)據(jù)至服務(wù)器,以供系統(tǒng)下載使用。原始降水?dāng)?shù)據(jù)為BUFR格式的二進(jìn)制數(shù)據(jù),因此在進(jìn)行水文模擬前需要對數(shù)據(jù)進(jìn)行預(yù)處理,將其轉(zhuǎn)換為格柵格式的數(shù)據(jù),以便后續(xù)程序直接讀取。在進(jìn)行降雨 — 徑流模擬前,需要計(jì)算每個子流域的降雨強(qiáng)度過程線。將雷達(dá)格柵與數(shù)字流域重疊,可計(jì)算每個子流域落在每個雷達(dá)格柵單元格的面積百分比,取面積百分比加權(quán)值作為子流域的降雨強(qiáng)度。
圖8 雷達(dá)定量降水估測圖(5.0 min累計(jì)降水量2015 -10 - 03 - 18∶30—18∶35)
系統(tǒng)的主程序用R語言編寫,負(fù)責(zé)“降雨 — 徑流 —影響”預(yù)報(bào)鏈中的數(shù)據(jù)處理及傳遞。模擬使用的CINECAR分布式水文模型用C++編寫,可在R主程序中調(diào)用該C++程序完成水文模擬。“降雨 — 徑流 — 影響”模擬每5.0 min進(jìn)行1次,與雷達(dá)定量降水估測的時間分辨率一致,以保證影響數(shù)據(jù)是最新的。
Shiny是R語言中的一種開發(fā)框架,使用戶可以不必了解網(wǎng)頁編程語言(如html、CSS、JS)也可以快速搭建動態(tài)交互、即時更新的網(wǎng)頁,以完成對數(shù)據(jù)的可視化。本系統(tǒng)的Shiny網(wǎng)頁主要包括2個功能,一為展示實(shí)時影響地圖;二為展示歷史山洪事件影響地圖。用戶可拖動時間進(jìn)度條以查看不同時間點(diǎn)的影響地圖。在查看歷史山洪事件功能中,用戶可以在列表中選擇想要查看的事件。
圖9 Shiny網(wǎng)頁應(yīng)用用戶界面圖
Leaflet包是R語言中的一個地圖集,可用于制作動態(tài)交互地圖,該地圖能夠嵌入Shiny網(wǎng)頁中。本系統(tǒng)的影響地圖即是用Leaflet繪制的,其底圖來自于CartoDB基礎(chǔ)地圖,河網(wǎng)、洪水重現(xiàn)期和受影響房屋數(shù)通過加入線圖層、面圖層,并進(jìn)行顏色映射操作的方式展現(xiàn)。
本文開發(fā)的系統(tǒng)針對山區(qū)小流域突發(fā)性洪水預(yù)見期短的特點(diǎn),使用高分辨率的雷達(dá)定量降水估測進(jìn)行“降雨 —徑流”模擬,實(shí)現(xiàn)山洪的即時預(yù)報(bào);系統(tǒng)還進(jìn)一步將受災(zāi)房屋數(shù)納入到預(yù)報(bào)內(nèi)容中,實(shí)現(xiàn)對山洪產(chǎn)生的社會經(jīng)濟(jì)影響的預(yù)報(bào)。系統(tǒng)將用于測試PICS項(xiàng)目開發(fā)的各個組件,不斷進(jìn)行完善,并最終應(yīng)用于法國地中海沿岸地區(qū)的山洪預(yù)警,為當(dāng)?shù)胤篮楣ぷ魈峁┘夹g(shù)支持。本系統(tǒng)對我國的山洪災(zāi)害預(yù)報(bào)預(yù)警研究具有參考價(jià)值:①提高水文資料信息化程度,加強(qiáng)中小河流水文頻率數(shù)據(jù)庫建設(shè);②針對山洪災(zāi)害高發(fā)地區(qū),完善洪水淹沒范圍地圖的編制工作;③在山洪災(zāi)害預(yù)警預(yù)報(bào)系統(tǒng)開發(fā)方面,應(yīng)考慮到非水文專業(yè)用戶的使用需求,更加直觀地展示水文預(yù)報(bào)預(yù)警信息。