1  Detect flood

The downloaded flood extent covered a larger area in Rio Grande do Sul. However, the area of interest containing the urban dense cities such as Porto Alegre was smaller. Therefore, a series of operations reduced the large area focusing on Porto Alegre, while improving the query performance.

A positive aspect of these set of operations in the methodology is to allow demanding analysis without expensive hardware (See computational environment), which indirectly means a lower consume of energy and environment impact. If a still larger area of interest is required, we carried out tests with duckdb and its extension spatial ¹ with sucess.

library(sf)
library(DBI)
# Load the flood extent and remove Z dimension
flooding_subdivided_porto <- st_zm(sf::st_read(connection,
                                        DBI::Id(schema="agile_gis_2025_rs",
                                            table="flooding_subdivided_porto")))
core_metropolitan_area <- st_zm(sf::st_read(connection,
                                        DBI::Id(schema="agile_gis_2025_rs",
                                            table="core_metropolitan_area")))
# Create leaflet product
library(leaflet.extras)
library(leaflet)
leaflet() |>
  addTiles()  |>
  addProviderTiles("OpenStreetMap", group = "OpenStreetMap")  |>
  addProviderTiles("Esri.WorldImagery", group = "Esri.WorldImagery")  |>                    addPolygons(data=flooding_subdivided_porto,
                    fillColor = "#dec8b7ff",
                    fillOpacity = 0.6,
                    weight = 1,
                    color = "#4b2609",
                    dashArray = "1",
                    group= "Flood extent") |> 
  addPolygons(data=core_metropolitan_area,
                    fillColor = "yellow",
                    fillOpacity = 0.05,
                    weight = 2,
                    color = "red",
                    dashArray = "1",
                    group= "Core Metropolitan Area") |> 
 addLayersControl(
    baseGroups = c("OpenStreetMap", "Esri.WorldImagery"),
    overlayGroups = c("Flood extent","Core Metropolitan Area"),
    options = layersControlOptions(collapsed = FALSE)
  )

These set of operations are grouped into the section subdivide the flood extent and simplify geometry.

1.1 Subdivide the flood extent

The code chunk “[AoI-6] Dect flood extent” includes one query with several CTE to create the table flooding_subdivided_porto. Firstly, we subdivided the extent of the flooding following the post published by the co-founder of PostGIS ², Paul Ramsay.

CREATE TABLE flooding_subdivided_porto AS
 WITH porto_alegre_ghs AS(
    SELECT
      ghs.*
    FROM
      urban_center_4326 AS ghs
    JOIN
      nuts
    ON
      st_intersects(nuts.wkb_geometry, ghs.geom)
  WHERE
    nuts.shapename = 'Porto Alegre'),
--- Select core metropolitan areas in the municipalities that intersects with ghs
core_metropolitan_area_ghs AS (
    SELECT 
        st_union(wkb_geometry) AS core_geom
    FROM 
        nuts
    JOIN 
        porto_alegre_ghs
    ON 
        st_intersects(nuts.wkb_geometry,
                      porto_alegre_ghs.geom)),
---Subdivide the flood extent to increase spatial indexes performance
flooding_sul_subdivided AS (
    SELECT
      st_subdivide(geom) as the_geom
    FROM
      flooding_raw)
/* Select the flooding subunits that intersects with
the previous bounding box */
    SELECT
      flooding_sul_subdivided.*
    FROM
      flooding_sul_subdivided,
      core_metropolitan_area_ghs
    WHERE
      st_intersects(flooding_sul_subdivided.the_geom,
                    core_metropolitan_area_ghs.core_geom);
--- time: 39099.137 ms

1.2 Simplify geometry

The reason for simplifying the geometry by removing the slivers polygons was that spatial operations such as st_difference, which are required to obtain the post-disaster network, are computationally expensive. Based on the Carto’s post, we removed the slivers polygons using the calculated area ³. Another technique in PostGIS to remove these slivers polygons is via dilation and erosion ^{4 5}, which could be used in future studies. The table flooding_cleaned_porto_union_simple is a lighter version of the flood extent

# Reduce geometry
CREATE INDEX idx_flooding_subdivided_porto ON flooding_subdivided_porto USING gist(the_geom);
EXPLAIN ANALYZE
CREATE TABLE flooding_cleaned_porto AS
  SELECT
    ST_Force2D((ST_Dump(the_geom)).geom)::geometry(Polygon, 4326) AS geom
FROM flooding_subdivided_porto;
--- time: 25.746 ms
    
# Unite
CREATE INDEX idx_flooding_cleaned_porto ON flooding_cleaned_porto USING gist(geom);
EXPLAIN ANALYZE
CREATE TABLE flooding_cleaned_porto_union AS
  SELECT
    ST_Union(geom)::geometry(multipolygon, 4326) geom
  FROM
    flooding_cleaned_porto;
--- time: 348.667 ms

# Reduce geometries again
CCREATE INDEX idx_flooding_cleaned_porto_union ON flooding_cleaned_porto_union USING gist(geom);
EXPLAIN ANALYZE
CREATE TABLE flooding_cleaned_porto_union_simple AS
  SELECT
    (ST_Dump(geom)).geom::geometry(polygon, 4326) geom
  FROM
    flooding_cleaned_porto_union;
---time:  16.893 ms
/*This allowed to calculate the area of each polygon finding
slivers that were removed using the following code.*/
SELECT COUNT(*)
FROM flooding_cleaned_porto_union_simple ; --- count= 99

DELETE FROM
flooding_cleaned_porto_union_simple
WHERE
ST_Area(geom) < 0.0001; ---count= 4.

1.3 Visualize flood extent

The following map shows the difference between the original flood extent subdivided and delimited by the AoI and the flood extent cleaned after removing the slivers polygons. The purpose is to check the results from a qualitative point of view.

# Load the flood extent and remove Z dimension
flooding_cleaned_porto_union_simple <- st_zm(sf::st_read(connection,
                                        DBI::Id(schema="agile_gis_2025_rs",
                                            table="flooding_cleaned_porto_union_simple")))
flooding_all <- st_zm(sf::st_read(connection,
                                        DBI::Id(schema="agile_gis_2025_rs",
                                            table="flooding_raw")))
# Visualize the raw flood extent and its simplified version
library(mapview)
library(leafsync)
m1_flood <- mapview(flooding_subdivided_porto, 
              color="#f1a340",
              col.region="#f1a340",
              layer.name ="Raw flooding extent")
m2_flood <- mapview(flooding_cleaned_porto_union_simple,
              color ="#998ec3",
              col.region="#998ec3",
              layer.name="Simplified flooding extent")
m1_flood | m2_flood 

Meanwhile, the following bar chart shows the difference on the size of the flood extent. Originally the size of the flood extent was nearly 20 MB covering an area larger than the AoI. After being reduced to the AoI, the flood extent reduced its size to nearly 664 kB Lastly, simplifying the geometry of the flood extent in the AoI reduced its size to 464 kB.

Fig.1 Flowchart

# Compare size of both versions
## Create dataframe
library(tidyverse)
df_flood <- data_frame(
    version = c("1.Original",
                "2.AoI",
                "3.AoI simplified"),
    size_bytes = c( object.size(flooding_all)/1000,
                    object.size(flooding_subdivided_porto)/1000,
                    object.size(flooding_cleaned_porto_union_simple)/1000)
) 

## Custom Y-axis labels 
labels <- function(x) {
  paste(x, "kB")
}
## Create column plot
library(ggplot2)
library(ggbreak)
flood_plot <- ggplot(data=df_flood) +
  ggplot2::geom_col(aes(x=version, y=size_bytes, fill=version)) +
  labs(x="Flood extent",
       y="Size in kBytes") +
  theme_minimal() +
  theme(legend.position = "bottom") +
  scale_fill_manual(values= c("#e66101",
                              "#f1a340",
                              "#998ec3")) +
  scale_y_break(c(664,20000), ticklabels = c(df_flood$size_bytes[3],
                                             df_flood$size_bytes[2],
                                             round(df_flood$size_bytes[1])))
flood_plot + scale_y_continuous(label = labels)

---

1. DuckDB spatial: https://duckdb.org/docs/stable/extensions/spatial/overview.html

2. https://blog.cleverelephant.ca/2019/11/subdivide.html

3. https://icarto.es/en/improve-performance-making-the-difference-between-two-layers-with-postgis/

4. https://trac.osgeo.org/postgis/wiki/UsersWikiBufferCoast

5. https://cdn2.hubspot.net/hubfs/2283855/PostGIS%20Day%202019%20-%20Overview.pdf