Simulating Flood Inundation with Python and Elevation Data: A Beginner’s Guide

have become more frequent and devastating across the globe, with the effects of climate change in recent decades. In this context, flood modeling has an important role in risk assessment and disaster-response operations, while also remaining a key focus of advanced research and academic studies.

In this article, we’ll build a basic flood inundation model using Python and a Digital Elevation Model (DEM). We’ll use a flood fill technique to incrementally simulate how rising water levels affect a landscape and animate the inundation process. It’s a visual and hands-on way to explore geospatial data and flood risks, even without a background in hydraulic modeling.

What you’ll Learn

1. What is a Digital Elevation Model (DEM)

A Digital Elevation Model (DEM) is a numerical representation of the Earth’s surface, where each cell (or pixel) in a regular grid (known as raster data) contains an elevation value. Unlike digital images that store color information, DEMs store height data, typically excluding surface features like vegetation, buildings, and other man-made structures.

DEMs are commonly used in fields such as mapping, hydrology, environmental monitoring, and earth sciences. They serve as a foundational dataset for any application that requires a detailed understanding of terrain and elevation.

Many sources offer free and reliable DEM data, including the USGS National Map, NASA Earthdata, and the Shuttle Radar Topography Mission (SRTM).

In this article, we’ll be using a DEM provided by the USGS National Geospatial Program, which is freely available and released into the public domain.

Note: The data provided by USGS has a spatial resolution of 1 arc second (approximately 30 meters at the equator).

The area of interest (AOI) in this study is located in the Northeast region of Brazil. The DEM file covers a 1° × 1° tile, extending from 6°S, 39°W to 5°S, 38°W, and utilizes the WGS84 coordinate system (EPSG: 4326), as illustrated below.

2. How to load and visualize elevation data with Python

Now we’ll use Python to set a viable environment to visualize and analyze some initial information about DEM data. First, let’s import the necessary libraries.

# import libraries
import rasterio
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.animation import FuncAnimation

• rasterio: Reads and writes geospatial raster data like DEMs.

• matplotlib.pyplot: Creates static and interactive visualizations.

• numpy: Handles numerical operations and array-based data.

• FuncAnimation: Generates animations by updating plots frame by frame.

Next, let’s use the rasterio library to open and visualize a DEM file of the AOI.

# Helper function to load DEM Files
def load_dem(path):
    with rasterio.open(path) as src:
        dem = src.read(1)
        transform = src.transform
        nodata = src.nodata

        if nodata is not None:
            # Mask no-data values
            dem = np.ma.masked_equal(dem, nodata)

        return dem, transform

The function above reads the elevation data and checks whether the file includes “no-data values”. No-data values are used to represent areas without valid elevation data (e.g., outside coverage or corrupted pixels). If a no-data value is present, the function replaces those pixels with np.nan, making it easier to handle or ignore them in later analysis and visualizations.

Visualizing DEM data

dem = load_dem("s06_w039_1arc_v3.tif")

plt.imshow(dem, cmap='terrain')
plt.title("Digital Elevation Model")
plt.colorbar(label="Elevation (m)")
plt.show()

• Using geographic coordinates in the visualization

As we can see, the axes are in pixel coordinates (columns and lines). To better understand flood inundation, it’s vital to know the geographic coordinates (latitude and longitude) associated with each pixel of the image.

To achieve that, we’ll use the coordinate reference system data of the DEM file. As said earlier, the DEM we are using utilizes the WGS84 coordinate system (EPSG: 4326).

We can adapt the helper function to load DEM files as follows:

def load_dem(path):
    with rasterio.open(path) as src:
        dem = src.read(1)
        transform = src.transform
        nodata = src.nodata

        if nodata is not None:
            # Mask nodata values
            dem = np.ma.masked_equal(dem, nodata)

        return dem, transform

The function retrieves the transform data from the DEM, which is an affine object that maps pixel positions (row, column) to geographic coordinates (latitude and longitude).

To represent the geographic coordinates on the axes of the plot, it’ll be necessary to explore the extent parameter from the imshow() function.

dem, transform = load_dem("s06_w039_1arc_v3.tif")

# Compute extent from transform
extent = [
    transform[2],                          # xmin (longitude)
    transform[2] + transform[0] * dem.shape[1],  # xmax
    transform[5] + transform[4] * dem.shape[0],  # ymin (latitude)
    transform[5]                          # ymax
]

# Plot with using geographic coordinates
fig, ax = plt.subplots()
img = ax.imshow(dem, cmap='terrain', extent=extent, origin='upper')
ax.set_xlabel('Longitude')
ax.set_ylabel('Latitude')
plt.colorbar(img, label='Elevation (m)')
plt.title('DEM Visualization')
plt.show()

The extent parameter will be used to define the spatial bounds of the DEM plot using values derived from the raster’s transform affine object. It sets the minimum and maximum longitude (xmin, xmax) and latitude (ymin, ymax) so that the plot shows coordinates on the axes instead of pixel indices.

Finally, we have the following results:

3. How to simulate flood scenarios with elevation thresholds

Now, we’ll demonstrate a simple yet useful method for visualizing flood scenarios and simulating inundation. It consists of defining a height threshold and generating a binary mask that identifies all areas with elevation below this level.

In this example, we simulate Flooding across all areas with elevations below 40 meters.

flood_threshold = 40  # meters
flood_mask = (dem <= flood_threshold).astype(int)

plt.imshow(flood_mask, extent=extent, cmap='Blues')
plt.title(f"Flooded Area (Threshold: {flood_threshold}m)")
plt.xlabel("Longitude")
plt.ylabel("Latitude")
plt.show()

With just a few lines of code, we can visualize the impact of different flood scenarios on the area of interest (AOI). However, because this visualization is static, it doesn’t show how the flood progresses over time. To deal with that, we’ll use matplotlib’s FuncAnimation to create a dynamic visualization.

4. How to animate flood progression with Python

We’ll now simulate a progressive flood scenario by increasing the water level incrementally and generating a new mask at each step. We’ll overlay this mask on the terrain image and animate it.

# flood_levels defines how high the flood rises per frame
flood_levels = np.arange(15, 100, 5)

# Set up figure and axes
fig, ax = plt.subplots()
img = ax.imshow(dem, cmap='terrain', extent=extent, origin='upper')
flood_overlay = ax.imshow(np.zeros_like(dem), cmap='Blues', alpha=0.4, extent=extent, origin='upper')
title = ax.set_title("")
ax.set_xlabel("Longitude")
ax.set_ylabel("Latitude")

# Animation function
def update(frame):
    level = flood_levels[frame]
    mask = np.where(dem <= level, 1, np.nan)
    flood_overlay.set_data(mask)
    title.set_text(f"Flood Level: {level} m")
    return flood_overlay, title

# Create animation
ani = FuncAnimation(fig, update, frames=len(flood_levels), interval=300, blit=True)
plt.tight_layout()
plt.show()

# save the output as a gif
ani.save("flood_simulation.gif", writer='pillow', fps=5)

If you’re interested in creating animations with Python, this step-by-step tutorial is a great place to start.

Conclusion and next steps

In this article, we created a basic workflow to perform a flood simulation in Python using elevation data from a DEM file. Of course, this model doesn’t implement the most advanced techniques on the subject, still for visualization and communication, this elevation threshold method offers a powerful and accessible entry point.

More advanced simulation techniques include:

• Connectivity-based flood propagation
• Flow direction and accumulation
• Time-based flow modeling

Nevertheless, this hands-on approach can be of great benefit for educators, students, and analysts exploring Geospatial Data in disaster response studies and environmental modeling.

The complete code is available here.

I strongly encourage readers to experiment with the code using their own elevation data, adapt it to their specific context, and explore ways to enhance or expand the approach.

References

[1] U.S. Geological Survey. National Map. U.S. Department of the Interior. Retrieved May 17, 2025, from https://www.usgs.gov/programs/national-geospatial-program/national-map

[2] U.S. Geological Survey. What is a digital elevation model (DEM)? U.S. Department of the Interior. Retrieved May 17, 2025, from https://www.usgs.gov/faqs/what-a-digital-elevation-model-dem

[3] Gillies, S. Georeferencing — Rasterio documentation (stable). Rasterio. Retrieved May 27, 2025, from https://rasterio.readthedocs.io/en/stable/topics/georeferencing.html

[4] Gillies, Sean. Affine Transforms — Rasterio Documentation (Latest). Accessed May 27, 2025. https://rasterio.readthedocs.io/en/latest/topics/transforms.html.

Data Source: DEM data used in this project is provided by the U.S. Geological Survey (USGS) through the National Map and is in the public domain.