WASR 4500/6500 · Quantitative Methods in Hydrology · PBA 1 · Evaristo Critical Zone Hydrology Lab · University of Georgia

SEBAL over the Flint — turning one thermal image into a map of water

A satellite measures temperature, not evaporation. This page walks the SEBAL algorithm step by step — the same five modules as your MATLAB exercise — over a center-pivot landscape in the Lower Flint River Basin of southwest Georgia. Every number here is computed live, with the equations from the handout.

What the satellite sees: surface temperature (K)
What SEBAL infers: daily ET (mm d⁻¹)

Cool pixels are wet, hot pixels are dry — SEBAL turns that single idea, calibrated by two anchor pixels, into a full evapotranspiration map. By the end of this page you will have run the entire chain yourself.

§0

Why evapotranspiration, and why here

The study region: the Lower Flint River Basin and the Dougherty Plain, southwest Georgia — the most intensively irrigated landscape in the state.

Evapotranspiration (ET) is the hardest term of the terrestrial water budget to measure. A rain gauge catches rain and a stream gauge integrates runoff, but the water that leaves as vapor — from soil, open water, and plant transpiration — crosses no weir. Yet in a Georgia summer it is the largest outgoing term: on a June day a well-watered corn pivot returns 5–7 mm of water to the atmosphere, roughly 55,000 liters per hectare, every day.

In the Lower Flint Basin that flux is the heart of a three-decade water conflict. Several thousand center pivots irrigate on the order of half a million acres of corn, cotton, and peanuts from the Upper Floridan aquifer — a karst limestone aquifer lying within tens of feet of the surface here, in direct hydraulic connection with the Flint River and its spring-fed creeks. Pumping shows up quickly as lost baseflow: in the 2011 drought the USGS mapped roughly 300 miles of losing streams and 600 miles of dry streambed in the lower basin. Those depleted flows sit at the center of the Apalachicola–Chattahoochee–Flint (ACF) dispute, argued all the way to the U.S. Supreme Court (Florida v. Georgia, decided 2021). Managing the system — the Flint River Drought Protection Act, irrigation permits, aquifer models — requires knowing how much water the landscape actually evaporates, field by field.

No network of flux towers can cover every field. But a satellite can see every field — in the thermal infrared. A transpiring canopy is evaporatively cooled; a dry fallow field bakes. SEBAL (Surface Energy Balance Algorithm for Land; Bastiaanssen et al., 1998) converts that temperature contrast into an ET map using the surface energy balance, needing only one satellite scene and one weather station.

~½ M
acres under center-pivot irrigation in the Lower Flint Basin — one of the densest pivot concentrations in the eastern U.S.
5–6 mm d⁻¹
typical June grass-reference ET in southwest Georgia (UGA weather network, FAO-56)
1–2 days
MODIS revisit time — Terra crosses Georgia mid-morning to local noon every day
2
calibration points SEBAL needs: one hot, dry pixel and one cold, wet pixel
The key principle, up front

SEBAL assumes the near-surface temperature difference dT = Ts − Ta varies linearly with surface temperature Ts across the image. Two extreme pixels — a hot, dry field where all available energy becomes sensible heat, and a cold, wet crop where nearly all of it evaporates water — pin down that line. Everything else on this page is bookkeeping around that one idea.

§1

The surface energy balance

One conservation law does all the work: the energy arriving at the land surface must go somewhere.

Eq. 1 λET = RnHG

Net radiation Rn is the energy supply. It is spent three ways: warming the air (H, sensible heat), warming the soil (G, ground heat), and evaporating water (λET, latent heat — λ ≈ 2.45 MJ kg⁻¹ is the energy to vaporize a kilogram of water). SEBAL measures nothing about water directly: it estimates Rn, G, and H from the satellite image, and whatever energy is left over must have evaporated water. λET is a residual.

Partition the energy yourself

A clear June midday in Mitchell County delivers about 620 W m⁻² of net radiation. Drag H and G and watch what is left for evaporation. Notice the residual logic — and what happens if you claim more heating than the energy supply allows.

Where Rn comes from

Eq. 2 Rn = (1 − α)RS↓ + RL↓RL↑ − (1 − ε0)RL↓

The radiation budget, live

Four streams cross the surface: incoming and reflected shortwave (sun), incoming and emitted longwave (thermal). The satellite supplies the two surface properties that control them — albedo α (how reflective) and emissivity ε₀ — plus the surface temperature that sets RL↑ = ε₀σTs⁴. Slide them and watch the budget rebalance.

Reading it like a hydrologist

Bright surfaces (high α) reject solar energy; hot surfaces (high Ts) shed longwave energy. That is why the dry fallow field, despite receiving the same sun, ends up with less net radiation than the irrigated pivot next door — and spends nearly all of what remains heating the air.

E1

Reference ET from a weather station

E1_ref_ET.mlx

Before touching the satellite image, SEBAL needs the atmosphere's demand for water: hourly FAO-56 Penman–Monteith reference ET from a ground station.

ETr is the ET of a hypothetical well-watered reference grass (0.12 m tall, surface resistance rs = 70 s m⁻¹, albedo 0.23) — the Penman–Monteith equation applied to a standard surface. In this Georgia edition the weather data come from the UGA Weather Network station at Camilla, Mitchell County (31.280°N, −84.292°E, 49 m) — fittingly, it sits at UGA's Stripling Irrigation Research Park, in the middle of the pivot country you will map below. The station records the same variables as the MATLAB exercise's meteodata.mat: air temperature, humidity, wind at 2 m, solar radiation, and pressure — and the network even publishes its own computed daily ET you can check against.

Eq. 8 ETr = Δ(RnG) + ρacp(esea)/rahΔ + γ(1 + rs/rah) energy term + demand term, throttled by two resistances

A June day at Camilla — build it, then break it

Defaults reproduce a clear 10 June day (DOY 161) at Camilla. The engine generates the diurnal weather, then computes hourly ETr exactly as E1_ref_ET does. Try: a windless muggy day (wind ↓, humidity ↑) vs. a dry cold-front day (humidity ↓, wind ↑). Which lever moves daily ETr most?

The dashed line marks the satellite overpass — the hour whose ETr becomes the scaling reference in module E5. Reset any time by reloading the page.

Code notes: what to watch for in E1_ref_ET.mlx
  • Ground heat under grass: both this page and E1_ref_ET use the FAO-56 hourly convention G = 0.1·Rn (daytime) and 0.5·Rn (night). Why is a fixed fraction defensible over a clipped reference grass — and for what surfaces would it fail?
  • Aerodynamic resistance: ra = 208/u2 is the FAO-56 closed form for the 0.12-m reference grass — it is not valid for forest or bare soil.
  • Units discipline: Rn arrives in W m⁻², vapor pressures in kPa, ET leaves in mm h⁻¹ — trace the 3600 and 10⁶ factors through Eq. 8 once by hand, and confirm the aerodynamic term ρacp(es−ea)/rah really is W m⁻² (cp in J kg⁻¹ °C⁻¹). Every SEBAL bug you will ever write is a unit bug.
E2

The scene, net radiation, and the two anchor pixels

E2_anchor_pixels.mlx

Now the satellite: a MODIS-style scene of the Lower Flint landscape at the moment of overpass, and the two pixels that calibrate everything.

The scene below is a synthetic training landscape — generated by this page, but built to the June statistics of Mitchell County: the Flint River and its cypress fringe on the east, loblolly-pine blocks, dryland fields, the town of Camilla, and the tell-tale circles of center pivots pumping from the Floridan aquifer. (Your MATLAB exercise does the same thing with a real MODIS scene; §9 shows where to download one for Georgia.) Hover any map — every pixel carries the full set of MODIS variables.

Choose the anchors

The sensible-heat equation (E3) has two unknowns per pixel. SEBAL closes the problem with two boundary conditions, taken at extreme cropland pixels in homogeneous neighborhoods:

Eq. 6–7 Hhot = RnG (dry bare field: λET = 0)     Hcold = RnG − 1.05·λETr (irrigated crop: ET = 105% of reference)

Pick a cold and a hot pixel interactive

Select a mode, then click the temperature map. The rules from E2_anchor_pixels: the pixel and its 3×3 neighborhood must be cropland (the river is colder than any field — and exactly the wrong choice; the widget will tell you why). Or let the algorithm choose, as the script does: coldest and hottest homogeneous cropland pixels. (Keyboard users: the algorithm's-picks button applies those same selection rules.)

Land surface temperature
K

Under the hood: Rn and G pixel by pixel

Shortwave in: Rs↓ = G_sc·d_r·τ_sw·cos(θ_sz) with transmissivity τ_sw = 0.75 + 2·10⁻⁵·z from the DEM. Longwave out: RL↑ = ε₀σTs⁴. Longwave in uses an air temperature tied to the cold pixel with a lapse correction. Ground heat uses the Bastiaanssen relation (handout Eq. 3):

Eq. 3 G/Rn = Ts [°C]α(0.0038α + 0.0074α²)(1 − 0.98·NDVI⁴)

Dense canopy (NDVI → 1) shuts G down; hot bare soil sends 20–30% of Rn into the ground. Check that your MATLAB implementation reproduces Eq. 3 exactly — including the division by α.

E3

The dT line and a first guess at H

E3_initial_H_calc.mlx

Two anchor points define a line; the line gives every pixel a temperature difference; the temperature difference gives H.

Eq. 5, 9 H = ρacp·dTrah   with   dT = a·Ts + b

We cannot measure air temperature at every pixel, so SEBAL never uses dT directly. Instead: at the anchors we know H (from E2), and we can compute the aerodynamic resistance rah from the wind profile — so we can invert Eq. 5 for dT at exactly two points. A line through them (Eq. 9) assigns dT to every other pixel from its surface temperature alone. That is the entire trick of SEBAL.

Wind enters through a blending height: the station wind at 2 m is projected up to 200 m — high enough that the wind no longer feels individual fields — then brought back down over each pixel using that pixel's roughness length zom (from land cover: tall pine is rough, water is glassy smooth).

The calibration line, drawn from your anchors

This scatter shows every pixel of the scene in (Ts, dT) space, with your anchors as the two fixed points. Go back to E2 and move an anchor — the line, and therefore every H in the scene, follows. (This is the first, neutral-atmosphere guess; E4 will bend it.)

Roughness lengths by land cover (zom)

Cropland is dynamic: z_om = 0.018·LAI — a leafy pivot is ~15× rougher than a freshly planted field. Note what this does to rah over the pine blocks.

E4

Stability corrections — the loop

E4_stability_correction.mlx · StabCorr.m

On a June midday the air over a hot field boils upward. Monin–Obukhov theory corrects the resistances for that buoyancy — but the correction depends on H, which depends on the correction. So: iterate.

M–O length L = −ρacp·u*³·Tsk g H the height at which buoyant turbulence overtakes wind-shear turbulence

Run the iteration interactive

Each pass: compute L per pixel → apply ψ corrections to u* and rah → re-invert the anchors for dT → new line → new H → repeat, until dThot and rah,hot stop moving (<1% change). The page arrives already converged (E5 needs the final state) — press Reset to neutral, then step through the passes to watch it happen.

Why the values move the way they do

Over the hot field L < 0 and |L| is small: buoyant plumes carry heat far more efficiently than the neutral log-profile assumed in E3, so rah drops on the first pass — and since Hhot is fixed by the energy balance, dThot = H·rah/(ρcp) drops with it. The first pass overshoots, then the adjustments shrink in a damped oscillation; on this scene the <1% rule lands at about seven passes — which is exactly why the MATLAB exercise hard-codes seven.

E5

Closing the balance: from one instant to a daily map

E5_ET_final.mlx

λET falls out as the residual at the overpass instant; the reference-ET fraction carries it through the day.

Eq. 10–12 ETins = 3600·λETλ·10⁶ mm h⁻¹     ETrF = ETinsETr,overpass     ET24 = ETr,24 · ETrF mm d⁻¹

The satellite sees the landscape for one instant (). The reference-ET fraction ETrF — this pixel's ET as a fraction of the reference crop's — is assumed constant through daylight, so multiplying by the day's total ETr (from E1) yields daily ET. ETrF is the satellite's version of a crop coefficient.

The full energy-balance suite

These update live with your anchors (E2), the weather (E1), and the iteration state (E4). Hover any pixel to read its complete budget — the inspector shows Rn − G − H = λET closing exactly.

Does it make sense? ET by land cover

Scene-average daily ET for each land cover, against the ranges we expect for southwest Georgia in June from lysimeter, eddy-covariance, and OpenET studies. This is validation thinking — an algorithm can run perfectly and still be wrong.

§7

Sensitivity & uncertainty — is the map trustworthy?

sensitivity_uncertainty.mlx

Learning goal 3: distinguish "the code ran" from "the result is valid." Two experiments.

Experiment 1 — the anchors are a choice

Every anchor pair you have tried on this page is logged here with the scene-mean ET it produced. Go back to E2, pick deliberately bad anchors (a mild "hot" pixel; a half-stressed "cold" pivot), and watch the whole map shift. In METRIC calibration, analyst anchor choice is one of the largest single sources of uncertainty.

TrialTs,cold (K)Ts,hot (K)slope a Scene ET̄₂₄ (mm/d)Pivot ET̄₂₄vs. algorithm pick

Experiment 2 — Monte Carlo error propagation

MODIS LST is good to about ±1 K; albedo to ±0.02; your station ETr to a few percent. Perturb the inputs within those errors, re-run all of SEBAL each time, and the spread of scene-mean ET₂₄ is your uncertainty band.

Failure modes to diagnose in your own runs
  • Anchor contamination: a "cold" pixel that is actually water or wetland — wrong roughness, wrong heat storage, and the 1.05·ETr boundary condition no longer holds.
  • Energy-balance non-closure: pixels where λET clamps at zero mean H > Rn − G — the dT line is over-extrapolated there.
  • Tall-canopy bias: the dT–Ts line is calibrated on cropland; over the pine blocks the roughness is 20× larger and the linear model is least defensible (this motivates the trapezoid models — Long & Singh 2012; Yang & Shang 2013 in your readings).
  • Advection: ETrF constant-through-the-day fails on days when dry air blows in mid-afternoon.
§8

Your turn: nine stations across Georgia

The handout's closing exercise: run the pipeline for UGA weather-network stations spanning Georgia's physiographic provinces — then explain the geography in the numbers.

Georgia climbs from sea level on Jekyll Island to nearly 600 m in the Blue Ridge at Blairsville. Along that climb, June radiation, humidity, and temperature all shift — and so does everything SEBAL computes. Before you reveal each row, predict: will daily ETr at this station be higher or lower than Camilla's? Then check yourself, and work through the discussion points below.

Ranges are June-typical values computed with this page's FAO-56 engine from each station's June-normal weather (UGA weather network climatology), ± the day-to-day spread. They are estimates for prediction-checking, not measurements.

Discussion points in class
  • Coastal Plain vs. Blue Ridge: which drives the ETr difference more — temperature, humidity, or radiation? Use the E1 sliders to test your claim.
  • Two stations at nearly the same latitude can differ in ETrF for the same crop. What surface properties could do that?
  • The satellite overpass time is essentially the same statewide. Why does that make the ETrF method portable across the state?
§9

Georgia data sources — build this with real data

Everything the MATLAB exercise ships as .mat files, you can pull for Georgia yourself. Updated product IDs replace the deprecated ones in older SEBAL guides.

Check yourself against the professionals

Two operational products publish satellite ET for every Georgia field: OpenET (ensemble of six models including the METRIC descendant of SEBAL, 30-m resolution, openetdata.org) and the USGS SSEBop actual-ET product. After you run SEBAL on a real Georgia scene, comparing your map against theirs is the single fastest sanity check you can do.

§10

Self-check

Answer before peeking. Each maps to a handout question or a step you just ran.