Is my reservoir / lake / river running low compared to historical normals?

Water-surface elevation is the height of a lake or river’s surface above sea level, like reading a giant ruler standing in the water. When a reservoir drops, that number falls. Satellites can now measure this height from orbit, and they also watch the water’s surface area (how much ground it covers). So you can tell whether your lake is running low compared to its normal years without ever visiting it.

Which data to use, and how

Use this dataset: SWOT lake heights (short name SWOT_L2_HR_RiverSP_reach_2.0), a satellite that measures the height of rivers (≥100 m wide) and lakes (≥250×250 m) every 21 days. It’s the most direct answer to “how high is the water.” There’s an easier route, too: get those SWOT heights through a free helper service called Hydrocron, which hands you a ready-made time series for a named lake instead of raw satellite pixels.

Then do four steps (the runnable code further down does all of this; this is just the idea):

1. Pick your water body. Draw a small latitude/longitude box around it and choose the years you care about.
2. Get the height series. Ask SWOT (via Hydrocron) for the water-surface elevation on each pass. That gives you a height-over-time chart.
3. Cross-check the area. Use HLS optical images to measure how much ground the water covers (a “surface-area” series). If both height and area are falling, you can trust the signal.
4. Compare to normal. Line up this year against the multi-year average for the same months. Below normal means drought conditions; near normal is fine.

That’s the whole method. Below, “How a scientist answers this” names the exact tools, and the Run it cell does it live on synthetic data so you can see each step work.

What you can find out

• Whether the water level is rising or falling over time. That’s a height series from SWOT for bigger bodies, or ICESat-2 (a laser-altimetry satellite) for smaller ones.
• Whether the water is shrinking in area, using HLS optical images on cloud-free days.
• The bigger regional water picture, from GRACE-FO, a pair of satellites that weigh changes in all the water in a 300 km region (surface + soil + groundwater combined). That gives a “total water storage” trend.
• How much rain is feeding in, from GPM IMERG precipitation, so you can see if a drop follows a dry spell.
• A drought flag: this year compared against the multi-year normal for each of those.

What it can’t tell you

• Underground water on its own. Groundwater is hidden below the surface. Isolating it needs GRACE-FO plus soil-moisture data plus a land-surface model. The satellites can’t separate it by themselves.
• Fast, same-day changes. During a storm surge or a sudden dam release, the water moves faster than the satellite revisits (every few days to weeks). For hour-by-hour levels you need ground gauges or hydraulic models.
• Water quality. Sediment, algae, and salinity don’t show up in a height or area measurement. That needs separate optical or hyperspectral data layers.

Gotchas to watch for

• SWOT can miss your lake on a pass. It revisits every 21 days and has a ~20 km blind strip directly beneath it (the “nadir gap”), so some passes simply don’t cover your water body. Gaps in your height series are normal rather than errors, so collect enough passes and don’t panic over a missing date.
• Clouds wreck the area estimate. The surface-area method (NDWI on optical images) can’t see through cloud, and a cloudy scene gives a garbage area number. Keep only low-cloud scenes; the code filters with cloud_cover=20.
• GRACE-FO is coarse. At ~300 km resolution it can’t isolate a single reservoir, seeing the region’s water as one blurry blob instead. Treat it as context for the big-picture trend, not as a direct measurement of your specific lake.
• Mixing image sensors isn’t perfect. HLS blends data from different Landsat satellites (the Landsat 7 era vs the 8/9 era), and they’re calibrated well but not identically. A small step in the area record can be the sensor rather than real change, so be cautious comparing across the sensor handover.
• Use Hydrocron, not raw pixels, to start. The raw SWOT product (the “L2 HR pixel-cloud”) is a firehose that’s slow and fiddly for a first look. Let Hydrocron pre-aggregate the heights to your lake for you.

The real-data code

The Run it cell above runs the method on synthetic data with no login. Below is the whole recipe against the real archive: pull the SWOT height series, cross-check it with an HLS surface-area series, compare this year to the multi-year normal, and plot both. It needs a free Earthdata Login.

import earthaccess, requests
import numpy as np, rasterio
import matplotlib.pyplot as plt
from datetime import datetime

earthaccess.login(strategy="netrc")

aoi = (-114.7, 36.0, -114.0, 36.4)        # Lake Mead, NV/AZ as (W, S, E, N)
window = ("2024-01-01", "2025-12-31")

# 1. SWOT water-surface elevation via Hydrocron (pre-aggregated; skip the raw pixel-cloud).
lake_id = "7250000000"   # Hydrocron PriorLake id for Lake Mead (look up in their catalog)
url = ("https://soto.podaac.earthdatacloud.nasa.gov/hydrocron/v1/timeseries"
       f"?feature=PriorLake&feature_id={lake_id}"
       f"&start_time={window[0]}T00:00:00Z&end_time={window[1]}T23:59:59Z"
       "&output=csv&fields=time_str,wse,area_total")
rows = [ln.split(",") for ln in requests.get(url).json()["results"]["csv"].splitlines()[1:]]
swot_t = np.array([datetime.fromisoformat(r[0].replace("Z", "+00:00")) for r in rows])
wse = np.array([float(r[1]) for r in rows])         # metres above geoid
wse = wse[(wse > -1e8)]; swot_t = swot_t[: len(wse)]  # drop Hydrocron fill values

# 2. Surface-area cross-check from HLS: open each low-cloud scene, build an NDWI water mask,
#    count water pixels and convert to km^2 (HLS is 30 m).
hls = earthaccess.search_data(short_name="HLSL30", bounding_box=aoi,
                              temporal=window, cloud_cover=20)
hls_t, area = [], []
for g in hls:
    urls = {u.split(".")[-2]: u for u in g.data_links()}   # ...B03.tif, ...B05.tif
    with rasterio.open(urls["B03"]) as gx, rasterio.open(urls["B05"]) as nx:
        green = gx.read(1).astype("f4"); nir = nx.read(1).astype("f4")
    ndwi = (green - nir) / (green + nir + 1e-6)
    water = ndwi > 0.0
    area.append(water.sum() * (30 ** 2) / 1e6)             # km^2
    hls_t.append(datetime.fromisoformat(g["umm"]["TemporalExtent"]
                 ["RangeDateTime"]["BeginningDateTime"].replace("Z", "+00:00")))
hls_t, area = np.array(hls_t), np.array(area)

# 3. Compare to normal: is the latest height below the period's typical level?
normal, latest = np.nanmedian(wse), wse[-1]
drop = latest - normal
verdict = "BELOW normal (running low)" if drop < -0.5 else "near normal"
print(f"Lake Mead WSE: latest {latest:.2f} m vs normal {normal:.2f} m "
      f"({drop:+.2f} m) -> {verdict}; "
      f"area {area[-1]:.0f} km2 (was {area[0]:.0f})" if len(area) else "")

# 4. Plot height and area on a shared time axis so the two signals corroborate.
fig, ax = plt.subplots(2, 1, sharex=True, figsize=(9, 6))
ax[0].plot(swot_t, wse, "o-"); ax[0].axhline(normal, ls="--", c="k", label="normal")
ax[0].set_ylabel("WSE (m)"); ax[0].legend()
ax[1].plot(hls_t, area, "s-", c="tab:green"); ax[1].set_ylabel("Area (km2)")
plt.tight_layout(); plt.show()

# Trust it only if height and area fall together; a coarse GRACE-FO (GRACEFO_L3_JPL_RL06.X_M)
# regional storage trend is the third, independent cross-check (see the gotchas on its 300 km blur).

Where the data comes from

Everything comes from NASA through one free login (Earthdata Login), even though it’s gathered by different teams. The water heights (SWOT) and the regional water-storage data (GRACE-FO) come from the PO.DAAC archive; the optical images (HLS) come from LP DAAC; and the rainfall (GPM IMERG) comes from GES DISC. Three different archives, one set of credentials.

Sources

• SWOT Hydrocron API: https://podaac.jpl.nasa.gov/hydrocron
• SWOT-galleries community: https://swot-community.github.io/SWOT-galleries/
• GRACE Tellus: https://grace.jpl.nasa.gov/