Why Is My Burned Area From GOES Fire Pixels Wrong?

On this page

Problem statement

GOES fire detections are tempting for burned area: they arrive every few minutes, cover the whole continent and are free. Turning them into square kilometres goes wrong in four common ways. Multiplying the number of detections by a pixel area counts the same pixel every scan. Using the nominal "2 km" pixel ignores how much larger pixels are away from the point under the satellite. Treating a pixel as burned ignores that a small hot fire anywhere in it is enough. And the product's own Area variable sounds like burned area but is something else.

This page puts numbers on each mistake for the first two and a half days of the 2024 Park Fire, using GOES-18 scans every 15 minutes and the fire's final CAL FIRE perimeter.

Quick answer

Count each pixel once, give it its real footprint, and treat the union as an outline of where the fire was, not as burned area:

foot = gpd.GeoSeries(list(CELL.values()), crs=3310)
union = foot.union_all()
print(f"{len(det):,} detections in {len(CELL)} distinct pixels; footprint union {union.area / 1e6:,.0f} km², "
      f"of which {union.difference(PARK).area / 1e6:,.0f} km² lies outside the final perimeter ({PARK.area / 1e6:,.0f} km²)")

6,535 detections in 201 distinct pixels; footprint union 1,361 km², of which 153 km² lies outside the final perimeter (1,739 km²)

Map of the Park Fire's final perimeter covered by 201 GOES-18 fire pixel footprints from 24 to 26 July 2024: 159 pixels lie inside the perimeter, 34 cross its edge, mostly along the southern half, and 8 lie wholly outside it; the parts of the footprints outside the perimeter, 153 square kilometres, are hatched. — Every fire pixel counted once, with its real footprint. Along the edge, footprints reach well past the burned area.

Never pixel counts times area, use real footprints

1. Do not multiply detections by pixel area

A burning pixel is flagged in every scan until the fire there dies down. Summing detections counts it each time:

med = foot.area.median() / 1e6
print(f"{len(det):,} detections x {med:.2f} km² = {len(det) * med:,.0f} km², {len(det) * med / (PARK.area / 1e6):.0f} times the final fire")
print(f"scans per pixel: median {det.groupby(['j', 'i']).size().median():.0f}, max {det.groupby(['j', 'i']).size().max()}")

6,535 detections x 6.77 km² = 44,252 km², 25 times the final fire
scans per pixel: median 30, max 103

2. Use the real pixel size

"2 km" is the pixel spacing at the sub-satellite point, on the equator at 137.2° W for GOES-18. Northern California is seen at an angle, and the footprint grows:

print(f"footprint of the {len(foot)} fire pixels: {foot.area.min() / 1e6:.2f} to {foot.area.max() / 1e6:.2f} km² (nominal 2 km x 2 km = 4 km²)")
print(f"distinct pixels x 4 km² = {len(CELL) * 4:,} km²; with real footprints = {foot.area.sum() / 1e6:,.0f} km²")
c = np.array(foot.iloc[0].exterior.coords)
print(f"one footprint: {np.hypot(*(c[1] - c[0])) / 1000:.2f} km east-west by {np.hypot(*(c[2] - c[1])) / 1000:.2f} km north-south")

footprint of the 201 fire pixels: 6.71 to 6.81 km² (nominal 2 km x 2 km = 4 km²)
distinct pixels x 4 km² = 804 km²; with real footprints = 1,361 km²
one footprint: 2.23 km east-west by 3.13 km north-south

Bar chart on a logarithmic scale of six areas from the same pixels: detections times pixel area 44,252 square kilometres; the Area variable summed over scans 665; distinct pixels times 4 square kilometres 804; distinct pixels with real footprints 1,361, of which 1,207 inside the perimeter; and the final CAL FIRE perimeter 1,739. — The same pixels give answers that differ by a factor of more than 60, depending on how they are counted.

3. Expect footprints to spill over the edge

A pixel is flagged when some part of it is hot enough, so pixels along the fire's edge include unburned land. Compare the footprints with the final perimeter:

edge = foot[foot.intersects(PARK) & ~foot.within(PARK)]
print(f"{len(edge)} of {len(foot)} pixels cross the final perimeter; on average {1 - edge.intersection(PARK).area.sum() / edge.area.sum():.0%} of their area lies outside it")
print(f"{(~foot.intersects(PARK)).sum()} pixels lie wholly outside it: flagged, though no part of them burned in this fire")

34 of 201 pixels cross the final perimeter; on average 43% of their area lies outside it
8 pixels lie wholly outside it: flagged, though no part of them burned in this fire

The same effect works the other way for small fires. A fire of a few hectares fills one pixel of almost 7 km², so the footprint can be a hundred times the burned area.

4. Do not read the Area variable as burned area

The FDC file's Area variable is the estimated instantaneous area of the flaming fire within a pixel, from the same sub-pixel fit as the fire temperature. It exists only for some fire categories:

print(f"Area present for {det.area_m2.notna().mean():.0%} of detections; median {det.area_m2.median() / 1e4:.1f} ha per pixel, "
      f"max {det.area_m2.max() / 1e6:.2f} km²")
print("detections with an Area value, by mask code:", det[det.area_m2.notna()].code.value_counts().to_dict())
print(f"summed over all scans: {det.area_m2.sum() / 1e6:,.0f} km², a number with no physical meaning")

Area present for 78% of detections; median 9.5 ha per pixel, max 1.18 km²
detections with an Area value, by mask code: {30: 5115, 10: 1}
summed over all scans: 665 km², a number with no physical meaning

The flaming area in each scan

The data used on this page

GOES-18 ABI L2 Fire - Hot Spot Characterization (ABI-L2-FDCC, CONUS sector) from NOAA's archive on AWS, one scan per 15 minutes from 21:00 UTC on 24 July to 23:45 UTC on 26 July 2024, and the Park Fire's final perimeter from CAL FIRE's last-five-years perimeter service. Pixel footprints are built from the corners of each pixel in the satellite's fixed-grid projection, then projected to California Albers (EPSG:3310).

import re
from pathlib import Path

import geopandas as gpd
import numpy as np
import pandas as pd
import requests
import xarray as xr
from pyproj import Proj
from shapely.geometry import Polygon

Path("data/fdc").mkdir(parents=True, exist_ok=True)
B = "https://noaa-goes18.s3.amazonaws.com/"
keys = []
for day in (206, 207, 208):
    for h in range(21 if day == 206 else 0, 24):
        text = requests.get(B, params={"list-type": 2, "prefix": f"ABI-L2-FDCC/2024/{day}/{h:02d}/"}, timeout=60).text
        keys += [k for k in re.findall(r"<Key>([^<]+)</Key>", text) if int(re.search(r"_s\d{9}(\d{2})", k).group(1)) % 15 < 5]
FILES = []
for k in keys:
    f = Path("data/fdc", k.rsplit("/", 1)[1])
    if not f.exists():
        f.write_bytes(requests.get(B + k, timeout=120).content)
    FILES.append(f)
with xr.open_dataset(FILES[0], engine="h5netcdf") as ds:
    p = ds.goes_imager_projection
    meanings = dict(zip(ds.Mask.attrs["flag_values"].tolist(), ds.Mask.attrs["flag_meanings"].split()))
    proj = Proj(proj="geos", h=p.perspective_point_height, lon_0=p.longitude_of_projection_origin, sweep=p.sweep_angle_axis,
                a=p.semi_major_axis, b=p.semi_minor_axis)
    x, y = ds.x.values * p.perspective_point_height, ds.y.values * p.perspective_point_height
fire_codes = [c for c, n in meanings.items() if n.endswith("fire_pixel") and "no_fire" not in n]
bx, by = proj([-122.3, -121.3, -122.3, -121.3], [39.6, 39.6, 40.5, 40.5])
IX, IY = (x >= min(bx)) & (x <= max(bx)), (y >= min(by)) & (y <= max(by))
rows = []
for f in FILES:
    with xr.open_dataset(f, engine="h5netcdf") as ds:
        m = ds.Mask.values[np.ix_(IY, IX)]
        a = ds.Area.values[np.ix_(IY, IX)]
        t = pd.Timestamp(ds.attrs["time_coverage_start"].rstrip("Z"), tz="UTC")
        for j, i in zip(*np.nonzero(np.isin(m, fire_codes))):
            rows.append((t, j, i, int(m[j, i]), float(a[j, i])))
det = pd.DataFrame(rows, columns=["time", "j", "i", "code", "area_m2"])
perims = gpd.read_file(requests.get("https://services1.arcgis.com/jUJYIo9tSA7EHvfZ/arcgis/rest/services/California_Fire_Perimeters_last_5_years_View/FeatureServer/0/query",
                                    params={"where": "FIRE_NAME='PARK' AND YEAR_=2024", "outFields": "FIRE_NAME", "outSR": 4326, "f": "geojson"}, timeout=120).text)
PARK = perims.to_crs(3310).geometry.make_valid().union_all()
xs, ys = x[IX], y[IY]
dx, dy = abs(xs[1] - xs[0]), abs(ys[1] - ys[0])
CELL = {}
for j, i in set(zip(det.j, det.i)):                                              # footprint of each fire pixel from its corners
    cx = [xs[i] - dx / 2, xs[i] + dx / 2, xs[i] + dx / 2, xs[i] - dx / 2]
    cy = [ys[j] - dy / 2, ys[j] - dy / 2, ys[j] + dy / 2, ys[j] + dy / 2]
    lon, lat = proj(cx, cy, inverse=True)
    CELL[(j, i)] = gpd.GeoSeries([Polygon(zip(lon, lat))], crs=4326).to_crs(3310).iloc[0]
print(f"{len(FILES)} scans; pixel spacing {dx:,.0f} m in the satellite's scan plane; {len(det):,} fire detections")

204 scans; pixel spacing 2,004 m in the satellite's scan plane; 6,535 fire detections

Example 1: the flaming area in each scan

Summed within one scan, Area estimates how much ground was actively flaming at that moment. That is a measure of intensity, closer to fire radiative power than to burned area:

per_scan = det.groupby("time").area_m2.sum() / 1e6
peak = per_scan.idxmax()
print(f"largest flaming area in one scan: {per_scan.max():.1f} km² at {peak.tz_convert('America/Los_Angeles'):%d %b %H:%M} Pacific; "
      f"median over scans with fire {per_scan[per_scan > 0].median():.1f} km²")

largest flaming area in one scan: 10.8 km² at 25 Jul 15:16 Pacific; median over scans with fire 3.0 km²

A fire pixel is a yes-or-no answer

A fire pixel is a yes-or-no answer for a cell several square kilometres in size: did something in this cell burn hot enough, in this scan, to stand out against the background? Every way of turning those answers into burned area has a bias. Counting detections multiplies by time. Counting distinct pixels is better, but the pixel size depends on where the fire is in the satellite's view, and pixels along the edge are partly unburned: here 43% of their area lies outside the final perimeter.

Which product answers which question.

GOES detections answer when and roughly where a fire is burning, and how hard (through fire radiative power). Burned area comes from perimeters, mapped from aircraft and on the ground, or from before-and-after imagery at 10 to 30 m, such as Sentinel-2 or Landsat burn severity. Polar-orbiting VIIRS detections at 375 m sit in between: fine enough to trace the fire front, too coarse and too infrequent for an area.

East versus West, parallax and fill values

GOES burned area problems and fixes.

GOES-East vs GOES-West: the same fire seen from 75.2° W has different, larger pixels than from 137.2° W; never mix their footprints.
Parallax: pixel positions assume a surface at ellipsoid height; on high terrain they shift by up to a few kilometres.
Missing scans shorten the time a pixel stays flagged; they do not change the distinct-pixel count much.
Full-disk vs CONUS products share the pixel grid; footprints are the same.
Area fill values: pixels without a sub-pixel fit (saturated, cloudy, low-probability) have no Area; do not fill them with zeros.

Internal links

How to track a wildfire with GOES fire detections — reading the product
How to map a wildfire's spread from GOES detections — the footprints used here, through time
How to map burn severity with dNBR — burned area from imagery

FAQ

How big is a GOES fire pixel?

2 km at the sub-satellite point. Over northern California from GOES-18 the footprint is about 2.2 by 3.1 km, roughly 6.8 km².

Can I estimate burned area from GOES fire detections?

Only roughly. Count each pixel once and use its real footprint, and you get an outline of where the fire was that includes unburned land along every edge. Use perimeters or burn severity imagery for area.

What is the Area variable in the GOES fire product?

The estimated area of the actively burning part of a pixel at the time of the scan, from a sub-pixel fit. It is not burned area.

Why does my GOES burned area exceed the fire size many times over?

You are probably summing detections over time. A pixel flagged in a hundred scans contributes its area a hundred times.