Floods through the clouds — radar flood mapping with Sentinel-1

Take-home notebook · Workshop 1 · IC-TWEP 2026 · QuestFeed

Floods come with clouds, and clouds blind optical satellites. Sentinel-1 carries a radar that sees through cloud, day and night. Calm open water reflects the radar pulse away from the satellite like a mirror, so it appears dark; land, vegetation and buildings scatter energy back and appear brighter.

The event. In September–October 2022 Nigeria suffered its worst floods in a decade: more than 600 people died and about 2.5 million were affected in 31 of the country's 36 states (government figures, October 2022). Months of heavier-than-usual rain across the Niger and Benue basins were intensified by water released from dams upstream, including Lagdo Dam on the Benue in Cameroon, which began releasing on 13 September. The two rivers meet at Lokoja, in Kogi State — one of the hardest-hit areas.

A transboundary river system. Above Lokoja, the Niger and the Benue drain about 2.1 million km² across 11 countries — from Guinea, where the Niger rises, to Cameroon, where the Benue begins. Mali (27 %), Nigeria (25 %), Niger (24 %) hold most of it; Algeria holds 8 %, in the dry Sahara.

Data: Sentinel-1 radiometrically terrain-corrected backscatter and JRC Global Surface Water (Microsoft Planetary Computer), Copernicus DEM (Earth Search on AWS), all open; the catchment outline comes from the workshop data store. No account is needed.

Ready.

# @title 1 · The river system above Lokoja
catchment = fetch_vector("flood_lokoja/catchment.geojson")
m = folium.Map(tiles=None, control_scale=True, height=430)
folium.TileLayer(tiles=ESRI, attr=ESRI_ATTR, name="Satellite imagery").add_to(m)
folium.GeoJson(catchment.__geo_interface__, name="Niger and Benue catchment above Lokoja",
               style_function=lambda f: {"color": "#FACC15", "weight": 2, "fillColor": "#FACC15", "fillOpacity": 0.08}).add_to(m)
folium.Marker([7.795, 6.745], tooltip="Lokoja: the Niger–Benue confluence").add_to(m)
w, s, e, n = catchment.total_bounds
m.fit_bounds([[s, w], [n, e]])
display(m)
print(f"Catchment above Lokoja: {int(catchment['area_km2'].iloc[0]):,} km² — rain anywhere inside the yellow line ends up at the confluence.")

Make this Notebook Trusted to load map: File -> Trust Notebook

Catchment above Lokoja: 2,086,091 km² — rain anywhere inside the yellow line ends up at the confluence.

1 — Area, period and one viewing geometry

We use one relative orbit so that every image is taken from the same viewing geometry — otherwise changes in viewing angle could look like changes on the ground. The notebook picks the orbit that covers the whole area most often.

# @title 2 · Choose the area and period
BBOX = [6.45, 7.55, 7.05, 8.15]          # the Niger–Benue confluence at Lokoja, Nigeria  (west, south, east, north)
START, END = "2022-06-01", "2022-11-30"
BASELINE_END = "2022-08-01"              # images before this date describe normal river levels
RES = 40                                  # metres; 20 for more detail (slower)

lon_c, lat_c = (BBOX[0] + BBOX[2]) / 2, (BBOX[1] + BBOX[3]) / 2
CRS = f"EPSG:{(32600 if lat_c >= 0 else 32700) + int((lon_c + 180) // 6) + 1}"
catalog = Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=pc.sign_inplace)
found = list(catalog.search(collections=["sentinel-1-rtc"], bbox=BBOX, datetime=f"{START}/{END}").items())
aoi = box(*BBOX)
days = defaultdict(list)
for it in found:
    days[(it.properties["sat:relative_orbit"], it.datetime.date())].append(shape(it.geometry))
cover = defaultdict(list)
for (orbit, day), geoms in days.items():
    cover[orbit].append(unary_union(geoms).intersection(aoi).area / aoi.area)
ORBIT = max(cover, key=lambda o: (sum(c >= 0.99 for c in cover[o]), len(cover[o])))
items = sorted([it for it in found if it.properties["sat:relative_orbit"] == ORBIT], key=lambda i: i.datetime)
others = [f"{o} ({len(c)} dates)" for o, c in cover.items() if o != ORBIT]
print(f"{len({i.datetime.date() for i in items})} Sentinel-1 acquisitions from relative orbit {ORBIT} "
      f"({items[0].datetime:%-d %b} – {items[-1].datetime:%-d %b %Y})" + (f"; other orbits: {', '.join(others)}" if others else ""))

15 Sentinel-1 acquisitions from relative orbit 30 (3 Jun – 30 Nov 2022)

2 — Where can a river flood? Height above the river

Floods spread over the low ground next to rivers. From the Copernicus elevation model we compute, for every pixel, its height above the nearest river (rivers taken from the 1984–2021 JRC water record). Ground more than 8 m above the river is left out: it cannot be flooded by the river, and leaving it out removes most false alarms from dark, dry surfaces such as sand or tarmac.

# @title 3 · Load the radar images, the elevation model and the river network (a few minutes)
ds = load(items, bands=["vv"], bbox=BBOX, crs=CRS, resolution=RES, groupby="solar_day", pool=16)
vv_db = 10 * np.log10(ds["vv"].where(ds["vv"] > 0))
occ = load(list(catalog.search(collections=["jrc-gsw"], bbox=BBOX).items()), bands=["occurrence"], like=ds,
           resampling="nearest")["occurrence"].isel(time=0).values
es = Client.open("https://earth-search.aws.element84.com/v1")
dem = load(list(es.search(collections=["cop-dem-glo-30"], bbox=BBOX).items()), bands=["data"], like=ds,
           resampling="bilinear")["data"].isel(time=0).values.astype(float)
water_ever = ndimage.binary_closing((occ >= 10) & (occ <= 100), iterations=2)
lab, n = ndimage.label(water_ever, structure=np.ones((3, 3)))
sizes = np.bincount(lab.ravel()); sizes[0] = 0
river = np.isin(lab, np.flatnonzero(sizes >= 0.05 * sizes.max()))          # the main river system(s)
dist, (iy, ix) = ndimage.distance_transform_edt(~river, return_indices=True)
hand = dem - dem[iy, ix]                                                    # height above the nearest river
permanent = (occ >= 50) & (occ <= 100)
floodplain = (hand < 8) & (dist * RES / 1000 < 15) & ~permanent
print(f"{vv_db.sizes['time']} radar images loaded · floodplain (≤ 8 m above the river): {floodplain.sum() * (RES / 1000) ** 2:,.0f} km² "
      f"· normal river water: {permanent.sum() * (RES / 1000) ** 2:,.0f} km²")

15 radar images loaded · floodplain (≤ 8 m above the river): 625 km² · normal river water: 101 km²

3 — From backscatter to flood water, date by date

Radar backscatter is expressed in decibels (dB). Open water is darker than about −18 dB in VV polarisation. A pixel counts as flood water on a date when, inside the floodplain and outside the normal river, it is darker than −18 dB and at least 3 dB darker than its own median before the flood. The second condition removes surfaces that are always dark.

# @title 4 · Flooded area on every date
THRESHOLD_DB, DROP_DB = -18, 3
t = pd.to_datetime(vv_db.time.values)
baseline = vv_db.isel(time=np.flatnonzero(t < BASELINE_END)).median("time").values
vv = vv_db.values
flood = (vv < THRESHOLD_DB) & ((vv - baseline) < -DROP_DB) & floodplain
px_km2 = (RES / 1000) ** 2
series = pd.Series(flood.sum(axis=(1, 2)) * px_km2, index=t)
fig = go.Figure(go.Scatter(x=series.index, y=series.values, mode="lines+markers", line=dict(color=COLORS["water"], width=3),
                           name="flooded land"))
fig.update_layout(title="Flooded land seen by Sentinel-1 radar (km², outside the normal river)", yaxis_title="km²",
                  height=420, **LAYOUT)
fig.show()
peak = series.idxmax()
print(f"Normal conditions (before {pd.Timestamp(BASELINE_END):%-d %B}): {series[t < BASELINE_END].median():,.0f} km² · "
      f"peak: {series.max():,.0f} km² on {peak:%-d %B %Y}")

Normal conditions (before 1 August): 1 km² · peak: 169 km² on 13 October 2022

# @title 5 · Why −18 dB? The backscatter of the floodplain before and at the peak
k = int(np.argmax(series.values))
edges = np.arange(-30, 5.5, 0.5)
fig = go.Figure()
for vals, name, color in [(baseline[floodplain], "before the flood (median)", COLORS["muted"]),
                          (vv[k][floodplain], f"peak, {peak:%-d %b %Y}", COLORS["water"])]:
    counts, _ = np.histogram(vals[np.isfinite(vals)], bins=edges)
    fig.add_trace(go.Bar(x=edges[:-1] + 0.25, y=100 * counts / max(counts.sum(), 1), name=name, marker_color=color,
                         opacity=0.65, width=0.5))
fig.add_vline(x=THRESHOLD_DB, line_dash="dash", line_color=COLORS["ink"], annotation_text="−18 dB")
fig.update_layout(barmode="overlay", title="Floodplain backscatter (VV): a dark population appears at the flood peak",
                  xaxis_title="dB", yaxis_title="% of floodplain pixels", height=400, **LAYOUT)
fig.show()

4 — Where did the flood reach?

# @title 6 · Flood map at the peak
fig, ax = plt.subplots(1, 2, figsize=(15, 7))
ax[0].imshow(baseline, cmap="gray", vmin=-25, vmax=0); ax[0].set_title(f"Radar backscatter (VV), normal conditions (median before {pd.Timestamp(BASELINE_END):%-d %b})")
ax[1].imshow(vv[k], cmap="gray", vmin=-25, vmax=0)
ax[1].imshow(np.ma.masked_where(~flood[k], flood[k]), cmap=ListedColormap(["#E4572E"]), alpha=0.8, interpolation="nearest")
ax[1].set_title(f"{peak:%-d %B %Y}: flooded land in red ({series.max():,.0f} km²)")
for a in ax:
    a.axis("off")
plt.tight_layout(); plt.show()

# @title 7 · Interactive flood map
from rasterio.transform import from_bounds
gb = ds.odc.geobox
west, south, east, north = transform_bounds(gb.crs.to_wkt(), "EPSG:4326", *gb.extent.boundingbox)
h, w = flood[k].shape
ll = np.zeros((h, w), np.uint8)
reproject(flood[k].astype(np.uint8), ll, src_transform=gb.transform, src_crs=gb.crs.to_wkt(),
          dst_transform=from_bounds(west, south, east, north, w, h), dst_crs="EPSG:4326", resampling=Resampling.nearest)
rgba = np.zeros((h, w, 4)); rgba[ll == 1] = [0.89, 0.34, 0.18, 0.8]
m = folium.Map(tiles=None, control_scale=True, height=520)
folium.TileLayer(tiles=ESRI, attr=ESRI_ATTR, name="Satellite imagery").add_to(m)
folium.raster_layers.ImageOverlay(rgba, bounds=[[south, west], [north, east]], name=f"Flooded land, {peak:%d %b %Y}",
                                  mercator_project=True).add_to(m)
folium.LayerControl(collapsed=False).add_to(m)
m.fit_bounds([[south, west], [north, east]])
display(m)

Make this Notebook Trusted to load map: File -> Trust Notebook

5 — Reading the result critically

• Missed flooding: water among trees, crops or buildings scatters the radar back strongly and can look like land, so radar flood maps underestimate flooding in forests, tall crops and towns. The map above is a conservative minimum.
• False alarms: very smooth dry surfaces can also look dark; the height-above-river mask and the comparison with each pixel's own normal level remove most of them.
• Timing: each orbit passes every 12 days in 2022, so the true peak can fall between two images.

Your turn

• The White Volta below the Bagre Dam (Burkina Faso → Ghana). Spillage from Bagre floods farmland along the river in Ghana's Upper East Region in wet years. Try BBOX = [-1.0, 10.45, -0.25, 11.15], START, END = "2020-06-01", "2020-11-30", BASELINE_END = "2020-07-15". The flooded band along the river is narrow: try RES = 20.
• The Lower Volta, Ghana, 2023 (BBOX = [0.25, 5.85, 0.70, 6.20], START, END = "2023-06-01", "2024-01-31", BASELINE_END = "2023-09-01"): after the Akosombo releases, most flooding was among vegetation and settlements — see how little of it the radar shows, and why field reports matter.
• A flood you know in your own country: change the area and the dates.

---

Data: Sentinel-1 RTC — contains modified Copernicus Sentinel data (ESA), processed by Microsoft Planetary Computer · Copernicus DEM GLO-30 (© DLR e.V. and © Airbus Defence and Space GmbH, provided under COPERNICUS by the European Union and ESA) · JRC Global Surface Water (Pekel et al. 2016) · Catchment: HydroBASINS (Lehner & Grill 2013) and Natural Earth · Basemap © Esri, Maxar, Earthstar Geographics.

QuestFeed · IC-TWEP 2026 take-home materials · questfeed.com