Burn severity explained: NBR, dNBR and what the classes mean
On this page
Problem statement
After a fire, the first map anyone wants is severity: not just where it burned, but how badly. The standard method is two satellite images and one index — the normalised burn ratio, differenced before and after — and it is simple enough that the difficult parts get skipped. The classes are not physical measurements, the thresholds are conventions from one study area, the index responds to moisture as well as to combustion, and the answer depends on how long after the fire you look.
Over a 555 km² window inside the 2024 Park Fire in northern California, dNBR from Sentinel-2 puts 58% of the window in the "low severity" class and 9% in "moderate-high". Shift the post-fire image from September to the following June and those numbers change, because vegetation has started to return.
Quick answer
NBR is a normalised difference of near-infrared and shortwave infrared; dNBR is pre minus post, scaled by 1000:
import numpy as np
def nbr(nir, swir2):
return (nir - swir2) / (nir + swir2)
def dnbr(pre_nir, pre_swir2, post_nir, post_swir2):
return (nbr(pre_nir, pre_swir2) - nbr(post_nir, post_swir2)) * 1000
USGS_CLASSES = [(-100, 99, "unburned / regrowth"), (99, 269, "low"),
(269, 439, "moderate-low"), (439, 659, "moderate-high"), (659, 1300, "high")]
d = np.load("dnbr.npy") # Park Fire, Sentinel-2 2 July vs 10 September 2024
for lo, hi, name in USGS_CLASSES:
m = (d >= lo) & (d < hi)
print(f"{name:20} {np.nanmean(m):6.1%} {m.sum() * 400 / 1e6:6.1f} km²")unburned / regrowth 12.6% 70.1 km²
low 58.4% 323.9 km²
moderate-low 19.3% 107.3 km²
moderate-high 9.4% 52.4 km²
high 0.2% 1.1 km²
Healthy vegetation is bright in near-infrared and dark in shortwave infrared, so NBR is high. Fire removes the vegetation and dries the surface, which reverses both — a large drop in NBR.
Two dates, matched phenology, the right bands, an offset
1. Choose the two dates deliberately
An initial assessment uses an image within days or weeks of containment and measures the fire's immediate effect. An extended assessment uses the following growing season and measures what survived. They are different questions with the same formula.
2. Match the phenology
Pre- and post-fire images should be at the same point in the season, or the difference includes ordinary green-up and senescence. Where that is impossible, the offset correction in step 4 matters more.
3. Use the right bands
Near-infrared around 850 nm and shortwave infrared around 2200 nm: Sentinel-2 bands 8A and 12, Landsat 8/9 bands 5 and 7.
4. Apply an offset if the unburned area moved
Compute the mean dNBR over unburned ground outside the perimeter and subtract it. It removes seasonal and atmospheric differences between the two dates.
5. Classify with stated thresholds
The widely used class breaks come from the USGS FIREMON programme and were derived in specific ecosystems. Use them, but say so, and check a few pixels against imagery.
6. Consider RdNBR in mixed vegetation
Relativised dNBR divides by the square root of pre-fire NBR, so a shrubland with little to lose is not automatically classed as low severity.
dNBR, RdNBR and the Park Fire
The data used on this page
Three real inputs for the 2024 Park Fire in northern California. The perimeter from the WFIGS interagency fire perimeters service (park_perimeter.geojson). The 2020 Census blocks of California with their POP20 and HOUSING20 counts, from the TIGER/Line block file (382 MB), kept for the bounding box of the main perimeter plus a 0.25° margin (blocks_ca_north.parquet) — only blocks touching the perimeter enter any result, so the margin does not change a number. And the burn severity from How to map burn severity with dNBR: its two Sentinel-2 scenes (2 July and 10 September 2024) on its 20 m window, saved as the uncorrected dNBR (dnbr.npy) and the pre-fire NBR (pre_nbr.npy).
import os
import geopandas as gpd
import numpy as np
import pystac_client
import rasterio
import requests
from rasterio.enums import Resampling
from rasterio.transform import from_origin
from rasterio.vrt import WarpedVRT
from rasterio.warp import transform_bounds
from shapely.geometry import box, shape
WFIGS = ("https://services3.arcgis.com/T4QMspbfLg3qTGWY/arcgis/rest/services/"
"WFIGS_Interagency_Perimeters/FeatureServer/0/query")
if not os.path.exists("park_perimeter.geojson"):
reply = requests.get(WFIGS, params={"where": "poly_IncidentName='PARK' AND attr_POOState='US-CA'",
"outFields": "*", "f": "geojson", "outSR": 4326}, timeout=300)
reply.raise_for_status()
open("park_perimeter.geojson", "w").write(reply.text)
if not os.path.exists("blocks_ca_north.parquet"):
fire = gpd.read_file("park_perimeter.geojson")
x0, y0, x1, y1 = fire[fire.poly_GISAcres > 1000].total_bounds # WFIGS has two small PARK fires too
blocks = gpd.read_file("https://www2.census.gov/geo/tiger/TIGER2020/TABBLOCK20/tl_2020_06_tabblock20.zip",
bbox=(x0 - 0.25, y0 - 0.25, x1 + 0.25, y1 + 0.25))
blocks[["GEOID20", "POP20", "HOUSING20", "ALAND20", "geometry"]].to_parquet("blocks_ca_north.parquet")
if not os.path.exists("dnbr.npy"):
BBOX = [-122.00, 40.00, -121.74, 40.22]
left, bottom, right, top = transform_bounds(4326, "EPSG:32610", *BBOX)
GRID = {"crs": "EPSG:32610", "transform": from_origin(left, top, 20, 20),
"width": int((right - left) // 20), "height": int((top - bottom) // 20)}
es = pystac_client.Client.open("https://earth-search.aws.element84.com/v1")
def scene(dates):
items = [i for i in es.search(collections=["sentinel-2-c1-l2a"], bbox=BBOX, datetime=dates,
query={"eo:cloud_cover": {"lt": 5}}).item_collection()
if shape(i.geometry).contains(box(*BBOX))]
return sorted(items, key=lambda i: i.datetime)[0]
def band(item, name):
with rasterio.open(item.assets[name].href) as src, WarpedVRT(
src, crs=GRID["crs"], transform=GRID["transform"], width=GRID["width"],
height=GRID["height"], resampling=Resampling.bilinear) as vrt:
a = vrt.read(1).astype("float32")
return np.where(a > 0, a * 0.0001, np.nan)
def nbr(item):
return (band(item, "nir08") - band(item, "swir22")) / (band(item, "nir08") + band(item, "swir22"))
pre, post = scene("2024-07-01/2024-07-10"), scene("2024-09-08/2024-09-15")
np.save("pre_nbr.npy", nbr(pre))
np.save("dnbr.npy", (nbr(pre) - nbr(post)) * 1000)
print(len(gpd.read_parquet("blocks_ca_north.parquet")), "blocks;", np.load("dnbr.npy").shape)8735 blocks; (1235, 1123)
Example 1 — the offset, measured on unburned ground
import geopandas as gpd, rasterio
from rasterio.features import geometry_mask
from rasterio.transform import from_origin
from rasterio.warp import transform_bounds
bbox = [-122.00, 40.00, -121.74, 40.22]
left, bottom, right, top = transform_bounds(4326, "EPSG:32610", *bbox)
transform = from_origin(left, top, 20, 20)
fire = gpd.read_file("park_perimeter.geojson")
fire = fire[fire.poly_GISAcres > 1000].to_crs("EPSG:32610")
inside = geometry_mask(fire.geometry, out_shape=d.shape, transform=transform, invert=True)
offset = np.nanmedian(d[~inside])
print(f"pixels inside the perimeter: {inside.mean():.0%}")
print(f"median dNBR inside {np.nanmedian(d[inside]):6.0f}, outside {offset:6.0f}")
corrected = d - offset
for lo, hi, name in USGS_CLASSES[1:]:
before = np.nanmean((d >= lo) & (d < hi) & inside)
after = np.nanmean((corrected >= lo) & (corrected < hi) & inside)
print(f" {name:14} {before:6.1%} -> {after:6.1%} after the offset")pixels inside the perimeter: 94%
median dNBR inside 205, outside -9
low 58.4% -> 58.0% after the offset
moderate-low 19.3% -> 20.3% after the offset
moderate-high 9.4% -> 10.1% after the offset
high 0.2% -> 0.3% after the offset
Example 2 — dNBR against RdNBR
pre_nbr = np.load("pre_nbr.npy")
rdnbr = d / np.sqrt(np.abs(pre_nbr).clip(0.001))
for name, arr, breaks in [("dNBR", d, (269, 439)), ("RdNBR", rdnbr, (500, 750))]:
lo, hi = breaks
print(f"{name:6} moderate-low {np.nanmean((arr >= lo) & (arr < hi) & inside):6.1%} "
f"moderate-high or worse {np.nanmean((arr >= hi) & inside):6.1%}")dNBR moderate-low 19.3% moderate-high or worse 9.6%
RdNBR moderate-low 34.9% moderate-high or worse 45.3%
RdNBR is not a rescaled dNBR with the same classes: it has its own thresholds, and on this window it assigns far more area to the higher classes. Pick one index and its own breaks, and say which.
Example 3 — severity is a distribution, not five buckets
q = np.nanpercentile(d[inside], [5, 25, 50, 75, 95])
print("dNBR percentiles inside the perimeter:", np.round(q).astype(int).tolist())
print(f"pixels within 20 of a class break: "
f"{np.nanmean(np.min(np.abs(d[inside, None] - np.array([99, 269, 439, 659])), axis=1) < 20):.1%}")dNBR percentiles inside the perimeter: [83, 154, 205, 303, 510]
pixels within 20 of a class break: 19.1%
Nearly a fifth of the burned pixels sit within 20 dNBR units of a class boundary, where the class is decided as much by speckle, shadow and co-registration as by fire.
Why NBR uses shortwave infrared
Why NBR uses shortwave infrared
Live vegetation contains water, which absorbs strongly around 2200 nm, so healthy canopy is dark there and bright in the near-infrared. Burning removes the leaf water and exposes char and soil, which raises shortwave reflectance and lowers near-infrared. Differencing the two bands captures both halves of that change at once.
Why the classes are conventions
The commonly quoted breaks — 99, 269, 439, 659 — come from field plots in specific forest types, matched to a composite burn index measured on the ground. Applied to chaparral, grassland or boreal forest they do not mean the same thing. They are useful because everyone uses them, not because they are physical constants.
Why timing changes the answer
An image taken days after the fire sees ash and scorched canopy. A year later it sees what regrew, so dNBR falls and the map moves towards lower severity. Neither is wrong; they answer "how hard was it hit" and "what survived". The date pair belongs in the title of the map.
Cloud, smoke, water and snow
- Cloud and smoke must be masked; smoke lowers NBR and looks like burning.
- Water has extreme NBR values and should be masked before classification.
- Snow after the fire breaks the comparison entirely.
- Sentinel-2 band 8A (narrow NIR) is preferred over band 8 for NBR because it matches Landsat's band better.
- Field validation with the composite burn index is what the thresholds were calibrated against.
- Fire perimeters and severity are different products; a perimeter says where, severity says how much.
Internal links
- How to map burn severity with dNBR in Python — the full workflow
- dNBR shows change outside the fire perimeter — the fix guide
- How to work with fire perimeter data in Python — the perimeter half
- Hazard, exposure and vulnerability explained — where severity fits into risk
- Vegetation indices explained — the same arithmetic for greenness
- How to detect change between satellite images in Python — change detection generally
FAQ
What is dNBR?
The difference between pre-fire and post-fire normalised burn ratio, scaled by 1000, used as a burn severity index.
Which bands does NBR use?
Near-infrared and shortwave infrared — Sentinel-2 bands 8A and 12, Landsat bands 5 and 7.
What do the severity classes mean?
They are conventional breaks from USGS FIREMON work, calibrated to field plots in particular ecosystems.
Should I use dNBR or RdNBR?
RdNBR where pre-fire vegetation varies a lot, because it scales by what there was to burn.
When should the post-fire image be taken?
Soon after containment for an initial assessment, or the next growing season for an extended one. Say which.
Why is there severity outside the perimeter?
Usually phenology, moisture or smoke rather than fire — apply an offset and check.