How much carbon is locked in this forest, and where is it densest?

A forest stores most of its carbon in wood: trunks and big branches. So the real question is “how much woody mass is standing here?” To measure that from space, NASA’s GEDI instrument fired laser pulses down from the International Space Station, timed how long they took to bounce back, and used that to map the height and structure of the tree canopy below. Taller, denser canopy means more wood, and roughly half of that wood is carbon. A map of woody mass is effectively a carbon map.

Which data to use, and how

Use this dataset: GEDI L4A biomass (short name GEDI_L4A_AGB_Density_V2_1_2056). It’s NASA’s ready-made answer: each laser shot already comes with an above-ground biomass density number in Mg/ha (megagrams, i.e. tonnes, of dry wood per hectare). You don’t have to convert raw laser data yourself; the hard physics is already done. The catch is that GEDI only measures along its orbit tracks: narrow ~25 m footprints scattered across the forest, not a continuous wall-to-wall image.

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

1. Pick your forest. A small latitude/longitude box around it, plus a multi-year window (biomass barely changes month to month, so you pool several years to gather enough shots).
2. Keep only the good shots. Each footprint has a quality flag (l4_quality_flag); keep the ones marked good (== 1) with a real positive biomass value, and drop the rest.
3. Find the typical value. Take the median biomass of all your good shots. That’s the representative density of the forest, in Mg/ha. Multiply by ~0.47 to estimate carbon.
4. Find the densest spots. Look at the highest shots (e.g. the top 10%) and their locations. Those mark the carbon-richest, tallest, least-disturbed patches of forest.

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.

Verified locally. For a 1°×1° block of Rondônia in the Brazilian Amazon (9.2–10.2 °S), the good GEDI shots had a median of 142.4 Mg/ha across 10,054 quality footprints. That’s dense tropical forest, with the top 10% of shots above ~317 Mg/ha where the canopy is tallest. Because GEDI samples along tracks, the honest read is “the typical and the densest spots we sampled,” not “every hectare.”

What you can find out

• The typical biomass of a forest. The median of the good shots over your area, in Mg/ha (multiply by ~0.47 for an above-ground carbon estimate).
• Where the densest, tallest forest is. The high-biomass footprints cluster in intact, undisturbed stands; map their locations to find the carbon-richest patches.
• How biomass varies across a landscape. Compare riverside, upland, and forest-edge zones by pooling the shots that fall in each.
• A disturbance story. Low-biomass shots near roads and clearings versus high-biomass interior forest sketch a degradation / regrowth gradient.
• A picture behind each shot. Pair the footprints with HLSL30 30 m satellite imagery to check whether a low reading is a clearing, a road, or genuinely short forest.

What it can’t tell you

• A total carbon figure for the whole area. GEDI samples along tracks, so you have scattered footprints, not a gap-free map. A grand total needs a gridded product or a model that fills in the gaps between shots.
• Below-ground or soil carbon. This biomass number is above-ground woody mass only. Roots and soil, often the bigger carbon store in some ecosystems, aren’t measured here.
• Change from one season to the next. Biomass moves slowly; GEDI’s strength is the multi-year structural snapshot, not a month-to-month signal.
• Exact carbon. The ~0.47 carbon-fraction is a rule of thumb; real wood density varies by species, and each shot carries its own measurement uncertainty (the agbd_se field).
• Steep or cloud-covered ground reliably. Laser returns get noisy on slopes, and the quality flag throws out many such shots, so coverage is uneven across rugged terrain.

Gotchas to watch for

• It’s samples, not a map. GEDI gives you dots along orbit lines, not a filled-in image. If your box has few tracks, you’ll have few shots. Pool a multi-year window so enough footprints accumulate to trust the median.
• You must filter for quality. Raw files include shots that failed (bad slope, clouds, weak return). Keep only those with the good flag (l4_quality_flag == 1) and a positive biomass value, or junk readings drag your numbers off. Any decent notebook does this automatically.
• The file format is fiddly. Each granule is an HDF5 file split into per-laser groups (BEAM0000…BEAM1011); the biomass lives in the agbd variable and each shot’s position in lat_lowestmode / lon_lowestmode. You loop over the beams and stitch the good shots together. The code below shows the pattern.
• “Carbon” is an estimate, not a measurement. Reporting “biomass × 0.47” is fine as a convention, but say so. Don’t present it as a measured carbon value.

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: search GEDI, open the HDF5 beam groups, filter for quality shots, take the median biomass, find the densest footprints, and map them. It needs a free Earthdata Login.

import earthaccess, h5py, numpy as np
import matplotlib.pyplot as plt

earthaccess.login(strategy="netrc")          # free Earthdata Login

W, S, E, N = -63.2, -10.2, -62.2, -9.2       # your forest as (W, S, E, N) — here, Rondônia, Amazon

# 1. Find GEDI L4A granules over the box, pooled across a multi-year window
#    (biomass moves slowly, so several years gather enough laser footprints).
granules = earthaccess.search_data(
    short_name="GEDI_L4A_AGB_Density_V2_1_2056",
    temporal=("2020-06-01", "2022-06-30"),
    bounding_box=(W, S, E, N),
)
print("granules found:", len(granules))

# 2. Open each HDF5 file, loop its per-laser beam groups, and keep only good shots:
#    quality flag == 1, positive biomass, and the footprint actually inside the box.
agbd_all, lat_all, lon_all = [], [], []
for g in granules[:12]:                       # a handful of granules is plenty for a median
    f = h5py.File(earthaccess.open([g])[0], "r")
    for beam in [k for k in f.keys() if k.startswith("BEAM")]:
        grp = f[beam]
        if "agbd" not in grp or "l4_quality_flag" not in grp:
            continue
        agbd = grp["agbd"][:]                 # above-ground biomass density, Mg/ha
        qf   = grp["l4_quality_flag"][:]      # 1 = good shot
        lat  = grp["lat_lowestmode"][:]
        lon  = grp["lon_lowestmode"][:]
        keep = (qf == 1) & (agbd > 0) & (lat >= S) & (lat <= N) & (lon >= W) & (lon <= E)
        agbd_all.append(agbd[keep]); lat_all.append(lat[keep]); lon_all.append(lon[keep])
    f.close()

agbd = np.concatenate(agbd_all)
lat  = np.concatenate(lat_all)
lon  = np.concatenate(lon_all)

# 3. The typical value is the median of the good shots; carbon ≈ biomass × 0.47.
med = float(np.median(agbd))
p90 = float(np.percentile(agbd, 90))          # threshold for the densest 10% of footprints
dense = agbd >= p90
print(f"{agbd.size} quality shots | median {med:.1f} Mg/ha "
      f"(~{med * 0.47:.1f} t C/ha) | densest 10% above {p90:.1f} Mg/ha")

# 4. Map the footprints, coloured by biomass; ring the densest stands.
fig, ax = plt.subplots()
sc = ax.scatter(lon, lat, c=agbd, s=6, cmap="YlGn", vmax=p90)
ax.scatter(lon[dense], lat[dense], s=30, facecolors="none", edgecolors="red",
           label="densest 10%")
plt.colorbar(sc, label="AGBD (Mg/ha)")
ax.set_xlabel("lon"); ax.set_ylabel("lat"); ax.legend()
ax.set_title(f"GEDI biomass, median {med:.0f} Mg/ha"); plt.tight_layout(); plt.show()

# Before you trust it: pair the low-biomass shots with HLSL30 30 m imagery to confirm they're real
# clearings/roads, and remember GEDI samples along tracks — this is "typical & densest sampled," not
# a wall-to-wall total (see the gotchas above).

Where the data comes from

All from NASA through one free login (Earthdata Login). The GEDI L4A biomass product comes from the ORNL DAAC (the archive for biosphere and carbon data), and the optional HLSL30 imagery you pair it with comes from the LP DAAC (the land-imagery archive). The free geoBoundaries state and country outlines (CC-BY) are separate and just help you draw and label your area.