How severe was this wildfire’s burn?

Burn severity is a map of how badly each patch of land was hit by a fire. It ranges from a light scorch where the grass browned to total destruction where every tree died. You can’t always tell from the ground, but you can tell from space. Charred, dried-out vegetation reflects infrared light very differently from healthy green plants. Compare a picture of the land before the fire with one after, and a satellite can grade the damage patch by patch.

Which data to use, and how

Use this dataset: HLS (Harmonized Landsat Sentinel-2) (short name HLSS30). HLS is NASA’s blended Landsat + Sentinel-2 product. It gives you sharp 30-metre pictures of the land surface in many colours of light, including the infrared bands you need to spot burned ground. It’s the easiest single source here because it’s already cleaned up and ready to use.

Then do four steps. The runnable code further down does all of this; this is just the idea.

1. Grab two clear pictures, one before the fire and one after, of the same burned area. “Clear” means no clouds and (for the after image) no fresh smoke hiding the ground.
2. Turn each picture into a burn index. For each image, compute the Normalized Burn Ratio (NBR). It’s a simple formula that contrasts near-infrared light (healthy plants reflect lots of it) with shortwave-infrared light (char and dry ground reflect lots of it): NBR = (NIR − SWIR) / (NIR + SWIR). Burning lowers NBR.
3. Subtract before from after. The change, dNBR = NBR_before − NBR_after, is the severity signal. A bigger drop means a worse burn.
4. Sort into severity classes and clip to the fire. Bin dNBR into unburned / low / moderate / high, then trim it to the fire’s outline so you’re only mapping the area that actually burned.

That’s the whole method. Below, “How you’d approach it” 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

• A burn-severity map showing where it burned and how hard, patch by patch, using the Normalized Burn Ratio described above.
• Severity, not just extent. By differencing before-and-after NBR (dNBR), you separate a light scorch from a stand-replacing burn and bin the land into unburned / low / moderate / high classes.
• Knock-on effects. Pair in MODIS/VIIRS LST and NDVI differences to see the change in land surface heat (LST = land surface temperature) and greenness (NDVI = a vegetation index) the fire left behind.

What it can’t tell you

• What caused the fire, or how it spread. This is an after-the-fact damage map, made from before/after snapshots. It’s not a record of the fire’s behaviour while it was burning. For that you’d need active-fire detections taken during the event.
• Anything hidden under fresh smoke or cloud. Smoke and cloud sit between the satellite and the ground and throw off the infrared bands. The fix: use the first clear post-fire scene you can find, even if it’s a few days later.
• Whether the damage is to the soil or the canopy, precisely. dNBR is a proxy. It correlates with severity but doesn’t measure it directly. To turn it into ground-truth severity, scientists calibrate it against field plots (a survey method called CBI, Composite Burn Index).

Gotchas to watch for

• Pick clear scenes, not just convenient ones. Clouds and lingering smoke corrupt the infrared bands that NBR depends on, and a smoky “after” image makes the burn look worse than it was. Reach for the first genuinely clear post-fire image, even if it means waiting a little.
• Before and after must line up. Compare a green spring “before” to a dry summer “after” and you inflate the apparent burn. Use two images that cover the same ground, in the same season if possible, so you’re measuring fire damage and not a normal seasonal change in greenness.
• dNBR is a proxy, not a verdict. The unburned/low/moderate/high thresholds are conventions, not laws of nature, and they vary by landscape, so the same dNBR number can mean different things in a forest vs grassland. The honest answer is to calibrate against field plots (CBI); short of that, treat the classes as relative, not absolute.

The real-data code

The Run it cell above runs the method on synthetic data (no login). Below is the whole recipe against the real archive: search a clear pre/post HLS pair, open the infrared bands, compute NBR, difference it, classify the severity, and map it. It needs a free Earthdata Login.

import earthaccess, numpy as np, rasterio
from rasterio.warp import transform_bounds
import matplotlib.pyplot as plt

earthaccess.login(strategy="netrc")

aoi = (-120.6, 38.5, -120.0, 38.9)   # burned area as (W, S, E, N) — here, the Caldor Fire, CA

# 1. Find one clear PRE-fire scene and one clear POST-fire scene over the box.
#    HLSS30 carries NIR in band B08 and SWIR in band B11 (Sentinel-2 naming).
def clearest(short_name, t0, t1):
    grans = earthaccess.search_data(short_name=short_name, bounding_box=aoi,
                                    temporal=(t0, t1), cloud_cover=(0, 20))
    grans.sort(key=lambda g: g["umm"].get("CloudCover", 100))   # least-cloudy first
    return grans[0]

pre  = clearest("HLSS30", "2021-07-01", "2021-08-13")   # before the fire
post = clearest("HLSS30", "2021-10-21", "2021-11-30")   # after the fire

# 2. Open the NIR + SWIR bands of one scene and compute NBR = (NIR - SWIR)/(NIR + SWIR),
#    cropped to the AOI. HLS COGs are one GeoTIFF per band; pick by the "_B08"/"_B11" suffix.
def nbr(granule):
    links = earthaccess.results.DataGranule(granule).data_links()
    pick  = lambda tag: next(u for u in links if u.endswith(f"_{tag}.tif"))
    bands = {}
    for tag in ("B08", "B11"):
        with rasterio.open(earthaccess.open([pick(tag)])[0]) as src:
            w, s, e, n = transform_bounds("EPSG:4326", src.crs, *aoi)
            win = src.window(w, s, e, n)
            arr = src.read(1, window=win).astype("float32") * 1e-4   # scale to reflectance
            arr[arr <= 0] = np.nan
            bands[tag] = arr
    nir, swir = bands["B08"], bands["B11"]
    return (nir - swir) / (nir + swir)

dnbr = nbr(pre) - nbr(post)   # 3. dNBR = NBR_pre - NBR_post; bigger drop = worse burn

# 4. Classify into the standard severity bins (USGS/Key & Benson dNBR thresholds).
edges  = [-np.inf, 0.10, 0.27, 0.66, np.inf]   # unburned | low | moderate | high
labels = ["unburned", "low", "moderate", "high"]
cls    = np.digitize(dnbr, edges[1:-1])
px_km2 = (30 * 30) / 1e6                         # HLS pixel ≈ 30 m
burned = np.isin(cls, [1, 2, 3]) & ~np.isnan(dnbr)
high   = (cls == 3) & ~np.isnan(dnbr)
print(f"Burned area: {burned.sum()*px_km2:.1f} km2; "
      f"high-severity: {high.sum()*px_km2:.1f} km2 "
      f"({100*high.sum()/max(burned.sum(),1):.0f}% of the burn)")

# Map the severity classes over the AOI.
plt.imshow(np.where(np.isnan(dnbr), np.nan, cls), cmap="YlOrRd", vmin=0, vmax=3)
cb = plt.colorbar(ticks=[0, 1, 2, 3]); cb.ax.set_yticklabels(labels)
plt.title("Caldor Fire burn severity (dNBR)"); plt.tight_layout(); plt.show()

# Before you trust it: the dNBR bins are conventions, not ground truth — calibrate against
# field CBI plots, and make sure the post-fire scene is genuinely smoke-free (see gotchas).

Where the data comes from

Everything here comes from NASA through one free login (Earthdata Login). HLS, the blended Landsat (USGS) and Sentinel-2 (ESA) product, is distributed by NASA’s LP DAAC, the land-data archive. The optional knock-on layers (MODIS/VIIRS land surface temperature and NDVI) come from the same NASA family of land products. This question backs the Recover phase of the NASA Disasters program, whose Disasters Portal publishes burn-severity as a headline fire product.

Sources

• Interactive NDVI guide — NBR is its close sibling.
• NASA Disasters program (Recover phase): https://disasters.openveda.cloud/