基于Sentinel-1/2数据的洪水淹没范围提取模型

摘要：遥感是监测洪水淹没范围、掌握洪涝灾情演变的重要手段，而光学影像在洪水发生时往往有较多缺失，全天候的SAR影像在提取水体时精度略低。为快速、精准提取洪水淹没范围，构建了一种综合利用Sentinel-2光学影像和Sentinel-1雷达影像数据的洪水淹没范围提取模型，采用一种自适应阈值分割算法即大津算法（OTSU）分别对两种数据以及该模型进行了水体范围提取试验，并以河北省保定市为例进行了应用分析。结果显示：云量较少的Sentinel-2影像水体提取效果最好，总体精度（OA）达到95.6%；所构建的模型在引入部分可用Sentinel-2数据后，OA达到95%，相比单独使用Sentinel-1数据OA和Kappa系数分别提升1.2%和2.4%。该模型搭载于Google Earth Engine平台，能实现快速、准确、低成本的地表水体空间范围连续输出，不受限于云雾且比单独使用Sentinel-1影像的提取精度更高，在云覆盖严重导致Sentinel-2数据缺少的情况下，该模型可作为洪水淹没范围提取方法的一种选择。

关 键 词：洪水淹没范围； Sentinel-1； Sentinel-2； 自适应阈值分割算法； Google Earth Engine； 保定市

中图法分类号： P954

文献标志码： A

DOI：10.16232/j.cnki.1001-4179.2024.09.010

0 引 言

由于地理位置和东亚季风的原因，中国一直饱受洪涝灾害的侵扰。据《中国水旱灾害防御公报2022》概要中统计［1］，2022年因洪涝灾害造成171人死亡失踪，直接经济损失达1 288.99亿元（占当年自然灾害总直接经济损失的54.01%）。及时、大范围地掌握洪涝灾情的演变情况，对灾情评估、灾后重建有重要的指导意义。

基于水的光谱特性，方便、快速的水体指数模型法可以在不同复杂情景的光学影像中取得较高的水体识别精度［2-4］，可运用于水体的动态监测［5］。但受限于云雾，洪灾发生时段实际可用的光学影像较少，因此在监测洪水面积时很难将光学影像作为唯一数据来源。SAR影像具有大范围、全天候、短周期、长时序的优点，不仅被广泛应用于城市矿山地面沉降［6］，在洪水动态连续监测方面也具备独特的优势。基于上述背景，国内外学者对洪水提取模型进行了广泛研究。例如贾诗超等［7］分析微波遥感水体信息的影像特点后，提出SDWI指数模型增强水体特征，提取效果显著；Amitrano等［8］提出了一种基于Sentinel-1数据的无监督洪水提取模型，能够快速提取洪水发生范围；Binh等［9］结合SAR数据（Sentinel-1）和光学数据（Landsat-8）的训练数据开发了一种神经网络分类模型，并运用于湄公河三角洲，证明该模型有较强鲁棒性和较高准确率。然而上述模型也有一定缺陷：Amitrano等所建模型虽然可以快速提取洪水发生范围，但数据源单一，精度较低；Binh等所建模型虽取得了较高的精度，但训练模型需要花费较长时间和较大人力，很难运用于强调时效性的抗洪救灾场景中。与此同时，数据的获取和必要的预处理也要尽可能方便，以提高计算效率，而Google Earth Engine（GEE）平台正好满足这一需求。GEE是向全球用户开放，可在线处理全球尺度卫星遥感数据并进行计算和可视化分析的平台，通过编码即可远程调用Landsat系列、Sentinel系列等免费数据进行处理而不占用计算机资源，极大提高了用户的效率，在水体监测方面得到了广泛应用［10-12］。

综上所述，本研究考虑综合利用光学（Sentinel-2）和SAR（Sentinel-1）两类数据，理由如下：一方面，覆盖范围广、全天候的Sentinel-1数据可以补足Sentinel-2数据的缺失；另一方面，部分存有的Sentinel-2数据可以提升模型的精度。在方法上选择不依赖前期训练、能够快速得到水体范围的自适应阈值分割法，搭载于GEE平台来满足洪水监测的时效性要求。最后，通过实验验证单独使用Sentinel-1/2数据以及两者综合利用构建的本文模型在不同场景的适用性和准确度并进行对比分析。研究成果以期实现洪水范围的长时序、全自动、高精度提取。

1 研究区概况及数据预处理

1.1 研究区概况

研究区位于河北省保定市东部，包含6个县级行政区，地理位置约在115°30′E～116°15′E，38°42′N～39°36′N之间。区域内地表永久性水体主要为位于涿州市的北拒马河与位于安新县的白洋淀，除了建成区以外，境内大部分土地类型为耕地，通过实地走访调查得知，种植的主要作物为玉米、小麦。2023年受到台风“杜苏芮”等多重因素影响，河北省多地发生特大暴雨从而造成洪水，研究区作为受洪灾影响最严重的区域，引起了社会广泛关注。

1.2 数据预处理

本研究SAR影像数据采用Sentinel-1A数据，数据分辨率为10 m，按工作模式为干涉宽幅模式（IW），数据类型为地距多视（GRD），极化方式为VV和VH双极化，成像轨道为升轨的条件筛选影像。随后按研究区范围对影像进行裁剪并用Refined Lee滤波［13］消除一定的SAR影像相干斑噪声。光学影像数据使用Sentinel-2数据，由多光谱成像仪（MSI）成像，包含可见光、近红外、短波红外等13个波段，分辨率包含10，20，60 m三种，本文选用10 m分辨率以获得更高的提取精度。光学影像去云处理采用GEE平台自带的云量分数S2_CLOUD_PROBABILITY数据集，该数据集通过s2cloudless算法计算每一景Sentinel-2影像像元的云量分数，像元云量分数越大，含云的概率越大。本研究所有Sentinel-2影像数据均按照云量分数小于50的准则进行去云处理。经过预处理，得到覆盖整个研究区的多幅Sentinel-1/2影像，每幅影像的数据类型及其对应的时间信息如表1所列。最后，SAR侧视成像还存在山体阴影和水体相互混淆的问题，因此利用美国航天飞机雷达地形测绘计划SRTM（Shuttle Radar Topography Mission）DEM数据计算坡度栅格［14］，设置坡度阈值对Sentinel-1数据水体范围进行掩膜以去除山体阴影的影响。

2 原理与方法

阈值法分割图像的关键在于阈值的准确度。在地表水体提取的情景中，当时间跨度较大时，同一地理位置的大气条件、地面作物的生长阶段等不同会引起最佳分割阈值的变化，此时利用固定的阈值分割会导致较为严重的误分。大津算法（OTSU）作为阈值算法的一种，具备计算快速、准确率较高的优点，相比于直方图双峰法受限于明显的双峰特性，其在阈值的选取上脱离人工分析，能够满足大批量遥感影像洪水淹没范围自动提取的需求。

2.1 自适应阈值分割原理

大津算法，又被称为最大类间方差法［15］，通过遍历图像的灰度值，寻找使得图像的前景和背景的类间方差最大的阈值，从而实现二值化阈值分割。其具体计算方式如下：

设图像的灰度空间为［0，T］，某个灰度级ti的像素个数为f（ti），N为图像像素的个数，则该灰度级出现的概率Pi为

Pi=ftiN（1）

设灰度级t为当前阈值，图像被分为A、B两类，两类像素的概率PA、PB及均值灰度MA、MB分别为

PA=ti=0Pi（2）

PB=Ti=tPi=1-PA（3）

MA=ti=0i（Pi/PA）（4）

MB=Ti=ti（Pi/PB）（5）

整幅影像的均值灰度M为

M=Ti=0tiPi=MAPA+MBPB（6）

两类的类间方差δ2为

δ2=PAMA-M2+PBMB-M2

=PAPBMA-MB2（7）

遍历图像的灰度空间使得δ2最大的灰度值即为需要求得的最佳阈值。

2.2 Sentinel-1数据的水体提取方法

SAR卫星传感器发射的微波对于地表目标的散射特性为：在粗糙表面（建筑、裸土等）的反射为漫反射，回波信号强，后向散射系数大，在影像中呈现亮白色；在光滑表面（水体等）的反射为镜面反射，回波信号弱，后向散射系数较小，在影像中呈现为暗黑色。因此可以选择合适的阈值来提取SAR影像中的水体。

Sentinel-1A影像由VV和VH两种极化方式成像，在双重反弹效应（电磁波在水面散射到部分被淹没的树干等处，然后直接反射回SAR传感器）下，VV极化相比VH极化对表面粗糙度更为敏感，并显示出更强烈的信号响应，有利于提取洪水中的半淹没区域，具有更好的洪水监测适用性［16］。因此，本研究采用OTSU算法自动提取VV极化SAR图像的分割阈值完成水体提取，即VV-OTSU的Sentinel-1A影像水体自动提取方法。

2.3 Sentinel-2数据的水体提取方法

光学影像的水体指数法是水体提取的一种常用方法。其原理是基于水体对绿光波段的高反射率以及近红外光波段的高吸收率特性，通过比值运算构建水体指数来增强水体在影像中的特征，而后使用简单的分割算法完成水体范围的提取。

水体指数中NDWI最为常见和流行［17］，其计算公式如下：

NDWI=（ρ3-ρ8）/（ρ3+ρ8）（8）

式中：ρ3、ρ8分别为Sentinel-2数据中绿光波段B3和近红外波段B8的反射率值。

Xu［18］通过研究改进NDWI指数提出了MNDWI，可以在水、建筑物以及植被混合的区域取得更高的水体提取精度。由于短波红外波段中水的吸收比近红外波段强，因此使用短波红外波段代替近红外波段。MNDWI的计算公式如下：

MNDWI=（ρ3-ρ11）/（ρ3+ρ11）（9）

式中：ρ3、ρ11分别为Sentinel-2数据中绿光波段B3和短波红外波段B11的反射率值。

需要注意的是，Sentinel-2数据中B3、B8波段的分辨率为10 m，而B11波段的分辨率为20 m，在计算MNDWI之前先将B11波段使用最邻近重采样至10 m。为了比较NDWI和MNDWI在研究区内的适用性，在同一时期的两个不同地点的Sentinel-2 NDWI/MNDWI图像上分别使用OTSU阈值对水体提取范围进行比较。图1显示了两者在典型区域的视觉结果，可以看出NDWI在建筑和水体混合的拒马河区域识别出了多余的噪声斑点，并且将白洋淀外围的沙丘误提为水体，整体表达效果不如MNDWI精细、简洁。其余部位也反映出类似的情况，因此选择MNDWI-OTSU作为Sentinel-2数据水体提取方法。

3 基于Sentinel-1/2数据的洪水淹没范围快速提取模型构建

为了得到准确的水体提取范围，本文首先开展了不同方法提取效果对比研究。根据对比结果，构建了综合利用Sentinel-1/2数据的洪水淹没范围快速提取模型，并对模型精度进行了评价。

3.1 MNDWI-OTSU与VV-OTSU水体提取效果对比

同时期的地面水体真实范围难以实际测量［18］，本研究采用目视解译方式对研究区的水体进行解译，作为精度验证的基准。Sentinel-2 FRB（Full Resolution Browse）数据是由短波红外、近红外、红波段进行合成，并经过压缩和拉伸（gamma = 2.0）生成的优化高分辨率图像。该图像在不剪切极值的情况下能够突出植被和水体的差异［19］，目视解译效果较好，适用于研究区耕地、湿地、水体混合的复杂情景，如图2所示为研究区2023年7月16日的Sentinel-2 FRB影像。因此本研究选择Sentinel-2 FRB图像进行目视解译并与水体提取结果进行直观对比，从而评价本文方法的水体提取精度。Sentinel-2 FRB的时间信息与本研究使用的Sentinel-2影像一一对应，并且同样经过去云处理，具体见表1。

由于研究区内永久性水体占比较小，为了减少多余计算量，提高计算效率，在研究区内圈定7个水体占比在5%～30%的样地（A，B，C，D，E，F，G）来代替整个研究区进行精度对比研究。

选择2023年7月12日的Sentinel-1影像和云量情况较好的7月16日Sentinel-2影像分别使用VV-OTSU和MNDWI-OTSU方法进行水体范围提取（7月12日至7月16日研究区内无降雨发生），并与相同时段的FRB图像进行一一比对，由于篇幅原因，仅展示具有代表性的样地A的比对效果，如图3所示，红圈揭示出两种方法提取水体的差异。总体来讲，两种方法都取得了较好的效果，但是VV-OTSU方法相对于MNDWI-OTSU方法在细小水体较多的场景下仍有不足，主要表现在：

①部分零碎水体漏提、误提;② 由于SAR图像滤波处理，弱化了部分水体边缘，将水体之间的一些空隙连接了起来，因此在此场景下使用VV极化提取的水体面积比光学提取的面积更大。

随后，为了排除时间上的偶然性，选择时间相近且时间段内无降雨发生的3组Sentinel-1/2影像对，每个影像对各自影像的成像时间如表2所列。

分别使用VV-OTSU和MNDWI-OTSU对7个样地的水体范围进行提取，并统计水体面积占比，共得到21组数据，对其进行拟合分析得到图4中的结果。由图4，可以看出两组数据的水体提取结果具有很强的线性相关性（R2=0.9，RMSE=2.02%），证明在Sentinel-2数据缺失的多云天气，利用VV-OTSU对Sentinel-1进行提取得到的水体范围也足够稳定、可靠，可在该方法的基础上，拓展出后续模型。

3.2 基于Sentinel-1/2数据的洪水淹没范围快速提取模型

从3.1节的实验结果来看，MNDWI-OTSU比VV-OTSU方法更优，因此本文构建了一种综合利用Sentinel-1/2数据的洪水淹没范围快速提取模型来优化单独使用Sentinel-1数据的提取精度，其流程如图5所示。

该模型首先筛选出区域内用户所需时间段的Sentinel-1/2影像；然后对Sentinel-1影像进行预处理（Lee滤波等），并采用VV-OTSU方法提取水体得到区域水体粗分布图；接着根据时间间隔不超过3 d、间隔期间不发生降雨的时间配对原则筛选Sentinel-2影像与对应Sentinel-1影像配对，然后对配对的Sentinel-2影像进行预处理，对s2cloudless去云处理（云量分数阈值设为50）和MNDWI-OTSU提取水体，得到部分区域的水体精分布，再对水体粗分布图进行掩膜得到初步结果。为了防止Sentinel-2影像缺失时，Sentinel-1影像山体阴影与水体相互混淆，采用SRTM数据计算坡度栅格，将坡度大于5°的水体像元去除。模型搭载于GEE平台，实现数据快速获取、影像批量调用，极大提高了模型运行效率。

量化验证MNDWI-OTSU、VV-OTSU、本文模型的精度如下：在7个样地内随机采集500个样本点，包含250个水体点和250个非水体点，然后根据FRB图像目视解译结果，得到3种方法的水体提取精度，如表3所列。云量较少的Sentinel-2影像配合MNDWI-OTSU方法取得了更好的水体提取精度（总体精度、Kappa系数均大于90%），而VV-OTSU方法提取精度略低一点，与3.1节中目视对比的结果保持一致；而本文模型相比VV-OTSU精度有一定提升（总体精度、Kappa分别提升1.2%和2.4%），在缺乏可用Sentinel-2影像的情况下，达到了优化VV-OTSU方法的效果。

x2p+p5ZpNzrVGoqqa4GXJw==

4 洪水淹没范围时间序列提取应用

利用本文提出的基于Sentinel-1/2数据的洪水淹没范围快速提取模型，得到研究区2023年5～9月的地表水体分布时间序列，将每幅影像的水体面积进行统计，绘制水体面积变化如图6所示。

2023年7月24日至8月5日水体面积陡然上升，由107.9 km2扩大到308.0 km2，水体大幅扩张了近2倍。根据保定日报报道：7月29日至8月1日保定市出现大暴雨和特大暴雨，总降水时长超过了80 h。全市平均降水量350 mm，最大雨强82.2 mm/h，250 mm以上（特大暴雨）降水面积占比辖区面积的90%。暴雨引起了严重的洪水灾情，研究区境内白沟河、小清河等河道水位一路走高。水利部信息中心（http：∥xxzx.mwr.gov.cn/）数据显示：8月4日08：00，涿州白沟河东茨村水文站水位达27.49 m，超警戒水位（25.50 m）1.99 m。图7为涿州市洪水现场图。

以2023年7月24日影像为灾前影像，连续绘制后续时间点相比灾前影像新增水体范围图。其中，离洪水发生时间最近、受灾面积最广的8月5日的研究区水体范围如图8所示。

淹没面积较广的行政区为高碑店市、雄县以及安新县，洪水淹没的面积分别约为70.1，52.3，40.6 km2。将新增水体与Google地图进行对比可以发现，洪水淹没首先发生在涨水的河道周围，如涿州市东北部的码头镇，该镇境内琉璃河、拒马河蜿蜒流淌，受暴雨影响出现严重积水；其次是天然的封闭洼地（图8中高碑店市中部以西、定兴东部）；最后是行洪河道（图8中雄县中部、安新县西部），由于地势差等因素，洪水通过这些河道被引至附近的湖泊和河流中。

2023年8月5日至9月22日，各县相比7月24日新增水体面积统计如表4所列。

退洪方面，至2023年8月17日，定兴县、涿州市、蓉城县、高碑店市新增水体面积仅剩约2.2，2.8，0.8，5.6 km2，相比上期分别缩减87.7%，82.1%，77.8%，92.0%；至9月22日，上述四县水体面积基本恢复到灾前水平。雄县境内主要的积水区域为地势较低的行洪道，而安新县内则有白洋淀蓄滞洪水，因此这两县的水体消退相比上述四县更为缓慢。

5 结 论

本研究以Sentinel-1/2数据为数据源，使用NDWI-OTSU、MNDWI-OTSU、VV-OTSU方法分别提取两类数据的地表水体范围。实验结果表明：提取云量较少（<10%）的Sentinel-2影像水体时，MNDWI-OTSU比NDWI-OTSU阈值模型表达更为精细、简洁，精度达到95.6%；VV-OTSU阈值模型提取Sentinel-1影像水体相比MNDWI-OTSU阈值模型存在零碎水体误提、边缘模糊的问题，但总体精度仍达到93.8%。两种阈值模型在时间上提取的水体范围拟合度较高（R2=0.9），不随时间改变发生较大变化，具有足够的可信度。本文在VV-OTSU方法的基础上引入Sentinel-2数据和MNDWI-OTSU方法，构建了一种洪水淹没范围快速提取模型，并以保定市东部6县为例，对该区域2023年5～9月时序地表水体进行了提取与分析。结果表明，该模型精度达到95%，不受限于云雾且比VV-OTSU方法的精度更高。该模型可通过编码搭载于GEE平台，从而能快速、准确、低成本地实现对一定时间范围内某一区域的地表水体空间演变图连续输出。

参考文献：

［1］ 《中国水旱灾害防御公报》编写组.《中国水旱灾害防御公报2022》概要［J］.中国防汛抗旱，2023，33（10）：78-82.

［2］ 胡荣明，姚晓宙，李朋飞，等.Landast-8影像的水体指数构建及稳定性研究［J］.测绘科学，2022，47（5）：150-155.

［3］ 吴庆双，汪明秀，申茜，等.Sentinel-2遥感图像的细小水体提取［J］.遥感学报，2022，26（4）：781-794.

［4］ 赵文举，李聪聪，马宏，等.基于Sentinel-2超分辨率影像的干旱区水体提取方法［J］.农业机械学报，2023，54（10）：316-328.

［5］ 刘浩，周万蓬，张宇佳，等.基于Landsat影像的1999～2019年鄱阳湖面积动态监测［J］.东华理工大学学报（自然科学版），2023，46（1）：68-76.

［6］ 赵蓓蓓，黄海峰，邓永煌，等.基于Sentinel-1A的三峡库区范家坪滑坡InSAR监测分析［J］.人民长江，2022，53（10）：103-107.

［7］ 贾诗超，薛东剑，李成绕，等.基于Sentinel-1数据的水体信息提取方法研究［J］.人民长江，2019，50（2）：213-217.

［8］ AMITRANO D，DI MARTINO G，IODICE A，et al.Unsupervised rapid flood mapping using Sentinel-1 GRD SAR images［J］.IEEE Transactions on Geoscience and Remote Sensing，2018，56（6）：3290-3299.

［9］ BINH P D，CATHERINE P，FILIPE A.Surface water monitoring within Cambodia and the Vietnamese Mekong Delta over a year，with Sentinel-1 SAR observations［J］.Water，2017，9（6）：366.

［10］刘小燕，崔耀平，史志方，等.GEE平台下多源遥感影像对洪灾的监测［J］.遥感学报，2023，27（9）：2179-2190.

［11］李健锋.基于多源遥感影像和Google Earth Engine的斯里兰卡地表水体及洪水信息快速提取研究［D］.西安：长安大学，2020.

［12］郝金虎，韦玮，李胜男，等.基于GEE平台的京津冀长时序水体时空格局及其影响因素［J］.生态环境学报，2023，32（3）：556-566.

［13］LEE J S.Speckle filtering of synthetic aperture radar images ：a review［J］.Remote Sensing Reviews，1994，8（4）：313-340.

［14］孙亚勇，黄诗峰，李纪人，等.Sentinel-1A SAR数据在缅甸伊洛瓦底江下游区洪水监测中的应用［J］.遥感技术与应用，2017，32（2）：282-288.

［15］OTSU N.A threshold selection method from Gray-Level Histograms［J］.IEEE Transactions on Systems Man & Cybernetics，2007，9（1）：62-66.

［16］陈赛楠，蒋弥.Sentinel-1 SAR在洪水范围提取与极化分析中的应用研究［J］.地球信息科学学报，2021，23（6）：1063-1070.

［17］MCFEETERS S K.The use of the normalized difference water index（NDWI） in the delineation of open water features［J］.International Journal of Remote Sensing，1996，17（7）：1425-1432.

［18］XU H.Modification of normalised difference water index（NDWI） to enhance open water features in remotely sensed imagery［J］.International Journal of Remote Sensing，2006，27（14）：3025-3033.

［19］KHUONG H T，MASSIMO M，LI J.Surface water mapping and flood monitoring in the Mekong Delta using Sentinel-1 SAR time series and Otsu threshold［J］.Remote Sensing，2022，14（22）：5721.

（编辑：谢玲娴）

Research on extraction model of flood inundation range based on Sentinel-1/2 data

DENG Qirui1，2，ZHANG Ying1，2，LIU Jia1，2，QIAO Qinghua1，2，ZHAI Liang1，2

（1.Chinese Academy of Surveying and Mapping，Beijing 100036，China； 2.Key Laboratory of Surveying and Mapping Science and Geospatial Information Technology of Ministry of Natural Resources，Beijing 100036，China）

Abstract：

Remote sensing is an important means of monitoring the extent of flood inundation and understanding the evolution of flood disasters.However，optical images often have many deficiencies during floods，and all-weather SAR images have slightly lower accuracy in extracting water bodies.A flood inundation range extraction model based on Sentinel-2 optical images and Sentinel-1 radar image data was constructed to extract the flood inundation range quickly and accurately.An adaptive threshold segmentation algorithm，the Otsu algorithm，was used to extract the water body range of two types of data and the proposed model，and an application analysis was conducted using Baoding City，Hebei Province as an example.The results showed that Sentinel-2 images with less cloud cover had the best water extraction effect，with an overall accuracy （OA） of 95.6%.After introducing some available Sentinel-2 data，the OA of the constructed model reached 95.0%，OA and Kappa coefficient were increased by 1.2% and 2.4% respectively compared to using Sentinel-1 data alone.This model is installed on the Google Earth Engine platform and can achieve fast，accurate，and low-cost continuous output of the spatial range of surface water bodies.Clouds and mist do not limit it and have higher extraction accuracy than Sentinel-1 images alone.In the case of severe cloud coverage leading to a lack of Sentinel-2 data，this model can be used as an alternative method for extracting flood inundation areas.

Key words：

flood inundation range； Sentinel-1； Sentinel-2； adaptive threshold segmentation algorithm； Google Earth Engine； Baoding City