Potential and Limitations of Open Satellite Data for Flood Mapping

Abstract

Satellite remote sensing is a powerful tool to map flooded areas. In recent years, the availability of free satellite data significantly increased in terms of type and frequency, allowing the production of flood maps at low cost around the world. In this work, we propose a semi-automatic method for flood mapping, based only on free satellite images and open-source software. The proposed methods are suitable to be applied by the community involved in flood hazard management, not necessarily experts in remote sensing processing. As case studies, we selected three flood events that recently occurred in Spain and Italy. Multispectral satellite data acquired by MODIS, Proba-V, Landsat, and Sentinel-2 and synthetic aperture radar (SAR) data collected by Sentinel-1 were used to detect flooded areas using different methodologies (e.g., Modified Normalized Difference Water Index, SAR backscattering variation, and supervised classification). Then, we improved and manually refined the automatic mapping using free ancillary data such as the digital elevation model-based water depth model and available ground truth data. We calculated flood detection performance (flood ratio) for the different datasets by comparing with flood maps made by official river authorities. The results show that it is necessary to consider different factors when selecting the best satellite data. Among these factors, the time of the satellite pass with respect to the flood peak is the most important. With co-flood multispectral images, more than 90% of the flooded area was detected in the 2015 Ebro flood (Spain) case study. With post-flood multispectral data, the flood ratio showed values under 50% a few weeks after the 2016 flood in Po and Tanaro plains (Italy), but it remained useful to map the inundated pattern. The SAR could detect flooding only at the co-flood stage, and the flood ratio showed values below 5% only a few days after the 2016 Po River inundation. Another result of the research was the creation of geomorphology-based inundation maps that matched up to 95% with official flood maps.

1. Introduction

In recent years, the increased availability of free-of-charge satellite data has allowed the study of many natural and human-made processes at low cost and has boosted research in many fields [1,2,3,4]. For instance, the Sentinel satellite constellation of the Copernicus program of the European Union [5] provides synthetic aperture radar (SAR) and multispectral data with global coverage, high-frequency pass, and high spatial resolution. Other examples of free remote sensing programs are Landsat, which has provided data since 1972 [6], and the MODIS daily satellites giving multispectral images [7]. These data are often available with the first level of atmospheric or radiometric calibration, allowing their use by different types of users and not only experts in remote sensing processing. An example of a user-friendly data portal is the Worldview [8] service for the visualization of MODIS products or the G-Pod service of European Space Agency (ESA), which allows the on-line processing of ENVISAT and Sentinel-1 SAR data [9,10]. In addition, free GIS plugins allow the downloading and processing of free multispectral satellite images [11]. The availability of these resources is useful for the management of natural hazard effects.

Every year, flood events cause great economic losses and victims [12]. For this reason, precise flood mapping and modelling are essential for flood hazard assessment [13], damage estimation [14] and sustainable urban planning to properly manage flood risk [15]. In such a context, satellite remote sensing is currently a low-cost tool that can be profitably exploited for flood mapping [16,17,18].

Extraction of flooded areas can be performed by using multispectral satellite data and their derived indexes [19,20,21,22,23,24], SAR images [25,26,27,28,29,30,31], or a combination of these data [32]. Other types of satellite data can be useful to improve flood mapping. For instance, the digital elevation models (DEMs) derived from satellite data, e.g., Shuttle Radar Topography Mission (SRTM) and ASTER, were also used to estimate flood-prone areas or to improve a SAR/multispectral data-based map [33,34,35,36]. In addition, water storage data from the GRACE satellite [37] or soil moisture data from ASCAT [38] were used to derive flood indicators. Each remote sensing technique for flood mapping presents advantages and drawbacks [39] that must be evaluated on a case-by-case basis.

The frequent passes of satellites and the availability of rapid processing chains allowed the development of services providing automatic and quasi-real time flood mapping such as, for example, the Copernicus Emergency Management Service (EMS) performed by the European Union [40,41], the Global Flood Detection System [42,43] and the NASA Global Flood Mapping System [44]. However, these services provide rapid mapping products that can be affected by uncertainty and are not always validated [45]. Maps of flooded areas produced by official authorities and based on bespoke aerial photos and field surveys are more accurate, although they are time-consuming and require higher costs to be generated. Consider also that flood mapping can be affected by uncertainty related, for instance, to hydraulic models or the resolution of DEMs [46].

In this work, we present a multi-sensor, low-cost and user-friendly approach for flood inundation mapping. Specifically, we combined semi-automatic and manual approaches for flooded area detection and created rapid and accurate maps by integrating satellite images with DEMs and ancillary data. The aim of our research is the definition of a procedure that can also be used by non-remote sensing processing experts to map flooded areas using free satellite data. We developed and tested the presented methodology on flood events that occurred in Spain in 2015 and 2017 and Italy in 2016. We used all available free satellite data services: Sentinel-1, Sentinel-2, Landsat-8, MODIS and Proba-V. We also focused on the factors that limit, or help, the capability of flooded area detection (e.g., the times of satellite passes with respect to the flood and the spatial resolution). The produced flood maps were compared and validated against official maps and Copernicus EMS.

2. Description of the Study Areas

The study aims to test the use of satellite data for flood mapping. The first thing that we have considered is the main flood characteristics: The extent of the flooded areas, the length of the river segment examined and the temporal evolution of the flood. This information drives the choice of the available satellite data and data processing to be used.

The three selected cases studies, located in the Ebro River Valley and Arahal in Spain and the Po River basin in Italy (Figure 1), are representative of different situations concerning remote sensing images, ancillary data, geomorphology and land use.

The cases of the Ebro valley and Po basin floods are representative of rivers with a wide floodplain in which the flood peak requires from many hours to several days to transit the river length, and the inundated area can reach several kilometers in width. This type of flood is more prone to be mapped by satellites because of its spatial and temporal scale. In addition, the availability of ancillary or field data is higher for the main river. The case of Arahal is representative of flash-flood processes common in semi-arid and hilly environments that have rapid temporal evolution and limited flood area extent. In this case, mapping with satellite data is less straightforward.

2.1. February–March 2015 Flood in Ebro Valley (Spain)

The area is located in the middle of the Ebro River course in the Aragon region in (NE Spain; Figure 1A). Here, the alluvial plain is located in a 5–10 km-wide valley. The land use of the floodplain is primarily cultivated land, urban areas (Zaragoza, with approximately 600,000 inhabitants and other minor towns) and natural vegetation along the riverbanks. The floodplain is affected by periodic floods, particularly in the winter and spring seasons with several extreme events [47]. In the February–March 2015 period, several events of heavy rainfall occurred in the upper and middle basin of the Ebro River [48]. Consequently, the middle Ebro valley near the city of Zaragoza was affected by four different inundations. We focused our analysis on a section of approximately 60 km along the Ebro River from the villages of Alagon to El Burgo de Ebro. The most severe flood occurred approximately 2–5 March 2015, when the Ebro River reached a maximum discharge of approximately 2600 m³/s in Zaragoza gauge (E_Za in Figure 1A). Here, the flooding threshold is approximately 1500 m³/s, and the mean discharge of the Ebro River for March is 400 m³/s. The flood severity was estimated at an approximately 10-year return interval [48,49,50]. Large cultivated areas, some small settlements and linear infrastructures were inundated, resulting in severe damages [51]. The activation of Copernicus EMS EMSR-120 [52] allowed the automatic mapping of the flooded areas (delineation maps) using RADARSAT-2 data. In our study, we used Landsat-8, MODIS-Aqua/Terra, Proba-V and Sentinel-1 data for the mapping of flooded areas with the support of a DEM-based water-depth. The flood maps derived from remote sensing data were validated by using inundation maps and river discharge data available on the geoportal of the Ebro River Basin Authority [53].

2.2. November 2016 Flood in Western Po River Basin (Italy)

In late November 2016, a severe flood hit the western sector of the Po River basin corresponding to the Piemonte region (NW of Italy). This area, periodically affected by floods, particularly during the autumn and spring seasons [54,55], has already been investigated with remote sensing technologies [28,32,36]. In the period of 21–25 November 2016, heavy rainfall of up to 50% of the mean annual precipitation (e.g., 700 mm of cumulated rainfall at Piaggia station [56]) occurred in the western basin of the Po River and its tributaries, the Tanaro River in particular. The Po River at Caringano gauge (P_Ca in Figure 1C) reached a maximum discharge of approximately 2000 m³/s in the late evening of 25 November 2016. However, the mean discharge of the monitoring station in November is 70 m³/s. The Po River reached a discharge of 9900 m³/s (November average is 680 m³/s) after the confluence with the Tanaro. The Tanaro River grew to an estimated discharge of 3800 m³/s at Montecastello (T_Mo in Figure 1B) station (average discharge of 220 m³/S) [56]). Thus, inundation of cultivated areas, linear infrastructures and some urban areas, such as the town of Moncalieri [57], occurred with estimated damages of 50 M€. In this case, the Copernicus EMS (EMSR-192 [58]) was also activated, allowing the production of automatic flood maps based on using RADARSAT-2 and Pleiades images in the most critical areas. In our work, we used data from MODIS-Aqua, Proba-V, and Sentinel-1/2 satellites with the support of a DEM of the area. In particular, we selected the two areas that were most heavily affected by the flood: (1) The plain located around the Alessandria town (Figure 1B) was inundated by the Bormida River on 25 November 2016 and by the Tanaro on 25–26 November 2016. Here, the main land use is cultivated land with small farms; the town of Alessandria (95,000 inhabitants) was not affected by inundation; (2) The plain located at south of Turin city (Figure 1C) was flooded by the Po River and the other tributaries (i.e., Chisola and Oitana streams) on 25 November 2016. Here, the main land use is cultivated land, some villages and the urban area of Moncalieri (60,000 inhabitants). To validate our results, we used the flooded area maps made available by the Regional Environmental Protection Agency (ARPA is its Italian acronym) of the Piemonte region.

2.3. November 2017 Small Flash Flood Near Sevilla (Spain)

This study area is located near the town of Arahal in southern Spain (Figure 1D). The area is characterized by smooth hills and valleys of minor streams; cultivated land and olive trees are the main land use of the area. In addition, the commercial area of Arahal town (20,000 inhabitants) was partly affected by floods. Intense storms that occurred on 29 November 2017 caused the rapid increase of some small streams. In this case, no river-gauge discharge data are available but, from ancillary data found on the web, it was possible to estimate that the flood peak occurred in the morning of 29 November 2017. The flood destroyed a segment of the Malaga-Sevilla railway, causing the derailing of a train and several injured [59]. The Copernicus EMS was not activated for this small and rapid flood. We used Sentinel-1, Sentinel-2 and Landsat-8 data jointly with the water-depth DEM to map the flooded areas. In this case, no official maps were available to validate our results.

3. Materials and Methods

To map flooded areas, different remote sensing data and methodologies were used considering the following criteria: (I) type of data (SAR, multispectral); (II) cost (i.e., we used only free-of-charge data); (III) availability of data in relation to the flood phase (pre-flood, co-flood and post-flood); (IV) spatial resolution; and (V) availability of support data (DEM or other ancillary data).

We collected (a) pre-flood data acquired before the inundation and used as a reference for the change detection analysis; (b) co-flood data, collected around the time of maximum inundation and used for mapping flooded areas with limited post-processing; (c) post-flood data, taken after the flood event. Co-flood data is relative to the segment of river considered; downstream or upstream time and magnitude of flood peak could be different, particularly for large river basins. For instance, the flood wave of the Po or Ebro River takes several days to transit from the upper to the lower part of the catchment.

In the case of the Ebro valley and Po basin floods, we used river stage data to select the satellite data most suitable for our analysis. Figure 2 shows the evolution of hydrometric levels during and after the flood events for Po River at Carignano gauge [60] (P_Ca in Figure 1C), Tanaro River at Alessandria gauge [61] (T_AL in Figure 1B) and Ebro River at Zaragoza gauge [53] (E_ZA in Figure 1A) versus the satellite passes.

In the case of the Po and Tanaro Rivers, the number of useful images acquired during the co-flood phase (approximately 12 h) is limited, and most of the analysis was performed by using post-flood data. In the case of the Ebro River valley, the long co-flood phase of March 2015 (2–3 days) allowed us to collect more co-flood data.

In the Supplementary Materials, a list of the acquisition dates used in this work is provided.

Table 1. Characteristics of the data used in our study.

Satellite	Band Used	Spatial Resolution (m)	Revisit Time (Day)	Data Availability
MODIS-Terra/Aqua	R; G; B; NIR; SWIR(2)	250 ÷ 500	1	[69]
Proba-V	R; B; NIR; SWIR	100 ÷ 300	5 ÷ 1	[71]
Landsat-8	R; G; B; NIR; SWIR(2)	30	16	[72]
Sentinel-2 (a/b)	R; G; B; NIR; RE(2); SWIR(2)	10 ÷ 20	5 ÷ 10	[62]
Sentinel-1 (a/b)	C	5 × 20 (SLC)	6 ÷ 12	[62]

and

Table 2. Bans and wavelength used for this study.

Band Name	Wavelength Interval
Blue (B)	~450–500 nm
Green (G)	~500–620 nm
Red (R)	~620–700 nm
Red Edge (RE)	~700–750 nm
Near Infrared (NIR)	~750–1400 nm
Shortwave Infrared (SWIR)	~1400–3000 nm
SAR-Band (C)	~5 cm

summarizes the revisit time of the satellite, the bandwidth, the spatial resolution and the data provider of the satellite images used for this study. Two types of satellite are exploited: Active SAR sensors operating in the microwave domain and passive multispectral sensors ranging from visible to thermal infrared. Among the multispectral sensors, the spatial resolution ranges from low (500 m/pixel) to medium-high (10 m/pixel).

3.1. SAR Data

The main advantage of SAR sensors is given by their capability to acquire images both at night and under all weather conditions, thus filling the gap often resulting from cloudiness-corrupted optical datasets.

Currently, the Sentinel-1 constellation of the Copernicus program is the only provider of SAR data, freely available through the Sentinel Scientific Data Hub [62]. The constellation is composed of two satellites, Sentinel-1A and Sentinel-1B, operating since 2014 and 2016, respectively. The twin satellites acquire C-band (central frequency of 5.404 GHz) SAR data all over the world, with a revisit time of 6 days over Europe and some hotspot areas, and 12 days in the rest of the world. Such a short repeat cycle increases the chance to collect free co-flood data that allow mapping inundated areas accurately.

For our study, we used the Interferometric Wide (IW) swath acquisition mode, which captures three sub-swaths by employing the Terrain Observation with Progressive Scans SAR ((TOPSAR) [63]). IW is the main acquisition mode for the systematic monitoring of surface deformation and land changes, providing data with a 250 km swath and a spatial resolution of 5 m × 20 m. In particular, we used the Single Look Complex (SLC) format, i.e., complex-valued images in the slant range by the azimuth imaging plane [64,65]. We also used the ESA SNAP open source software for geocoding SAR data [66] and converting the pixel values from digital numbers (DNs) into measurements of the sigma naught values (in dB). To this end, the 3-arcsec SRTM DEMs have been used on the investigated areas.

3.2. Multispectral Data

Multispectral satellites acquire images using several different wavelength bands. This technique enables the derivation of different spectral indexes, e.g., Normalized Difference Vegetation Index (NDVI) and Normalized Difference Water Index (NDWI), allowing the identification of flooded areas (including after the water has already withdrawn) by analyzing the flood’s effects on the soil [67]. Considering the spatial resolution, it is possible to distinguish two classes of multispectral satellites: The medium-low resolution (100–500 m) satellites, such as MODIS, Proba-V or Sentinel-3, and the medium-high resolution (10–30 m) satellites, such as Sentinel-2 or the Landsat series.

3.2.1. Medium-Low Resolution Multispectral Data

For our analysis, we used MODIS (Moderate Resolution Imaging Spectroradiometer) and PROBA-V (PRoject for On-Board Autonomy-Vegetation) data. MODIS is a system of two sun-synchronous, near-polar orbiting satellites called Aqua and Terra that daily acquire images all over the world [7] that can be visualized on the NASA Worldview portal [8]. Terra collects images in the late morning and Aqua in the early afternoon; they also have a night-time pass when they acquire in thermal bands. We used the MYD09 (Aqua) and MOD09 (Terra) Atmospherically Corrected Surface Reflectance 5-Min L2 Swath 500 m [68], downloaded from the Level-1 and Atmosphere Archive and Distribution System (LAADS [69]). These products include seven bands (RGB; NIR and SWIR) with a spatial resolution of 500 m. We also downloaded 250 m spatial resolution Red and NIR band MOD09/MYD09 through the SCP plug-in of QGIS [11] to calculate the NDVI.

MODIS images were used for the Po basin and Ebro areas, for which co-flood data not affected by extensive cloud cover and post-flood data were available. In the Ebro valley, MODIS images were available for all three of the flood events that occurred in February–March 2015.

For the Ebro valley and Po basin floods, we also used data acquired by the ESA PROBA-V, mission launched in 2013 for vegetation monitoring [70]. This satellite has four bands (Blue, Red, NIR, and SWIR) with spatial resolution ranging from 100 m (revisit time of 5–6 days) to 300 m (daily revisit time). We used the Level 2A—radiometrically corrected segments, projected in plate carrée, Top-of-Atmosphere with 100 m of spatial resolution downloaded from the ESA Proba-V portal [71]. The small area of Arahal is not suitable for analysis with the coarse spatial resolution of MODIS and Proba-V.

3.2.2. Medium-High Resolution Multispectral Data

For our analysis, we used Landsat-8 and Sentinel-2 data.

Landsat-8, launched in 2012, acquires images with 11 spectral bands; it has a spatial resolution ranging from 15 to 100 m and a revisit time of 16 days. In particular, we used six bands (visible, NIR and SWIR) with 30 m spatial resolution. Level-1 data were downloaded from the Earthexplorer portal of USGS [72] and used for the Arahal and Ebro case studies. For the Po and Tanaro flooded areas, no cloud-free images are available close to the flood event, and we did not use Landsat-7 images because they are affected by a no-data strip.

For the Arahal and Po basin areas, we used the Sentinel-2 constellation composed of two satellites, S-2A and S-2B, launched in June 2015 and on March 2017, respectively. These satellites are characterized by a revisit time of 5 days in Europe and 5 to 10 days in the rest of the world, and they acquire data in 12 bands with spatial resolution ranging from 10 to 60 m. For our study, we focused on the RGB, NIR and SWIR bands (10 to 20 m spatial resolution; data were downloaded from the Sentinel scientific data hub [62]).

3.3. Data Processing, Analysis and Validation

Figure 3 shows the methods and the procedures used for the flood mapping.

The analysis is divided into three main steps:

• The first step is the detection of the flooded area, which can be performed using a manual or a semi-automatic mapping approach:
  • The manual mapping consists of the direct visual interpretation of the images (SAR amplitude or color combinations of multispectral bands). In this case, the flooded areas were drawn manually directly on the georeferenced satellite images in QGIS software.
  • With the semi-automatic approach, we initially produced an automatic flooded area map in raster format. The map is extracted from SAR or multispectral satellite data using different methodologies such as band index, supervised classification or backscattering difference. In this step, we used an empirical threshold to detect flooded areas; for this reason, it is not a fully automatic approach.
• A possible improvement of manual and automatic detection could be made using a cloud mask and permanent water body (from ancillary data or pre-flood images). The second step is map improvement and refinement. In this step, we also consider additional information such as (a) water depth model derived from DEM, (b) hillshade and aerial photos to detect the geomorphological features, and (c) ancillary data such as georeferenced photos or documents found on the web to have ground information about the flooded area extent. These data allow the creation of an improved final version of flooded area maps, manually drawn, both for the semi-automatic and manual approaches.
• The third step is the flood map validation. This step is performed only when official flood maps or field survey maps are available. We used these maps to evaluate the quality of the flooded area maps and in particular the performance of semi-automatic mapping (flood ratio and not flood ratio).

3.3.1. Manual Mapping

The manual mapping of flooded areas is made by drawing on QGIS software the polygon of the flooded area directly on the georeferenced satellite image. The image is obtained with the composition of multispectral bands in true or false color. In the case of co-flood imagery, the identification is simple; however, for post-flood data, the creation of an enhanced color image is necessary to identify the traces of flooding. For the manual mapping of the flooded area, we also try to use the SAR amplitude images of the Sentinel-1, geocoded with ESA SNAP software.

3.3.2. Semi-Automatic Mapping Based on SAR Data Processing

• SAR Amplitude classification (SAR_AC). The calibrated SAR (amplitude) image was classified to identify the water-covered area according to the criterion that soil covered by a quiet water table shows low amplitude and can be easily extracted from the image. We applied the classification to Sentinel-1 SLC for the Zaragoza, Po basin and Arahal areas. We also implemented a simple filter, available in the processing raster tool of QGIS software, to smooth local noise effects. Using empirical thresholds based on a visual approach, we classified the SAR data into three classes: (a) Low values of sigma naught (σ°) (< −25/−20 dB), i.e., water-covered areas (lakes or flooded areas), flat surface or shadow; (b) middle values (−20/−10 dB), i.e., cultivated land or natural vegetation; and (c) high values (> −12 dB), i.e., dense urban or forested areas or area with bare rocks. We used optical images to validate the accuracy of classification qualitatively. During the refinement step, the operator can separate the flooded area from another type of low-dB area with the help of ancillary data (e.g., permanent water body mask or optical images).
• Backscattering variation (Δσ^o). Two SAR images (one post-/co-flood and one pre-flood) were used for change detection analyses. SAR data, acquired in the SLC format, were radiometrically calibrated [66] to obtain the relevant σ^o maps, i.e., the surface backscattering maps. Calibrated SAR data were averaged through a multi-look operation, with one look in the azimuth direction and four looks in the range one, accurately co-registered [73] and, finally, geocoded [74]. Surface changes due to flooding were detected by calculating the log ratio between the post- and the pre-flood images. Thus, the Δσ^o map showing the temporal variation of surface backscattering from the pre- to the post-flooding phase was produced. For the classification of such a map, we applied an empirical threshold of the Δσ^o value that allows identifying as many flooded areas as possible and minimizes false detection errors.

3.3.3. Semi-Automatic Mapping Based on Multispectral Image Processing

The first presented method is based on the analysis of a single image. The others are based on a comparison of two images (post or co-flood against pre-flood images). The thresholds used to discriminate flooded from not-flooded pixels are empirically based on a visual approach; we set the threshold that best defines a geomorphology- based flood pattern.

• Supervised classification (SC). We applied this technique to classify MODIS and Landsat-8 co-flood data. We initially detected the most representative land-cover types (i.e., water-covered area, vegetation, cloud, snow, and urban area/bare soil) in training areas selected over a composite band image on QGIS to create a spectral signature with the available bands. Then, we performed a supervised classification of the images using the SAGA-GIS. Specifically, we tested different classification methods and selected the maximum likelihood and spectral angle methods as the most appropriate for our study. Finally, using a raster query in QGIS, we extracted the “area covered by water or wetland” category. In this case, the use of a permanent water body mask also allowed us to separate the flooded area from the permanent rivers and lakes.
• Normalized Difference Vegetation Index variation (NDVI_var). We calculated the NDVI variation between the pre- and post-flood conditions (Equation (1)). The aim is to identify areas characterized by a decrease in NDVI values due to vegetation activity decreasing (area covered by sediments or damaged vegetation) or the presence of water covering the area [75,76]. We computed NDVI using the NIR and the red band of Sentinel-2 (10 m of spatial resolution), Landsat 8 (30 m of S.R.), Proba-V (100 m of S.R.) and MODIS data (250 m of S.R.):

  $$\mathrm{NDVIvar}=\mathrm{NDVIpost}-\mathrm{NDVIpre}\mathrm{where}\begin{array}{c}\mathrm{NNDVI}=\frac{\left(\mathrm{NIR}-\mathrm{Red}\right)}{\left(\mathrm{NIR}+\mathrm{Red}\right)}\end{array}$$

  (1)
  Modified Normalized Difference Water Index variation (MNDWI_var) and Normalized Difference Moisture Index variation (NDMI_var). These indexes allow the detection of water bodies or wetlands and were widely used to map flooded areas [77]. In the literature, different combinations for calculating MNDWI are presented and discussed [78,79,80].
• In our study, we used the Red Edge band—Short Wavelength Infrared bands (Equation (2a)) for Sentinel-2 data (20 m of S.R.), and the Red and SWIR bands for Landsat-8 (30 m of S.R.), MODIS (500 m of S. R.) and Proba-V data. In a pre-flood situation, these bands better supported noticing changes in soil moisture/water covered areas. The NDMI [81] is calculated by exploiting the NIR and SWIR bands and was used for the Arahal area (Equation (3)).

  $$\mathrm{MNDWIvar}=\mathrm{MNDWIpost}-\mathrm{MNDWIpre}$$

  (2)

  where

  $$\mathrm{MNDWI}=\frac{\left(\mathrm{RedEdge}1-\mathrm{SWIR}\right)}{\left(\mathrm{RedEdge}1+\mathrm{SWIR}\right)}$$

  (2a)

  or

  $$\mathrm{MNDWI}=\frac{\left(\mathrm{Red}-\mathrm{SWIR}\right)}{\left(\mathrm{Red}+\mathrm{SWIR}\right)}$$

  (2b)

  $$\mathrm{NDMIvar}=\mathrm{NDMIpost}-\mathrm{NDMIpre}\mathrm{where}\mathrm{MNDWI}=\frac{\left(\mathrm{NIR}-\mathrm{SWIR}\right)}{\left(\mathrm{NIR}+\mathrm{SWIR}\right)}$$

  (3)
• Variation of averaged visible bands (VIS_var). This index allows the detection of the reflectance variation induced by the presence of silt deposits on crop fields. The index is based on the averaging of the RGB bands (Equation (4)). We derived this index from the Sentinel-2 data (S.R. 10 m) acquired over the Tanaro and Po areas, allowing noting of the sediment deposits and mapping the flooded areas indirectly.

  $$\mathrm{VISvar}=\frac{\left(\mathrm{VISpost}-\mathrm{VISpre}\right)}{\left(\mathrm{VISpost}+\mathrm{VISpre}\right)}\mathrm{where}\mathrm{VIS}=\frac{\left(\mathrm{R}+\mathrm{G}+\mathrm{B}\right)}{3}$$

  (4)
• SAR-Optical combination: The extraction of the flooded area for the Arahal case study was also obtained by integrating the processed SAR (Sentinel-1) and multispectral (Sentinel-2 and Landsat-8) data. Specifically, we created a Boolean raster map based on the relationship AND/OR between SAR- and optical-based flood maps.

3.4. Improvement and Refinement of Flood Maps

The semi-automatic or manually drawn flood maps based on remote sensing data were improved and refined by using DEM modelling and ancillary data. For such support data, we also followed the criterion of using free-of-charge data. As a final step, manually drawn flood maps were produced using QGIS software.

3.4.1. DEM Support: Water Depth Model and Shaded Relief

We used 5 m and 25 m resolution DEMs produced using LiDAR data and photogrammetry techniques, freely available for downloading from the geoportal service of Regione Piemonte [82], and the National Geographic Institute of Spain (IGN is the Spanish acronym) [83]. We used the methodology described in [57,84] With this method, water height information based on ground truth data flooded area polygons based on the intersection of remote sensing data overlapped with DEM are used to compute a raster of water level (WL). Then, the difference between the WL and the DEM allows obtaining of the simulated water depth model (WD) map. The water depth is also an essential parameter for damage assessment. We also used the shade relief model derived from DEMs and aerial photos from the WMS services of Regione Piemonte or IGN to detect geomorphological features such as terraces or embankments that could constrain inundation areas.

3.4.2. Ancillary Ground Truth Data

We used ancillary data such as field observations, civil protection reports, geolocated photos or news found on reliable websites. This type of data is currently an essential source of information for flood map assessment [85,86]. For the Po basin, we also accessed volunteer geographic information that provided us photos that were georeferenced by using Google Streetview. These data provided spot information that allowed confirming of the reliability of flood maps. As mentioned previously, other ancillary data (such as permanent water body masks) were used to improve flood mapping.

3.5. Flood Map Validation and Quality Statistics

The final step of the work is the validation and quality assessment of all the produced flood area maps (semi-automatic and manual). To do this task, we used the flood maps made by official authorities and published online as a benchmark. These maps are based on aerial photos taken a few days after the floods, field surveys, and remote sensing data. For the 2016 Po River basin flood, we used a flood map available from ARPA (Regional Environment Protection Agency) Piemonte geoportal [87], whereas for the 2015 Ebro flood, we used the maps of the Ebro River Basin Authority [53]. We could not apply this step to the Arahal study case because no official maps are available.

The validation procedure is described in the following:

(a)

Both the produced flooded area maps and the official reference maps were converted into a Boolean raster format (1—flooded or 0—not flooded)

(b)

The two rasters were crossed, and a raster map with four possible values was generated: True positive (TP), true negative (TN), false negative (FN), false positive (FP). TP corresponds to areas correctly classified as flooded, TN are the areas correctly classified as not flooded, FN indicates the undetected flooded areas, and FP are the not-flooded areas erroneously classified as inundated. To make a more-uniform comparison between different types of semi-automatic flood detection approaches, it is possible to use a permanent water mask to exclude river and lake areas from FR statistics.

(c)

We evaluated the efficacy of used datasets and methods by computing the flood-mapping ratio (Equation (5)) and the not-flood ratio (Equation (6)) represented as a percentage format.

FR (%) = TP/(TP + FN)

(5)

NFR (%) = TN/(TN + FP)

(6)

4. Results

In this section, we present the results achieved in the selected study areas by applying the abovementioned methodologies. In the following figures, we use the polygon of definitive flooded areas as a benchmark.

4.1. Manual Mapping of the Flooded Area

We tested the possibility of mapping the flooded areas of all our case studies manually. In some cases, this mapping was possible, mainly if a co-flood image was available. In general, areas covered by water are clearly visible on multispectral or SAR images, and the manual drawing can be done. In particular, false color images (Figure 4) made by using SWIR-NIR-BLUE bands of multispectral data allow better definition of the inundated areas. The main advantage of visual mapping is the possibility of minimizing errors (e.g., false positives and false negatives) of the semi-automatic classification methods, using geomorphological criteria during drawing. Conversely, the manual mapping is more time-consuming and requires the operator to possess experience in geomorphology and the dynamics of floods.

Figure 4 shows examples in our study areas of visual interpretation using remote sensing data.

• In the Ebro valley near Zaragoza (Figure 4A), we used a Landsat-8 co-flood image acquired on 2 March 2018; in this case, the false-color image allowed mapping the entire flooded area.
• In the Po River basin, around Alessandria town (Figure 4C), we used the Proba-V image acquired in the morning of 26 November 2016. In this case, it was possible to map most of the area flooded by Tanaro River, but the area flooded by Bormida River was not clearly detectable.
• In the Po River basin, at the South of Turin city (Figure 4B), the Sentinel-2 post-flood image shows only weak traces, and the semi-automatic classifications better identified the inundated area.
• In the Arahal area (Figure 4D), the Sentinel-2 post-flood image shows flood traces that are very difficult to detect, and the manual mapping fails. Here, we needed to use the automatic extraction of the flooded area for more reliable results.

In most of the considered cases, the polygon of flooded areas drowned manually needs refinement via DEM and other ancillary data.

4.2. Semi-Automatic Mapping of the Flooded Area

The applied semi-automatic approach consists of initially elaborating the co-flood images or the post- and pre-flood images. For each type of elaboration, we defined an empirical threshold that allowed detecting of the flooded area. All of this processing was done in QGIS software, and the automatic flood maps are in raster format.

4.2.1. SAR Data Analysis

In this section, we show the results of the SAR data-based analysis performed by using SAR amplitude classifications (SAR_AC) (Figure 5) and Δσ^o backscattering variation in Figure 6:

• SAR_AC. In the case of the Ebro River flood, the classification of the Sentinel-1 image acquired on 5 March 2015, i.e., two-three days after the flood peak (when most of the areas remain inundated according to MODIS data), allowed good detection of the flooded areas (Figure 5A). It is also possible to observe false positives, i.e., areas wrongly classified as flooded, most likely related to topographic artefacts or flat surfaces (e.g., airport ramp). Figure 5B shows the SAR_AC based on the Sentinel-1 image acquired on 28 November 2016 over the areas inundated by the Tanaro and Bormida Rivers near Alessandria. The flooded areas were detectable more easily for the Tanaro River, where the soil remained submerged two days after the flood peak, than for the Bormida River. The analysis performed over the area flooded by the Po River did not show significant results. In the Arahal area, the flooded area was partially detected with images acquired in the early morning of 29 November 2017, a few hours before the flood peak (Figure 5C).
• Δσ^o backscattering variation-based analysis. Compared with the SAR_AR-based analysis, such an approach allows achieving more-accurate results by exploiting calibrated Sentinel-1 data and image log ratios between the post- and pre-flood images (i.e., Δσ^o = log σ^o post-log σ^o pre). In addition, the comparison with pre-flood data allowed discriminating flooded areas from permanent water bodies. For most of the cases, we used VH polarization and Δσ^o > 0.6 as a threshold. In the case of Zaragoza (Figure 6A), we used the post-flood image of 5 March 2015 (acquired three days after the discharge peak), and the pre-flood image of 9 February 2015. The flooded area was partially detected with the semi-automatic methodology. However, this result allowed mapping of the entire inundated area with a visual interpretation. The results are also in agreement with MODIS data acquired on the same date. The comparison with a water-depth model showed that the areas remaining classified as flooded on the Sentinel-1 image approximately correspond to regions in which the water depth during the flood peak was higher than 1.5 m. In the case of Alessandria and Turin, we used the post-flood image of 28 November 2016 and the pre-flood image of 22 November 2016. Results show residual inundated areas. In particular, as already mentioned, the flooding pattern remains visible in the case of Tanaro River but not for the southern region of Turin, in which the remaining flooded areas are too small to be detected with Sentinel-1 data. In the case of Arahal (Figure 6D), we used VV polarization, which showed better performance concerning VH. We compared pre-flood (23 November 2017) with near co-flood images of 29 November 2017. In this case, we found that the threshold of 0.3 was the best value to detect flooded areas. The Sentinel-1 image was acquired in the early morning of 29 November 2017, some hours before the flood peak; most likely, the inundated areas are underestimated.

In all of our study areas, the capability to map flooded areas by SAR data is also limited by the presence of dense arboreal vegetation (small areas in the Ebro valley) or dense urban settlements (area of Moncalieri [57]). This limitation is already known in the literature [26]; more-complex SAR processing was developed to detect flooded areas in forested or urban areas [88,89]; moreover commercial L-Band satellites (e.g., ALOS-2) have better performance on forested areas [90].

4.2.2. Multispectral Data Analysis

In this section, we present the results of the multispectral analyses relying on single data methods (e.g., supervised classification) and on methods based on two-image comparisons.

(a) Supervised classification.Figure 7 shows the supervised classification based on MODIS and Landsat-8 multispectral images. We have defined four primary land-use typologies: Vegetation, urban area and bare soils, water body/wetland (that mostly correspond to the flooded areas) and cloud cover. The supervised classification showed good results for both the 2 March 2015 MODIS-Terra image (Figure 7A) and the 3 March 2015 Landsat-8 (Figure 7B) collected over the Ebro River valley, for which the scenes were almost cloud-free and the river was at the maximum discharge. The Landsat-8-based classification showed errors in the areas where the vegetation emerged from the water. The details were provided by the finest resolution. In the case of the Po basin flood, the supervised classification of the 26 November 2016 MODIS-aqua image (12–24 h after the flood peak) showed an underestimation of the flooded areas in the sectors inundated by the Bormida (Figure 7C) and by the Oitana and Chisola streams (Figure 7D). Classification errors due to partial cloud cover and shadow classified as water can be seen in the area around Alessandria (Figure 7C), but they do not significantly affect flood detection.

(b) Post-Pre flood image comparison. The comparison of co-flood or post-flood images with pre-flood images allowed mapping the flooded areas using different band-ratio variations.

(I) The MNDWI_var shows different performances depending upon the time of the satellite pass and local area setting. The threshold used to detect flooded areas ranged from 0.1 to 0.2. There was a decrease in the efficacy of this index with time elapsed from the flood event.

In the Zaragoza case study (Figure 8A), the flooded area was almost entirely detected by the co-flood Landsat-8 image acquired on 3 March 2015. Some false positive variation at SW is related to cloud cover. The MNDWI_var computed with the 2 March 2015 MODIS image provided the same results.

In the case of the November 2016 Po basin flood, only post-flood data are available at high resolution. At south of Turin, the flooded areas were satisfactorily detected with the 26 November 100 m spatial resolution Proba-V image (Figure 8B), whereas by using Sentinel-2 imagery acquired (Figure 8C) five days after the water was withdrawn, it remains possible to detect a pattern of positive MNDWI_var inside the inundated areas, but with more false positive areas. Over the Alessandria area, the MNDWI based on 100 m spatial resolution Proba-V (Figure 8E) and 500 m spatial resolution Proba-V (Figure 8E) data acquired approximately 8–12 h after the flood peak detected most of the areas inundated by Tanaro, despite cloud cover and shadows. Sentinel-2 images (Figure 8F) acquired approximately 10 days after the flood peak only show residual patterns. Many pixels located out of flooded areas are classified as inundated, most likely due to a general increase of soil moisture. We also applied a cloud mask to the images of 26 November of MODIS and Proba-V; fortuitously, the cloud cover and shadow did not affect the detection of the flooded area.

In the case of the Arahal flood (Figure 8F), the MNDWI_var based on the Sentinel-2 image acquired 10 days after the flood event allowed the detection of only some residual parts of the original flood. The NDMI_var based on Landsat images showed a similar result.

(II) In the case of the 2016 Po basin flood (Figure 9), we also calculated the NDVI_var and the VIS_var to extract the flooded area automatically. The areas affected by floods show a negative variation of NDVI most likely related to a decrease of vegetation activity due to the deposition of silt layers or damages (e.g., broken stems of the crop). This effect is particularly evident in the wheat fields outside the flooded areas, showing an increase of NDVI. Moreover, for the NDVI_var method, the capability to detect flooded areas decreases with the time elapsed from the flood peak. This decreasing trend is independent of the spatial resolution of the satellite. The NDVI_var made by using the MODIS-Aqua images collected on 30 November 2016 and 12 November 2016 at a spatial resolution of 250 m/pixel shows a clear pattern of the flooded area for both the Po (Figure 9E) and Tanaro areas (Figure 9D).

Moreover, the high-resolution NDVI_var based on Sentinel-2 data continues to shows a pattern of negative variation both in the South Turin area with an image acquired on 1 December (Figure 9B) and in the Alessandria area with the 8 December images (Figure 9A,C). The VIS_var shows a clear pattern of positive variation only for the area flooded by the Tanaro River nearby Alessandria (Figure 9B) due to the widespread deposition of silt sediments covering the area two weeks after the inundation. The areas flooded by the Po and other rivers showed less evidence of sediment; accordingly, the VIS_var analysis is less performant.

4.3. Flood, Mapping Refinement Using Water Depth DEM and Ancillary Data

As already mentioned, the second step for the flood map creation is the refinement and improvement based on supporting data (i.e., DEM and ancillary data), the result of which is particularly useful when we have flood maps based only on post-flood images. The final result was the production of flooded area manually drawn polygons. Figure 10 shows an example of the integration of remote sensing data, DEM modelling and ancillary data for the mapping of the areas flooded by the Tanaro River near Alessandria. The map in Figure 10A represents the inundation area derived from the semi-automatic classification of remote sensing images based, in this case, on VIS_var of Sentinel-2 data.

Figure 10B shows the inundation maps based on a water depth (WD) model in which the pixel is classified as inundated when WD > 0 m,

Figure 10C represents the synthesis of the two presented maps, in which a pixel is classified as inundated if flooded in both satellite data and the WD model. The WD model allowed reducing pixels incorrectly classified as flooded (false positive). Furthermore, the geolocated photos (available on websites or from volunteers) allowed us to have some ground truth points.

Figure 10D represents the final version of the flooded area that was manually drawn considering the remote sensing data, the WD model, and the ancillary data.

We also used the same methods to create the final maps of the area flooded by the Po and Ebro Rivers. The refinement step was handy for the case of the Po River, whereas for the Ebro, for which the high-resolution co-flood map is available, the support data do not provide significant improvement.

For the Arahal case study, in which the flooded area was more difficult to detect, we made a more in-depth effort to combine satellite and support data. Figure 11 shows the steps used to create the refined flood map.

We initially created a semi-automatic map of the flooded area (Figure 11A) by crossing all of the available remote sensing data (i.e., Landsat-8, Sentinel-2 and Sentinel-1). Then, we created a WD model map for the two small streams affected by the flood, around the intersection with the Seville-Malaga railway (Figure 11B). The combination of such a model with remote sensing data allowed us to produce the map of the estimated flooded area (Figure 11C). To assess the reliability of the flood map, we exploited geolocated photos and news about the flood. Available images allowed an evaluation of the reliability of the inundation map locally, e.g., it was possible to locate the effect of embankment erosion on a railway line (photo 3).

4.4. Flood Map Validation and Statistics

The last step of the presented methodology consists of validating the automatic flood maps. We applied the validation to the flood maps of Ebro valley and the Po basin using the official flood maps of the Ebro River authority (CHE) and ARPA Piemonte.

The validation process, as described earlier, is based on the definition of a raster validation map. In these raster maps, each pixel can assume four values: (0) True negative (TN) case represented in brown; (1) false positive (FP) case represented in cyan; (2) false negative (FN) case represented in yellow; and (3) true positive (TP) case represented in blue. To make a more homogeneous comparison, the area covered by a permanent water body was not calculated in FR or NFR.

In Figure 12, examples of validation maps for semi-automatic flooded area maps are reported. For the area of Zaragoza (Figure 12A), the MNDWI_var derived from co-flood Landsat-8 images allowed correctly classifying both flooded and not-flooded areas (>90%); the Δσ^o Sentinel-1-based analysis (Figure 13B) also shows satisfactory performance (FR 50%; NFR > 90%).

For the area of Alessandria, the NDVI_var computed from the MODIS-Terra data shows better results (FR 70%, Figure 12C) compared with the Δσ^o-based analysis (FR of 30%, Figure 12D).

In the area at the south of Turin, the MNDWI_var (Figure 13E) made with a post-flood image (1 December 2016) from Sentinel-2 data detected approximately 45% of the flooded area.

We also applied the validation to the final refined maps of flooded areas. The results show that our maps mostly coincide (FR > 90%) with official maps of ARPA Piemonte and the Ebro River authority in the case of Zaragoza and Alessandria, but more differences can be found in the maps flooded by the Po (FR 70%).

Figure 13 shows a correlation between the river levels for the Ebro (Figure 13A), Tanaro (Figure 13B) and Po (Figure 13C) with the flood ratio (FR) for semi-automatic flood mapping based on different methodologies and satellites. The results give an idea of the various parameters that influence flooded area detection. In general, the quality of semi-automatic detection of flooded areas decreases with the time after the peak of the flood for all of the cases considered. The SAR data show a much more rapid decrease in the capacity to detect flooded areas.

For each area, the primary results are as follows:

• The area flooded by the Ebro River near Zaragoza in 2015 (Figure 13A). In this case, co-flood data at low resolution (MODIS-Terra/Aqua and Proba-V) and medium-high resolution (Landsat-8) allowed mapping the entire flooded area (FR > 95%) with few false positive values. In this case, spatial resolution has little influence on the accuracy of flooded area detection. The Δσ^o map made with Sentinel-1 data acquired two days after the maximum flood shows that approximately 50% of the area is detectable with an NFR (≈98%). With the map based on Sentinel-1 data, it was nonetheless possible to detect a clear pattern that identifies entire flooded areas.
• The area flooded by Tanaro near Alessandria in 2016 (Figure 13B). Here, there are high values of FR for near co-flood data of MNDWI based on MODIS-Aqua (FR ≈ 90%) and Proba-V (FR ≈ 82%). For a post-flood image, the relationship “time of flood peak vs. time of satellite pass” is more important than spatial resolution; 250 m NDVI based on MODIS data of 30 November 2016 show a better result (FR ≈ 80%) than do 10 m NDVI based on Sentinel-2 of 8 December 2016 (FR ≈ 60%). SAR data based on 28 November imagery from Sentinel-1 show low values of FR (≈20%), making detecting the entire flooded area more complicated. Another critical factor influencing flood detection with post-flood data is the local conditions. In particular, the presence of a thin layer of silt deposits on crop fields helped the identification of the inundated area. Two weeks after the flood, the VIS_var based on Sentinel-2 data showed an FR of 70% and NFR ≈ 84%
• The area flooded by the Po at the south of Turin, compared with Carignano gauge (Figure 13C). Here, the flood detection presented results similar to those of the Alessandria area. In particular, Sentinel-2 data acquired on 1 December 2016 show better results for MNDWI_var (≈50% FR) with respect to NDVI_var. (≈35% FR). It was not possible to the map the flooded area with Δσ^o of Sentinel-1 because three days after the flood peak, only 4% of the inundated area was detected. For this area, the accuracy in flood detection is also unrelated to spatial resolution; the best performance comes from MNDWI based on the MODIS image of 26 November 2016 (FR ≈ 80%).

4.5. Comparing Flood-Mapping Approaches: Automatic Emergency, Traditional and Combined Approach

For the case of the 2015 Ebro and 2016 NW Italy floods, Figure 14 shows a comparison between the flooded areas mapped by EMS-Copernicus [52,58], the flooded area made with manual interpretation for this work and the official flood maps produced by river authorities.

The EMS-Copernicus mapping is automatic mapping in quasi-real time that identifies critical areas during the emergency but does not represent the real flood extent for further post-emergency analyses. For instance, the EMS flood mapping in the case of the Alessandria area (Figure 14A) is based on Radarsat-2 data (acquired on 27 November 2016), and the map shows only the remaining flooded areas (dark blue). The EMS underestimates the flooded area with respect to the results of this work and official maps.

The area at the south of Turin (Figure 14B) was partially mapped by EMS only from Carignano downstream to Turin. The flood map for the area near Moncalieri, investigated with Pleiades satellite at 26 November 2016, detected most of the flooded areas and fitted well with our map. The sector near Carignano town, mapped with an image of 27 November 2016, shows only residual flooded areas.

In the case of Ebro valley, the EMS flood map was made with a Radarsat-2 image acquired on 2 March 2015 (Figure 14C). The map shows an almost complete fit with our flood mapping and the official flood mapping of CHE.

The comparison with official mapping showed good results. The flood maps made for this study based on remote sensing data match well with the official maps of the Ebro basin authority (CHE) with a 95% overlap (Figure 14C). The area inundated by Tanaro (Figure 14A) proposed in this study also shows a good match, approximately 90%, with ARPA Piemonte flood maps. There are more differences in the area at the south of Turin (approximately 85% matching, Figure 14B); this sector is mapped mostly with post-flood data, and the areas affected by shallow water depth have fewer signs of the flood. Flood mapping based on satellite data presented more challenge in the urban areas, in which a more detailed survey is often required, as shown in a previous work [47].

In the case of the Arahal flood, it is not possible to estimate the accuracy of flood mapping because official maps are not available. However, the combination of Sentinel-1, Landsat-8, and Sentinel-2 data implemented with a water depth model allowed creating a reliable map that was locally validated.

5. Discussion

The results of our study show that choosing the best satellite data or the most appropriate processing aimed at detecting flooded areas is a complicated process. It is necessary to consider the complex relationships between the flood and the satellite characteristics. In the following, we introduce a more detailed description of the parameters that play a fundamental role:

• The time factor is the most critical parameter. If a co-flood image is available, the detection and the refinement of the inundated area are more straightforward and require a short time. To know the time of maximum inundation, data from river gauge stations (e.g., for the Ebro, Po, and Tanaro case studies) can be considered. If river gauge station data are not available (e.g., Arahal case study), it is often possible to estimate the time of inundation from ancillary data found on the web such as news, photos, or videos. Medium-low resolution multispectral satellites (MODIS, Proba-V, Sentinel-3), with their daily revisit time, are more likely to acquire co-flood images for many flood events. However, it is less probable to have co-flood SAR images acquired by Sentinel-1, which in 2018 are the only SAR satellites providing free images. The time elapsed from the flood peak also influences the band ratios to use; MNDWI_var shows better performance in the short term, whereas NDVI_var and VIS_var have better performances over the long term because MNDWI is more sensitive to open water; changes in NDVI and VIS are primarily related to vegetation damage or silt deposits in flooded areas. Concerning the pre-flood images used as a benchmark, the one acquired as close as possible to the flood event should be selected to have similar conditions of land use and sun illumination.
• The use of a single acquisition date with respect to image differences with pre-flood must be considered in the semi-automatic mapping: In the first case, permanent water bodies outside flooded areas are wrongly classified as flooded. In the second case, the water bodies inside the flooded area are not classified as flooded. This issue can be solved using a mask of the permanent water body.
• The sky-condition is another relevant parameter that limits the availability of co-flood images. Indeed, only SAR satellites can acquire images under cloudy conditions. In our cases, the area flooded around Zaragoza was cloud free during the 2015 flooding because the flood peak arrived several days after the meteorological event, whereas cloudy conditions did limit the availability of images for the 2016 Po basin flood.
• The spatial features of the flooded area represent a constraint on satellite spatial resolution. Small, inundated areas (e.g., Arahal case study) cannot be mapped by low-resolution satellites. In addition, an area flooded by a small stream has a very short co-flood time interval.
• The extent of the affected area. For the main river, it might require more than one image to map the whole river section affected by the flood. In addition, the time of maximum flood changes along the stream and is not the same for tributaries. For instance, the flood waves of the Po and Ebro Rivers transited over several days from upstream to downstream.
• Land use and morphology of the flooded area. These factors are important in choosing the satellite image. Inundations occurring on cultivated crops in floodplains can be easily detected by most satellites. In addition, several days after the flood, with multi-spectral data, it is possible to recognize flooded areas easily, particularly if sediments cover crop fields. In urban or densely forested regions, detection of the flooded areas is more complicated. In the arid region, the effect on NDVI can be opposite, with an increase of vegetation index in the flooded area [67]. In many cases, high-resolution data and field surveys remain necessary for a reliable mapping.
• With the semi-automatic approach, we used empirical thresholds to detect flooded areas. The semi-automatic approach allows selecting the best threshold that defines a geomorphological based pattern of the flooded area. Conversely, the semi-automatic approach requires specific operator experience in flood dynamics.
• Results also show that the direct visual mapping of flooding areas can be performed only when the flooded area is visible in co-flood images. Otherwise, it is better to use semi-automatic detection followed by the refinement step.

In the presented study, we used only cost-free satellite data that allowed us to create reliable flood maps. One of the most important features of the used satellite data is their continuous temporal sampling and the rapid revisiting time. For example, MODIS satellites provide multispectral images on a daily basis, which allows having free cloud images close to the flood peak for the Ebro and Po flood areas. The Sentinel-1 system continually provides SAR data with short revisit times and global coverage, and Sentinel-2 satellites allowed the creation of flood maps with very good spatial resolution.

However, if we must map floods affecting dense urban areas with high accuracy or produce an exact flood water depth map for damage assessment, we must exploit high-resolution (satellite, aerial or terrestrial) data. Very often, the required high-resolution data are not free, for example, the 2016 Chisola stream, which flooded the Moncalieri town, as described in previous work [57].

In this work, we also exploited free ancillary data such as DEM, river gauge data, and geolocated information. Our areas of study were covered by free and high-resolution DEMs provided by local authorities and, for the Ebro and Po basins, were covered by gauge station networks providing river level data with a good density. As described earlier, these ancillary datasets are fundamental for a better definition of the limits of the flooded area.

6. Conclusions

In this study, we presented a methodology for low-cost and user-friendly flood mapping based on free satellite data jointly exploited with free ancillary data. In particular, we used data from Sentinel-1, Sentinel-2, Landsat-8, Proba-V and MODIS Terra/Aqua satellites and processed them with open-source software based-methods easily carrying out also by non-remote sensing expert community. We tested our methodology over three areas recently hit by floods: (i) The Ebro valley near Zaragoza, Spain, flooded in spring 2015; (ii) two areas near Alessandria and Turin (Piemonte region, NW Italy) inundated by Po, Tanaro, and their tributaries in November 2016; and (iii) a small area around the Arahal village, southern Spain, affected by a flash flood in November 2017.

The proposed flood mapping strategy is composed of progressive steps. Firstly, the detection of flooded area is performed, through semi-automatic extraction or by visual interpretation of satellite images. The flood map is later manually refined with the help of ancillary data such as digital elevation models (DEMs), water depth (WD) models or ground photos and web-based information. The final outcome is an accurate and morphologically based flood map in vector format.

In the presented case studies, we validate the methodology results with available flooded area maps provided by local authorities. For the case of the Ebro and Po basin flood, we considered a performance assessment of the different automatic maps, and we calculated the flood ratio (FR), which is the percentage of the flooded area detected by satellites with respect to official maps.

The study points out that the capability of spaceborne sensors to map a flood event depends upon several factors: (a) The time of the satellite pass with respect to the time of the flood peak, that is the most important factor; (b) the spatial resolution of satellite; (c) the sensors type (e.g., SAR or multispectral); (d) the sky conditions; (e) the land-use and morphology of the flooded area; (f) the type of data processing.

The results show that SAR-based flood mapping is more efficient when co-flood acquisitions are available, allowing the extraction of information on areas fully covered by water also in cloudy or night conditions. However, the only SAR system currently providing free and open data is Sentinel-1; for this reason, the availability of co-flood data is limited. In our study, a co-flood Sentinel-1 image was available only for the 2015 Ebro valley flood, allowing a satisfactory detection of the flooded area (≈50% FR). In the other study areas, only post-flood SAR data were available; these data, although acquired a few hours or days after the water withdrawal, were not suitable for flood mapping: FR showed values of about 20% for the area flooded by Tanaro to less than 5% for the area inundated by the Po River.

In the case of multispectral data-based mapping instead, the analysis shows that post-flood acquisitions can be profitably exploited, allowing to detect changes inside the flooded areas within a few weeks after water withdrawal. The results show that the FR decreases with the time elapsed from the flood peak but slower respect to SAR: For instance, the computed MNDWI_var decreases from values >90% with co-flood imagery to less than 50% if the available imagery is acquired two weeks later in all the flooded areas. It is worth noting that the local conditions, such as the depositions of silt sediments, contribute to mapping the flooded areas for a more extended period: e.g., VIS_var showed a FR≈ 70% in the area flooded by the Tanaro using Sentinel-2 data acquired 15 days after the floods.

The MODIS and Proba-V satellite sensors, with a daily revisit time, increase the probability of acquiring co-flood data and allow achievement of very good results despite the medium-low spatial resolution. In the Zaragoza area, the MODIS and Proba-V co-flood images allowed a very accurate mapping (FR >90%). For the area of Alessandria, quasi co-flood Proba-V data permitted achieving excellent results, with a FR of nearly 100% for the area flooded by the Tanaro.

The analysis that relied on medium-high resolution Landsat-8 and Sentinel-2 data shows the same flood ratio performance as the medium-low resolution satellites (FR > 90% for co-flood Landsat images in Zaragoza areas) but with a spatial precision several times greater. These data indeed also allowed the detection of small, flooded areas (up to ≈500 m²). In particular, the high-resolution was very useful for the Arahal flash-flood case study. In this case, the availability and the combination of high-resolution Sentinel-2, Sentinel-1 and Landsat-8 data joined with ancillary data and geomorphological interpretation was fundamental to producing a reliable flood map.

The study points out some drawbacks of flood mapping based on the use of SAR and multispectral satellite data. In the latter case, the main limitation is represented by cloud cover affecting co-flood multispectral acquisitions. The presence of clouds has indeed limited the flooded areas detection in the case of the Po flood, while has not affected the Ebro study because the inundation occurred several days after the rainfall event.

The presence of dense urban areas and forests affects both SAR and multispectral based flood mapping and requires a more-complex data processing which is not straightforward to accomplish with a user-friendly approach. High spatial resolution is also a key factor when mapping floods in dense urban areas, and it is one of the limitations of the free of charge satellite data approach. Based on our experience, on-demand high costs, high resolution data and field surveys are often necessary [57].

However, despite such drawbacks the integration of freely available satellite data with ancillary data allowed us to map about 450 km² of inundated areas with high accuracy for all studied areas.

Our methodology proves to be a good trade-off between emergency real-time automatic flood mapping and the traditional flood mapping performed with on-demand data. Real-time automatic flood mapping (e.g., Copernicus EMS) are usually fast, but not very accurate. In contrast, traditional flood mapping performed with on-demand data provides very accurate results, but it is a more resource-intensive process.

Furthermore, we propose a user-friendly approach for satellite data processing that can also be exploited by users with no remote sensing expertise, contributing to hydrological studies, urban planning or flood hazard/damage assessment. In such a context, having a network of operational satellites providing free- and open-access data with continuity is of great importance for flood mapping and hazard/risk management in floodplain areas and disadvantaged regions worldwide.