基于多源遥感数据的蓄洪区洪涝遥感监测与分析

朱 辉，何祺胜，李金阳，陈 丽

(1.河海大学水文水资源学院，江苏 南京 210098；2. 江苏省水利科学研究院，江苏 南京 210017)

我国地域辽阔、自然地理和气候条件复杂、各地自然灾害频繁发生，并且在全球变暖以及人为干扰的影响下，极端降水频发[1]，导致洪涝灾害的发生。洪涝灾害被认为是最具破坏性、最广泛和最频繁的灾害之一[2-3]，是自然灾害链中造成巨大损失的重要一环[4-6]，已成为制约经济社会持续健康发展的突出瓶颈[7]。为了科学有效地应对洪涝灾害，降低灾害造成的经济财产损失，快速准确地提取洪涝灾害范围，尤其是近实时、高频次的洪涝动态监测制图，以及灾后对洪涝灾害受灾区域的时空变化特征进行详细准确的分析至关重要，是洪涝灾害治理领域研究的重要内容。

由于遥感技术的时效性和低成本性，卫星数据在洪涝灾害淹没范围提取中发挥了重要的作用[8-11]。汤玲英等[12]基于Sentinel-1A SAR数据，采用面向对象的方法提取广西壮族自治区桂林市临桂县会仙岩溶湿地区域的洪涝灾害淹没范围，绘制灾前、灾中、灾后水体变化监测图。苏亚丽等[13]利用多源卫星遥感数据，采用一种基于NDVI变化特征值的方法对安徽省巢湖地区进行洪涝灾害动态信息的提取。Liu等[14]采用时空上下文学习的方法和Modest AdaBoost分类器，实现了在HJ-1A和GF-4影像上提取洪涝范围。Wan等[15]根据SAR影像，采用结合阈值自动选择、面向像元和面向对象分类以及双向区域增长的一种洪水自动提取混合方法对吉林省的洪涝灾害淹没范围进行提取。Abdelhakim等[16]采用SAR和sentinel-2训练监督分类器，实现洪水范围的自动提取。Son等[17]利用含水率模型(WFM)，在MODIS影像中提取越南南方的洪涝灾害动态淹没范围。Xie等[18]引入大津算法，实现基于Landsat 8影像的洪水无监督自动提取。高伟等[19]利用MODIS数据和达特茅斯洪水数据库，通过数据协同的方法对举水流域洪水范围进行时间序列分析。李思慧[20]基于GF-1和MODIS遥感数据，结合地面气象观测数据和基础地理信息数据，对通辽市洪水范围进行提取及分析。上述研究注重于洪涝灾害信息提取的精确性，在洪涝灾害信息的时效性和高频次方面有所欠缺，而时效性和高频次对于快速制定响应计划和减小损失非常重要。在时效性方面，需要考虑信息提取方法的效率问题；在高频次方面，由于卫星的重访周期以及洪涝灾害期间云量较大，光学影像质量差，需要考虑数据源的选取问题。

蓄滞洪区是指河堤外洪水临时贮存的低洼地区及湖泊等，包括行洪区、分洪区、蓄洪区和滞洪区，是流域防洪工程体系中的重要组成部分，是减轻灾害损失、保障人民生命财产安全的重大工程。如鄱阳湖区设有康山、珠湖、黄湖、方洲斜塘4座蓄滞洪区，极大地提升了鄱阳湖区防洪体系的灾害防御能力[21]；大名泛区是承担着漳卫河系滞沥和漳河超量洪水分洪任务的蓄滞洪区[22]；位于安徽省阜阳市的蒙洼蓄洪区是淮河流域极其重要的防汛工程，从1954—2020年期间已开闸蓄洪16次，在淮河流域抗洪救灾中具有不可替代的地位。开闸蓄洪对蒙洼蓄洪区造成相当大的经济损失，据统计，2003年蒙洼蓄洪区开闸蓄洪直接经济损失达3.8亿元，2007年直接经济损失达6亿元，2020年7—8月直接经济损失达6.7亿元。

基于此，本文借助遥感和GIS技术，以蒙洼蓄洪区为研究区，选取哨兵一号、二号和高分系列卫星影像作为数据源，采用区域生长算法对覆盖研究区的9景高质量遥感影像进行洪涝信息快速提取，分析洪涝期间蒙洼蓄洪区受灾区域的时空变化，并结合研究区的土地利用类型覆盖图，分析蒙洼蓄洪区各土地利用类型受灾情况的时空变化，以期为蓄洪期间减灾救灾、降低蒙洼蓄洪区经济损失等提供参考。

1 研究区概况

蒙洼蓄洪区位于安徽省阜阳市阜南县境内(115°35′20″E～115°56′38″E、32°25′1″N～32°36′47″N)，淮河干流洪河口至南照集间，南临淮河，北临蒙河分洪道，汛期四面环水[23]。蒙洼蓄洪区东西长约40 km，南北宽2～10 km，呈西南—东北走向的狭长地带，地势由上游王家坝向下游曹台孜倾斜，由淮堤渐向蒙河分洪道倾斜[24]。其总面积约180.4 km2，设计蓄洪量7.5亿m3，涉及4个乡镇，75个行政村、131座庄台、17万余人，耕地面积约13 160 hm2。图1为研究区Landsat 8 OLI 7、3、2波段遥感影像。

图1 蒙洼蓄洪区Landsat 8 OLI遥感影像(7、3、2波段)Fig.1 Landsat 8 OLI remote sensing image of Mengwa flood storage area (7, 3 and 2 bands)

2020年7月20日，由于淮河上游汛情紧张，王家坝闸开闸蓄洪，蒙洼蓄洪区第16次启用。7月23日，王家坝闸关闸。8月1日，经过12 d的蓄滞洪后，蒙洼蓄洪区退水闸正式开闸退水。

2 研究数据与研究方法

2.1 数据源

选取灾前、灾中质量较高的哨兵一号(Sentinel-1A)、二号(Sentinel-2A)数据和高分系列卫星(GF-1、GF-3和GF-4)2020年7—8月数据进行洪涝灾害淹没范围的提取，包括7月8日(灾前)、7月20日、7月21日、7月26日、8月1日、8月4日、8月17日、8月18日和8月24日的影像，各影像数据的成像时间及级别见表1，其他辅助数据包括Google高分辨率影像和研究区矢量边界。

表1 遥感卫星2020年7—8月数据信息

Sentinel-1A卫星于2014年4月3日发射，是欧洲空间局哥白尼计划发射的首颗环境监测卫星，在近极地太阳同步轨道上运行，轨道高度约700 km，重访周期为12 d。它搭载一台C波段合成孔径雷达，成像系统采用4种成像模式，分别为WV(wave mode)、EW(extra-wide swath mode)、IW(interferometric wide swath)和SW(strip map mode)，极化方式为VV/VH，其中WV模式为单极化成像，EW、IW和SW模式均分为单极化和双极化成像。

Sentinel-2A卫星是2015年发射的一颗新星高分辨率多光谱成像卫星，在可见光、近红外(VNIR)和短波红外(SWIR)光谱波段上设置有13个光谱波段(4个10 m波段，6个20 m波段，3个60 m波段)，主要用于陆地植被、土壤及水资源等在内的全球陆地观测[25]。Sentinel-2A卫星搭载一台多光谱成像仪，幅宽为290 km，辐射分辨率为12 bit，是第一个包含3条红边波段的陆地观测光学卫星。

GF-1卫星是2013年发射的我国第一颗高分辨率对地观测卫星，搭载了2台高分辨率相机，分别是1台2 m空间分辨率的全色相机、1台8 m空间分辨率的多光谱相机，以及4台16 m空间分辨率的中低分辨率多光谱宽幅相机。全色相机的重访周期为4 d，多光谱相机的重访周期为2 d，其中全色相机的波长范围为0.45～0.90 μm，包含蓝色波段(0.45～0.52 μm)、绿色波段(0.52～0.59 μm)、红色波段(0.63～0.69 μm)和近红外波段(0.77～0.89 μm)，幅宽大于60 km；多光谱宽幅相机幅宽大于800 km，因此GF-1影像具有探测范围大、空间分辨率高的优点。

GF-3卫星是我国于2016年8月10日发射的首颗全极化雷达卫星，轨道高度为755 km，轨道类型为太阳同步回归晨昏轨道，天线类型为波导缝隙相控阵，重访周期为 2～3 d[26]。GF-3卫星具有12种常规成像工作模式，4种在轨试验工作模式，分辨率为1～500 m，成像幅宽10～650 km，具有单极化、双极化、全极化成像能力，最长单次成像时间50 min。

GF-4卫星是我国于2015年12月29日发射的第一颗地球同步轨道对地观测遥感卫星，搭载了一台可见光和近红外50 m/中波红外400 m分辨率、幅宽大于400 km的凝视相机，通过面阵凝视方式成像，具备可见光、近红外和中红外同时成像的能力[27]。

2.2 数据预处理

Sentinel-1 GRD产品是经过多视和地理编码矫正处理的SAR数据，无需进行地理配准。使用欧洲航天局(European Space Agency)研发的SNAP(Sentinel Application Platform)软件对其进行预处理，主要流程包括辐射校正、几何地形校正、去除噪声和后向散射系数影像转换，并重采样为30 m分辨率。

Sentinel-2 L1C产品是已经过系统几何精校正的光学影像，无需再进行几何精校正。使用SNAP软件对其进行预处理，调用SNAP软件的Sen2Cor模块进行大气校正，然后对校正好的数据进行重采样。

GF-1和GF-4影像属于光学影像，以Sentinel-2为基准进行几何校正，误差控制在1个像素以内。GF-3影像属于雷达影像，使用航天宏图信息技术股份有限公司研发的PIE-SAR(Pixel Information Expert for SAR)软件进行预处理，主要流程包括辐射定标、极化滤波、正射校正和地理编码，最后利用ArcMap将影像重采样为30 m分辨率。

2.3 土地利用分类

为了得到蒙洼蓄洪区土地利用类型分布图，用于分析研究区内各土地利用类型受灾情况及其时空变化特征，同时考虑到蒙洼蓄洪区面积较小，人类活动因素对其土地利用类型分布影响较大，故基于谷歌地球引擎(Google Earth Engine)云平台，参照2020年4月的Google影像进行目视解译，以确保解译结果的准确性。蒙洼蓄洪区的土地利用类型分为农田、养殖塘、庄台、保庄圩、河流、道路和其他，其土地利用类型分布图如图2所示。

图2 蒙洼蓄洪区土地利用类型Fig.2 Land-use type of Mengwa flood storage area

2.4 水体提取

利用遥感技术获取大范围水体的空间信息具有很大的优势，目前水体提取较为常用的方法主要有单波段阈值法、水体指数法、谱间关系法、监督分类法以及区域生长法。①单波段阈值法是利用水体在近红外波段反射率较低、易与其他地物区分的特点，根据确定的分类阈值进行水体提取[28]。②水体指数法是根据水体和其他地物类型的光谱特征，选取相应波段进行代数运算，从而在影像上突出水体信息，通过阈值分析提取水体。该方法常用的水体指数有归一化差异水体指数(NDWI)、改进的归一化差异水体指数(MNDWI)和修订型归一化差异水体指数(RNDWI)等。③谱间关系法是通过分析地物的光谱曲线特征，利用逻辑判别运算进行水体提取[29]。上述3种方法操作流程简单，但易受到阴影、建筑物等地物的影响，且水体指数法和谱间关系法不能用于雷达影像。④监督分类法是根据选取的训练区样本，建立和训练分类器，利用分类器进行图像分类。该方法虽然精度较高，但需要选取样本，耗时过长，且同样不能用于雷达影像。⑤区域生长法是目前应用较为广泛的图像分割方法之一，思路简单，只需提供若干种子点即可完成，能够提供良好的边界信息，产生较好的分割结果及较快的分割速度[30-31]。该方法侧重于集中连片水体的提取，适合洪涝水体这类大范围水体的快速提取[32]。

本文选择区域生长法进行水体提取。区域生长法适合于光学和雷达影像，其基本思想是根据事先定义的准则将像素点聚合在一起，即以一组种子点作为生长的起点，根据定义的准则，将种子点周围的像素点与种子点进行对比，合并具有相似性质(光谱特征或纹理特征等)的像素点，并继续向外生长，直到没有满足同质准则的相邻像素点为止，然后进行区域合并以减少区域的数量[33]。

选择合适的种子点是区域生长法水体提取精度的关键因素，合适的种子点应该符合以下条件：①必须位于水体之内，最好位于水体的中央区域；②应选择在光谱反射或后向散射较为均匀的水体区域内，并且种子点的像素值在该区域内具有一定的代表性，否则会导致聚类结果的区域过小，影响水体提取的效率。若影像上水体区域的光谱反射或后向散射不均匀，需要将影像进行分块后再进行水体提取。

合适的阈值是水体提取的另一关键因素。区域生长法的阈值代表在同一景影像中用于判别像素点是否同质的光谱信息或后向散射信息的差值。阈值过小将导致聚类结果的区域太小，影响水体提取的效率；阈值过大将导致其他地物被划分在水体之中，降低水体提取的精度。图3(a)～(c)分别为基于2020年8月1日Sentinel-1雷达影像，选择同一种子点(红色点)，阈值设为10、20和40时，区域生长法提取的水体(图中蓝色边框)，图3(d)(e)分别为基于2020年8月4日GF-1光学影像，选择同一种子点，阈值设为5、10和15时，区域生长法提取的水体。由图3(a)～(c)可以看出，对于Sentienl-1影像，20是水体提取较为合适的阈值，阈值为10时聚类结果的区域偏小，阈值为40时其他地物被划分在水体之中，提取精度不高；而对于图3(d)(e)的GF-1影像，10是较为合适的阈值。

图3 区域生长法水体提取结果Fig.3 Water extraction results of region growing algorithm

采用区域生长法初步提取水体之后，参照对应影像，通过目视解译进行验证，并使用ArcGIS软件对提取的水体进行人工修改，得到高精度的水体矢量数据。

3 结果与分析

3.1 水体提取结果

采用区域生长法提取水体，并结合目视解译进行人工修改，得到研究区的水体分布数据。图4(a)～ (j)分别为蒙洼蓄洪区2020年7月8日、7月20日、7月21日、7月23日、7月26日、8月1日、8月4日、8月17日、8月18日和8月24日的水体范围。其中7月23日缺少高质量的遥感影像提取水体，但由于该日王家坝闸关闸，理论上7月23日是受灾最严重的一天，因此根据7月21日和26日的水体，使用ArcGIS进行叠加操作，以此作为7月23日当天的水体。

图4 蒙洼蓄洪区水体提取分布Fig.4 Extraction distribution of water body in Mengwa flood storage area

3.2 受灾区域的空间变化特征

根据提取的蒙洼蓄洪区水体分布数据，分析受灾区域的空间变化特征。使用ArcGIS软件的空间分析模块，对除7月8日(灾前)外的所有水体矢量数据进行空间分析，得到蒙洼蓄洪区各个区域淹没持续时间分布图，以此判断各个区域的受灾程度。蒙洼蓄洪区淹没持续时间分布如图5所示。

图5 蒙洼蓄洪区淹没时间分布Fig.5 Distribution of inundation time in Mengwa flood storage area

根据研究区的最大淹没历时，设淹没历时大于20 d的区域为重度受灾区，5～20 d为中度受灾区，小于5 d为轻度受灾区。从图5可以看出，蒙洼蓄洪区中部和东北部地区受灾情况较为严重，淹没历时均达20 d以上，为重度受灾区；轻度受灾区和中度受灾区大部分位于蒙洼蓄洪区的西南部地区；无受灾区，即未淹没的区域，大部分分布在蒙洼蓄洪区的圈堤附近，少部分分布比较广泛，无规律。

结合蒙洼蓄洪区的地域特征和土地利用类型分布图，分析蓄洪区受灾区域空间变化特征的原因为：①蒙洼蓄洪区的地势由上游王家坝向下游曹台孜倾斜，西南高，东北低，蓄洪期间洪水汇聚在地势低的区域，因此蓄洪区轻度受灾区和中度受灾区大部分位于西南部，重度受灾区位于中部和东北部；②保庄圩是由新建圩堤和淮河大堤连接起来形成的圩区，庄台是截面为正梯形的土台，地势较高，皆可有效预防洪涝灾害。无受灾区的地理位置与蓄洪区内的保庄圩和庄台的分布基本一致，因此蒙洼蓄洪区内的无受灾区基本是保庄圩和庄台。

3.3 受灾区域的时间变化特征

根据提取的蒙洼蓄洪区水体分布数据，统计出蒙洼蓄洪区2020年7月20日、7月21日、7月23日、7月26日、8月1日、8月4日、8月17日、8月18日和8月24日的淹没面积分别为38.08 km2、102.33 km2、150.29 km2、133.09 km2、125.02 km2、119.61 km2、91.80 km2、76.99 km2和40.14 km2。可以看出，受灾期间以7月23日的淹没面积最大，20日水体面积最小； 23日之前淹没面积呈上升趋势，23日之后呈下降趋势。

综上，总结研究区受灾区域时间变化特征为：①7月20日，西南部最先受灾，并逐渐向东北部延伸，淹没面积增长明显，至7月23日达到最大值；②7月26日，西南部开始退洪，总体淹没面积逐渐下降，下降趋势不大；③8月4日开始，淹没面积明显下降，至8月18日，淹没面积下降至最大淹没面积的一半以下。

蒙洼蓄洪区启用各措施的时间是7月20日开始蓄洪、23日王家坝闸关闸、8月1日开始退洪，上述淹没面积时间变化特征与蒙洼蓄洪区开展蓄洪措施的时间基本吻合。

3.4 各土地利用类型受灾情况的时空变化特征

根据图2，使用ArcGIS计算出各土地利用类型的面积为：农田120.7 km2，养殖塘37.7 km2，保庄圩9.6 km2，庄台5.7 km2，河流2.1 km2，道路1.4 km2，其他6.1 km2。可以看出，农田是蒙洼蓄洪区面积最大的土地利用类型，养殖塘其次，而其余土地利用类型的面积相对于农田和养殖塘都较小，蓄洪前后其受灾面积不会发生明显变化。此外，农田和养殖塘是蒙洼蓄洪区的主要经济来源，其洪涝灾害受灾情况与造成的经济损失密切相关，而保庄圩和庄台不会受洪涝灾害的影响，因此主要分析农田和养殖塘这两种土地利用类型的淹没面积时间变化。

表2为2020年7月8日、7月20日、7月21日、7月23日、7月26日、8月1日、8月4日、8月17日、8月18日和8月24日蒙洼蓄洪区内农田、养殖塘淹没面积的时间变化。由表2可以看出，蒙洼蓄洪区农田和养殖塘的淹没面积在7月23日达到最大值，分别为109.01 km2和35.30 km2，23日之前淹没面积呈上升趋势，23日之后淹没面积呈下降趋势；7月20日至8月18日农田和养殖塘的平均淹没面积分别为74.77 km2和24.19 km2，占各自总面积的61.95%和64.16%，受灾情况均较为严重。

表2 蒙洼蓄洪区农田、养殖塘淹没面积时间变化

由图2可以看出，农田和养殖塘在蓄洪区内的分布相对均匀，因此其淹没面积的时间变化特征与总区域淹没面积的时间变化特征基本一致。此外，农田在蓄洪区的东北部地区分布较多，养殖塘在蓄洪区的中部地区分布较多，而东北部和中部地区均为重度受灾区，因此农田和养殖塘受灾情况均较为严重。

4 结 论

a.区域生长算法能够高效精确地识别大面积连片的水体，并且适用于光学和雷达影像；该算法结合多源卫星遥感数据能够实现洪涝灾害动态监测的时效性和高频次。

b.2020年7月20日至8月18日蒙洼蓄洪区平均受灾面积约105 km2。蒙洼蓄洪区中部区域为受灾面积最大的区域，其次是东北部区域，最小是西南部区域。基于蒙洼蓄洪区西南部地势高、东北部地势低的地域特征，结合蒙洼蓄洪区受灾区域的空间变化特征，风险较高的区域主要位于蒙洼蓄洪区的中部和东北部；蒙洼蓄洪区内的农田和养殖塘受灾均较为严重，原因在于农田在蓄洪区的东北部分布较多，养殖塘在蓄洪区的中部分布较多，均为高风险、受灾严重地区。

c.蒙洼蓄洪区受灾情况的时间变化特征与蓄洪措施存在有很强的相关性。洪涝灾害期间连续降雨，云量较多，7月23日没有高质量的光学影像，而由于重访周期的原因同样没有相应雷达影像，因此7月23日的水体覆盖数据是叠加合成的数据，在后续研究中可增加更多的卫星影像，确保数据的获取。