Volumetric Quantification of Flash Flood Using Microwave Data on a Watershed Scale in Arid Environments, Saudi Arabia

Abstract

Actual flood mapping and quantification in an area provide valuable information for the stakeholder to prevent future losses. This study presents the actual flash flood quantification in Al-Lith Watershed, Saudi Arabia. The study is divided into two steps: first is actual flood mapping using remote sensing data, and the second is the flood volume calculation. Two Sentinel-1 images are processed to map the actual flood, i.e., image from 25 May 2018 (dry condition), and 24 November 2018 (peak flood condition). SNAP software is used for the flood mapping step. During SNAP processing, selecting the backscatter data representing the actual flood in an arid region is challenging. The dB range value from 7.23–14.22 is believed to represent the flood. In GIS software, the flood map result is converted into polygon to define the flood boundary. The flood boundary that is overlaid with Digital Elevation Map (DEM) is filled with the same elevation value. The Focal Statistics neighborhood method with three iterations is used to generate the flood surface elevation inside the flood boundary. The raster contains depth information is derived by subtraction of the flood surface elevation with DEM. Several steps are carried out to minimize the overcalculation outside the flood boundary. The flood volume can be derived by the multiplication of flood depth points with each cell size area. The flash flood volume in Al-Lith Watershed on 24 November 2018 is 155,507,439 m³. Validity checks are performed by comparing it with other studies, and the result shows that the number is reliable.

1. Introduction

In the recent century, floods have become one of the most devastating disasters and their occurrence has increased frequently due to climate change [1,2,3]. Floods can happen because they are influenced by the chain of hydrological and meteorological conditions [4].

Saudi Arabia is an arid region, characterized by low rainfall annually and no permanent natural water bodies, such as rivers and lakes [5]. One common disaster in the arid region is flash floods, which is influenced by less vegetation and an abundance of loose materials [6]. The loose sediment materials and lack of cohesivity on the soils can cause extreme sediment fluxes during a flash flood in an arid region [7,8]. In the arid region, sometimes water caused by extreme rainfall cannot penetrate easily into the ground due to impermeable rock (e.g., caliche) and causes high runoff water [9,10,11].

Flash floods in arid regions have caused numerous fatalities and property losses [12,13,14]. In Saudi Arabia, flash floods are included in the top 10 natural catastrophes. For example, during a devastating flood in Jeddah that occurred on November 25, 2009, 161 fatalities happened and required USD 900 million to reconstruct Jeddah city [15,16]. The flood risk in the arid region is increasing due to urbanization expansion [17,18,19]. Therefore, it is crucial to calculate the authentic flood volume to prevent more losses caused by the same disaster. This information will be beneficial for the stakeholder to make a proper recommendation and prevention program. In this study, the area we focused on is the flash flood disaster in Al-Lith City on 24 November 2018. Al-Lith City is a coastal area where it is known that most of its floods occur at the river and coastal area [20,21]. The recent flash flood hazard map in Al-Lith Watershed shows that topography and slope have an important role. The higher the elevation and slope, the higher the hazard will occur [22,23].

Remote sensing is a useful technique to assess Earth’s environment on a large or regional scale [24,25]. Sentinel-1 is a remote sensing image generated by the Sentinel-1 mission belonging to the European Space Agency (ESA). The Sentinel-1 image is now commonly used by many researchers in many countries to map the actual flood. The inland water bodies could be detected after thorough processing of certain algorithms using infrared (IR), short-wave infrared (SWIR) bands [26], and global navigation satellite system reflectometry (GNSS-R) [27,28,29]. Kouassi et al. [30] mapped the devastating flood in South-west Cote d’Ivoire during the high rainfall season in 2017. Tavus et al. [31] used the Sentinel-1 data to map the flood in Ankara, Turkey, and provide high accuracy. Twele et al. [32] mapped the flood in two locations (i.e., Greece and Turkey) using Sentinel-1 data and concluded that VV data give more accuracy than VH polarization data. When trying to map the flood in the arid region, there is a challenge because the flat sand will be identified as water by the satellite. Flood mapping in the arid region has been carried out by Elhag and Abdurahman [33] using the Sentinel-1 image, and it is concluded that the sand feature has a similar backscatter with the flooded area. Therefore, it required effort in a trial and error fashion when selecting the most representative backscatter value.

A flood depth quantification study has been performed by Cham et al. [34] andCohen et al. [35]. The idea is to use the flood map data generated from the remote sensing technique to determine the flood depth distribution within the inundation boundary. These two studies utilized different methods to define flood surface elevation inside of the flood boundary. Cham, Mitani, Fujii and Ikemi [34] used the Triangulated Irregular Network (TIN) in the 3D Analyst Tool, and Cohen, Brakenridge, Kettner, Bates, Nelson, McDonald, Huang, Munasinghe and Zhang [35] used the Focal Statistics method in the Spatial Analyst Tool. In this study, we decided to use the Focal Statistics method because Focal Statistics is a neighborhood method and it fills the value of the empty cells based on the adjacent cells’ value.

GIS software is used for this analysis. By using flood depth distribution data, the flood volume can be calculated. The concept of flood volume calculation is the multiplication of flood depth distribution with the raster cell size area. In this study, the Digital Elevation Map (DEM) will be used to assist the calculation [36]. The DEM source is from the American–Japanese satellite, ASTER GDEM, with a resolution value of ±30 m.

The aims of this study are to quantify the flood map area (based on remote sensing technique) and the flood volume, which happened in Al-Lith Watershed, on 24 November 2018. In the Study Area section, details on the location and watershed parameter will be explained. The step-by-step processes from remote sensing analysis (actual flood mapping) to the flood volume calculation will be in the Methodological Framework section. In the Result and Discussion section, there are more complex scientific discussions, such as: how to select dB values to represent the flood in arid region; how to use the focal statistics to interpolate water surface elevation; how to calculate the flood distribution raster; how to calculate flood volume. In the end, statistical analysis is carried out to remove the outliers and several validation analyses are conducted to ensure that the result is reliable; therefore, this study can be followed by any studies in another area.

2. Study Area

Al-Lith Watershed is located in the west of Saudi Arabia (Figure 1). It is part of Mecca Province and about 200 km from Jeddah City. It lies between 40°10′ E and 40°50′ E and 20°00 N and 21°15′ N. By using ArcGIS software, several morphometric parameters are calculated (

Table 1. Al-Lith Watershed morphometric parameters.

Parameter	Result	References
Watershed Area (km2)	3209	Schumm [37]
Watershed Perimeter (km)	624.6	Schumm [37]
Streams Length (km)	2613	Horton [38]
Elevation Range (m asl)	0–2663

). The watershed covers 3209 km² and 624.6 km perimeter. This watershed’s elevation range is between 0 and 2663 m above sea level (asl), and the total stream’s length is 2613 km.

The geomorphology condition (Figure 1) shows that the stream’s directions range from the eastern high mountainous zone to the western flat sediment area [39]. At the end of this watershed (i.e., coastal area), there is a city called Al-Lith City, in which water streams from the entire Al-Lith Watershed will be flowing to the city [40].

Inside the watershed, Al-Lith Dam is collecting water streams from the north before entering Al-Lith City. The study will be more focused on quantifying the flood below the dam or downstream area because it is believed that the water above the dam is conservative.

3. Methodological Framework

This study’s methodology is divided into two steps: flash flood mapping and flood volume quantification (Figure 2). The first step is to derive the actual flood map in an area using the remote sensing method. The second step explains how to use the existing flood map for flood volume calculation using GIS software.

3.1. Flash Flood Mapping

This step’s material requirements are two Sentinel-1 images representing the condition before the flood (archive image) and the condition during the peak of the flood (crisis image) [41]. Sentinel-1 image is free remote sensing data generated by the European Space Agency (ESA). In this study, the archive image and the crisis image acquisition are 28 May 2018 and 24 November 2018, respectively. The software utilized for image processing is the Sentinel Application Platform (SNAP) software. SNAP software is an open-access software that is also from ESA and able to generate the flood map. After the Sentinel-1 images are derived, these two Sentinel-1 images are cropped into a smaller area, representing the study area sufficiently; therefore, image processing will be faster (Figure 3).

• Calibration and radiometric correction

Image calibration and radiometric correction are essential for comparing two images and transforming the digital numbers into physical quantities [42,43]. The backscatter data with VV polarization are used for this process, because they are more suitable for inland water detection than VH polarization [44,45], and then the σ₀_VV is generated. The σ₀_VV data (linear scale) are converted into dB (log scale) for color manipulation purposes [46]. Thus, the image will have better visualization, and the histogram is more comfortable to be manipulated (Figure 4). The equation for this conversion is given by:

$${\mathsf{\sigma }}_0\left(\mathrm{dB}\right)=10\times log10\left[abs\left({\mathsf{\sigma }}_0\right)\right]$$

(1)

• Speckle filtering

Lee filter is applied for the σ₀_VV_db image to minimize the speckle [47]. The result will be smoother compared to the non-filtered image (Figure 5).

• Terrain correction

Terrain correction is performed to project the images into a map system’s coordinate and correct the terrain’s distortion.

• Create stack and RGB composite

Combine the images and create an RGB composite to distinguish the flood area and non-flood area. The RGB composite image should be used as a guide to identify the flood area [48,49]. The red color represents the flood near the permanent water body, and the blue color represents the flood in the non-permanent water body (e.g., seasonal streams). The actual flood can be predicted by selecting the dB value representing the actual flood. After a representative dB value is selected, the flood map should be exported into a raster format (e.g., GeoTIFF) as a requirement for the flood quantification process.

3.2. Flood Volume Quantification

The materials required for this step are Flood Map Raster and DEM. The DEM data must be processed with Sinks Filled Tools. The concept of flood volume quantification is to multiply the flood depth raster with the raster’s cell size area. Therefore, it is required to derive the raster image containing the flood depth data on each cell. The software used for this step is ArcGIS.

• Define the flood boundary

By using Conversion Tools, the flood map raster is converted into a polygon to acquire the flood boundary, and then the flood boundary polygon is converted into a polyline. This flood boundary polyline is converted into a raster to process further analysis. By using Raster Calculator, the value of the flood boundary raster is added with the DEM data.

• Define the flood surface elevation

After the flood boundary raster is generated, the next step is to find the flood surface elevation inside each boundary. The Focal Statistics method is used with statistical type using the mean data, circle geometry, and mean value. The circle with a 1-cell radius is used and believed to give a better result. Some iterations of Focal Statistics are required until all empty cells inside the flood boundary are filled. Then, the result is cropped by the flood boundary to remove overcalculation outside the flood boundary.

• Define the flood depth distribution

The flood surface raster is subtracted by the DEM using Raster Calculator to derive the flood depth distribution. The subtraction result raster’s negative value must be eliminated by using Con Tools in Spatial Analyst Tools. Then, Filter Tools with a smoothing method (low pass) is used to minimize the sharp change in flood depth spatial distribution [35]. The result is clipped by the flood boundary to ensure that the Filter’s result does not exceed the flood boundary. At this step, a proper flood depth raster can be generated.

• Flood volume calculation

Using Conversion Tools, the flood depth raster converted into a point (shapefile), in which this point will have a flood depth elevation value. Inside the attribute table of flood depth point shapefile, this value is multiplied with the raster size (e.g., 30.1595 m × 30.1595 m). The result represents the flood volume per-point or per-cell. In the end, by summation of all volumes, the flood volume is derived.

4. Result and Discussion

4.1. Selecting Representative dB Value from Crisis Image

After creating a stack and RGB composite (Figure 5A), it is required to select the dB value representing the flood. In an arid region, selecting this dB value is a challenge according to the following reasons:

• The flat sand morphology will be detected as water because it has similar backscatter data.
• There is almost no permanent river or lake; therefore, the condition before the flood occurred is dry. Water content inside the soil will cause stronger backscattering because the SAR signal penetration ability will become stronger than dryer sand. Due to this reason, the backscatter during rainy conditions will tend to increase [50].

It is believed that the dB value, which represents the flood is a range between 7.23 and 14.22, which is selected from the crisis image by using the following equation:

255 * (if σ₀_VV_db_slv1_24Nov2018 <−7.23 and σ₀_VV_db_slv1_24Nov2018 >−14.22 then 1 else 0)

(2)

This formula will make the desired dB value become 1 leaving the other values to be 0.

Another problem is that the backscatter dB value within the range 7.23–14.22 is not only reflecting the flood in the downstream area but also reflecting the mountainous area. Due to this reason, the dB range within 7.23–14.22 is only applied for the area below the dam and the downstream area. This value is not applied to the mountainous and upstream areas. The other consideration is because the water in the area above the dam is considered to be conservative. The result can be seen in Figure 6, where the total flood area that occurred in Al-Lith City on 24 November 2018 is 125,169,332 m² or 125.17 km².

4.2. Flood Surface Interpolation

Figure 7 shows the conversion result from the flood map raster into shapefile (polygon or polyline). Then, the flood boundary polyline is converted into a raster. Using DEM (sinks filled), the value of flood boundary raster is added with elevation unit meter above sea level (m asl). It is crucial to ensure that both raster cell sizes have the same values (i.e., 30.1595 m × 30.1595 m). To do this step, the resample method from Data Management Tools is used.

Focal Statistics in Neighborhood Tools is used to determine the inundation surface elevation inside the flood boundary. In this step, circle geometry with a 1 m radius and 3 iterations is required to fill all empty cells inside the boundary (Figure 8).

The Focal Statistics result is clipped by the flood boundary to eliminate the overcalculation outside the flood boundary (Figure 9).

4.3. Flood Depth Distribution

In determining the flood depth raster, the raster of flood surface elevation should be subtracted by the DEM (sinks filled). Figure 10A shows that the calculation results in a negative value. The overcalculation results cause a negative value during the interpolation step. Due to this circumstance, the negative value should be neglected using Con in Spatial Analyst Tools (Figure 10B). Then, in Spatial Analyst Tools, a Filter for smoothing (low pass) is used to minimize the sharp change within flood depth spatial distribution (Figure 10C). The Filter is one of the interpolation methods that possibly calculate the value outside the flood boundary; therefore, it must ensure that the result is not exceeding the flood boundary. Thus, the Filter result is cropped by the flood boundary (Figure 10D).

The flood distribution raster is converted into points (shapefile), which contain the flood depth information (meter). The statistical information can be seen in Figure 11A, where it has a range from 0.001 to 11.02 m, an average (mean) flood depth of 1.17 m, and a median of 1.03 m. This information is gathered from 161,064 points, where each point represents 30.1595 m × 30.1595 m cell size. This result is be analyzed statistically to remove the outliers. One of the common methods that is used to eliminate the outliers is the Interquartile Range (IQR) technique [51,52]. The equation for the upper boundary is given by:

$$\begin{gathered} Upper Boundary=Q3+1.5\left(Q3-Q1\right) \\ Upper Boundary=1.57+\left(1.5\times 0.97\right)=3.015 \end{gathered}$$

(3)

The result shows that values higher than 3.015 are justified as outliers and will be eliminated by using Con in Spatial Analyst Tools. The final statistical information after the outlier’s removal can be seen in Figure 11B, where now it has a range of 0.001–3.014 m, an average (mean) flood depth of 1.09 m, and a median of 1.01 m, gathered from 156,363 points. Figure 12 shows the final flood depth distribution raster in Al-Lith Watershed. Most of the highest flood depths are occurred surround the streams, especially on the upstream areas.

4.4. Flood Volume Quantification

The flood volume can be calculated on the attribute table of flood depth points using a field calculator. The flood depth is multiplied by the flood depth raster’s cell size area (30.1595 m × 30.1595 m). This multiplication summation is 155,507,439 m³, representing the total volume of the flood in Al-Lith Watershed. This result is compared with a study performed by other researchers for validation purposes.

Ewea et al. [53] analyzed the maximum flood volume of every watershed in Saudi Arabia. They made the envelope curves of flood volume versus watershed area (km²). The equation for maximum flood volume is given by:

$$V_{max}=0.0004A^{1.7}$$

(4)

where V_max is the maximum flood volume in 10⁶ m³, and A is a watershed area in km². By using the Al-Lith Watershed area (i.e., 3209 km²), the maximum flood volume is 365,500,950.2 m³. Based on this information, the result of this study (155,507,439 m³) does not exceed the value of maximum flood volume.

Another study was carried out by Ekraim et al. [23]. They calculated the flood volume in Al-Lith Watershed using the HEC-HMS model. Their study revealed that in a 20-year and 100-year storm return period, Al-Lith Watershed could produce 151,132,255.20 m³ and 208,212,192.9 m³ flood volume, respectively. In this watershed, the average rainfall that occurred on 23 November 2018 (which caused the flash flood on 24 November 2018) was 104 mm/day [54], which can be classified as a 20-year storm return period. This piece of information also emphasizes that this study’s result is matched and reliable.

5. Conclusions

In this study, there are two steps to quantify the actual flood in such an area; the first step is flash flood mapping using remote sensing technique, and the second step is calculating the flood volume. During Sentinel-1 image processing in SNAP software, there is a challenge in selecting the backscatter value representing the flood in an arid region. By using the RGB composite maps as a guide to select the backscatter value, the flood occurrences can be estimated. In this study, the dB range values that represent the flood are 7.23–14.22 based on the crisis image. The result shows that the flood area is 125,169,332 m² or 125.17 km². The interpolation technique used to define the flood surface within flood boundary is Focal Statistics, where it is believed to give a reliable result. Three iterations are required for the Focal Statistics processes until all empty cells are filled. Several steps are carried out to minimize the overcalculation outside the flood boundary, i.e., clipped by the flood boundary polygon and negative values removal using Spatial Analyst Tools. The Filter for smoothing purposes is used to decrease the spike change in the flood depth distribution. The flood depth distribution shows that the range of the flood depth is between 0.001 and 3.014 m, with an average depth is 1.09 m and a median is 1.01 m. The total volume of the flash flood in Al-Lith Watershed happened on 24 November 2020: 155,507,439 m³. This number agrees with the result of other studies in the same area, which used different methods. It is concluded that the methodology used in this study is reliable and can be used in another area.