HERE Technologies and Wherobots announce partnership Learn More

El Niño 2026, atmospheric rivers, and California’s burn scars: one SQL engine, two data models

Authors

Governor Gavin Newsom proclaimed a state of emergency on 21 September 2026 to prepare California for El Niño 2026. It names recent wildfire areas as particularly susceptible to mud and debris flows. Above Altadena, the Eaton fire left 42 basins that the US Geological Survey (USGS) rates high hazard for debris flows. USGS computes the rating for a 15-minute burst of 40 mm/h of rain (1.57 in/hour), and an atmospheric river can deliver that kind of deluge. The rating covers the basins on the hillside. The homes and roads downslope of the basins fall outside the rating, and the queries here count them.

Key takeaways

  • A raster vector join matches raster pixels to vector geometry by location. It is the case of a spatial join where one side is raster (elevation) and the other is vector (building and road geometry).
  • WherobotsDB runs both joins in SQL. It reads Overture buildings and roads, the ASTER GDEM elevation model, Regrid parcels and the USGS hazard basins in place, as vector and raster, with no data movement.
  • 7,713 building footprints sit inside the 42 Eaton high-hazard basins or within 500 m (1,640 ft).
  • 696 sit inside a basin, and 16,748 sit within 1 km (0.6 mi).
  • 388 km (241 mi) of road and path sit in the 16.8 km² (6.5 sq mi) of ground below the scar’s lowest point. 305 km (190 mi) is drivable and 24 km (15 mi) is motorway.
  • Joined to Regrid parcels, the footprints near the basins sit on 5,934 parcels assessed at $4.8 billion on the 2023 tax roll.
  • A bounding-box filter before the spatial test cut the first building count to 22.8 seconds. Without it, the run was cancelled at 7% after three minutes.
  • The same queries extend to every recent burn scar in California and across the western United States, where USGS publishes debris-flow hazard assessments after major fires.

WherobotsDB does the counting in SQL, reading every dataset in place. A join to the USGS basins returns 7,713 mapped building footprints inside the high-hazard basins or within 500 m (1,640 ft) of them. A raster vector join over elevation returns 388 km (241 mi) of mapped road and path.

Start with elevation, buildings and roads in the catalog

All four datasets sit in wherobots_open_data: the ASTER GDEM 30 m (98 ft) elevation model, Overture’s buildings and roads, and Regrid’s parcel sample. The companion notebook fetches the CAL FIRE perimeter and USGS basins. It also registers every shared view, so each query below runs as written, with no setup of its own. Each elevation row stores a footprint polygon, so WherobotsDB selects rows by location before it reads any pixel. Ten rasters cover the Eaton area.

Filter on the bounding box before the spatial test

Overture features carry a bbox struct beside the geometry. With a range filter on it before ST_Intersects, the first building count finished in 22.8 seconds. Without it, the run was cancelled at 7% after three minutes. The range check on four numbers is cheap, so the costly geometry test runs only on the rows that pass.

Show the code
-- scar holds the Eaton perimeter; the notebook fetches it from CAL FIRE and registers it as a view
WITH aoi AS (
  SELECT ST_Transform(
           ST_Buffer(ST_Transform(g, 'EPSG:4326', 'EPSG:32611'), 3000),
           'EPSG:32611', 'EPSG:4326') AS g
  FROM scar
),
bld AS (
  SELECT t.id, t.geometry AS geom
  FROM wherobots_open_data.overture_maps_foundation.buildings_building t
  WHERE t.bbox.xmin BETWEEN -118.20 AND -117.97   -- range scan first
    AND t.bbox.ymin BETWEEN  34.13  AND  34.27
)
SELECT count(*) FROM bld b JOIN aoi ON ST_Intersects(b.geom, aoi.g);

118,746 footprints pass the bounding box, 83,930 fall within 3 km (1.9 mi) of the perimeter, and 10,652 sit inside the scar.

Each of two questions below gets its own raster vector join. For buildings the answer is one elevation per footprint, so the first join reads the raster under every footprint. For roads, however, the answer is a length. So the second join turns the raster into a zone polygon first, then clips each road segment against it.

Buildings read the elevation under each footprint a b c elevation + building footprints RS_ZonalStats building elevation a 222 m b 285 m c 201 m one value per building, kept if below 230 m Roads turn the elevation into a zone, then clip each road 1 2 elevation below 230 m zone + roads 1 RS_MapAlgebra marks the cells below the line 2 RS_Polygonize turns them into polygons, ST_Union_Aggr merges them into one zone, and ST_Intersection measures the part of each road inside it
The two joins used below, drawn on a toy elevation grid with a 230 m (755 ft) line.

Buildings: read the elevation under each footprint

Which building footprints sit below the burn scar’s lowest point?

15,233 of the 83,930 footprints sit below that point, 223 m (732 ft). RS_ZonalStats runs the raster vector join per building, so every footprint gets one elevation value of its own. That value becomes one more column beside the building’s other attributes.

Show the code
-- scar, dem and inaoi (the 83,930 footprints within 3 km) are views from the notebook's setup cell
WITH scar_floor AS (
  SELECT min(RS_ZonalStats(d.rast, s.g, 1, 'min', true)) AS m FROM dem d, scar s
),
be AS (
  -- pool pixels across every raster a footprint touches: total value over pixel count
  SELECT i.id,
         sum(RS_ZonalStats(d.rast, i.geom, 1, 'sum', true))
           / NULLIF(sum(RS_ZonalStats(d.rast, i.geom, 1, 'count', true)), 0) AS elev
  FROM inaoi i JOIN dem d ON RS_Intersects(d.rast, i.geom)
  GROUP BY i.id
)
SELECT count(*) AS below_scar_floor
FROM be, scar_floor
WHERE be.elev < scar_floor.m;

Footprints that straddle two rasters have their pixels pooled before averaging, since the rasters meet under them. 380 of the 384 stay on the same side of the line. A 30 m (98 ft) pixel is wider than many houses, but every footprint still got an elevation value from the query.

2026-09-24T20:53:29.930329 image/svg+xml Matplotlib v3.11.2, https://matplotlib.org/ 500 750 1,000 1,250 1,500 1,750 2,000 2,250 feet 200 300 400 500 600 700 mean elevation under the footprint (m) 0 2,000 4,000 6,000 8,000 10,000 buildings per 20 m (66 ft) band scar’s lowest point, 223 m (732 ft) 15,233 below it 146 more sit above 700 m (2,297 ft), up to 1,809 m (5,935 ft), off the chart
Mean elevation under each of the 83,930 footprints within 3 km (1.9 mi) of the perimeter, in 20 m (66 ft) bands. Blue marks the buildings below the scar’s lowest point.

Roads: a raster vector join through a low-ground zone

Why do roads need a different join than buildings?

Roads need a different join because length, not a single point, decides the answer. A building gets one average elevation, and that single value places it above or below the line. A 1 km (0.6 mi) road dips below the line for a short stretch and stays above it for the rest. An average across the whole road would place it on one side, and the low stretch would go uncounted. The query takes three steps instead. RS_MapAlgebra marks every elevation pixel below 223 m (732 ft). RS_Polygonize turns the marked pixels into polygons, and ST_Union_Aggr merges them into one shape of low ground. ST_Intersection then cuts each road at the polygon’s edge, and the query adds up the length inside it. WherobotsDB builds the polygon once and tests every road against that same shape.

Show the code
-- aoi, dem and roads are views from the notebook's setup cell
WITH masked AS (
  SELECT RS_MapAlgebra(rast, 'D', 'out[0] = rast[0] < 223 ? 1 : 0;') AS m FROM dem
),
polys AS (
  SELECT t.p.geom AS g
  FROM masked LATERAL VIEW explode(RS_Polygonize(m, 1)) t AS p
  WHERE t.p.value = 1.0
),
low AS (
  SELECT ST_Intersection(u.g, a.g) AS g
  FROM (SELECT ST_Union_Aggr(ST_MakeValid(g)) AS g FROM polys) u, aoi a
)
SELECT r.class,
       sum(ST_Length(ST_Transform(ST_Intersection(r.geom, low.g),
                                  'EPSG:4326', 'EPSG:32611'))) / 1000 AS km
FROM roads r JOIN low ON ST_Intersects(r.geom, low.g)
GROUP BY r.class
ORDER BY km DESC;
2026-09-24T20:53:30.030379 image/svg+xml Matplotlib v3.11.2, https://matplotlib.org/ 0 20 40 60 80 100 120 140 160 kilometres inside the low zone Residential Service Footway Motorway Primary Tertiary Secondary Other 121.0 km (75.2 mi) 97.6 km (60.6 mi) 60.2 km (37.4 mi) 24.0 km (14.9 mi) 22.1 km (13.7 mi) 20.1 km (12.5 mi) 19.8 km (12.3 mi) 23.0 km (14.3 mi)
Kilometres of Overture road and path inside the 16.8 km² (6.5 sq mi) of ground below the scar’s lowest point, within 3 km (1.9 mi) of the perimeter.

The low ground covers 16.8 km² (6.5 sq mi) and holds 388 km (241 mi) of road and path within 3 km (1.9 mi) of the scar. Of that, 305 km (190 mi) is drivable and 24 km (15 mi) is motorway. Lengths follow each mapped segment, so both carriageways of every divided road are counted in full.

Joining to what USGS already modelled

An elevation line is only a screen. USGS publishes the basins as a live feature service, so they join in as one more vector layer. Each building counts once per distance, however many basins it touches.

Show the code
-- high holds the 42 dissolved high-hazard basins; the notebook fetches them from USGS
WITH z AS (
  SELECT g AS inside,
         ST_Transform(ST_Buffer(ST_Transform(g, 'EPSG:4326', 'EPSG:32611'), 500),  'EPSG:32611', 'EPSG:4326') AS b500,
         ST_Transform(ST_Buffer(ST_Transform(g, 'EPSG:4326', 'EPSG:32611'), 1000), 'EPSG:32611', 'EPSG:4326') AS b1k
  FROM high
)
SELECT sum(CASE WHEN ST_Intersects(b.geom, z.inside) THEN 1 ELSE 0 END) AS inside,
       sum(CASE WHEN ST_Intersects(b.geom, z.b500)   THEN 1 ELSE 0 END) AS inside_or_500m,
       sum(CASE WHEN ST_Intersects(b.geom, z.b1k)    THEN 1 ELSE 0 END) AS inside_or_1km
FROM bld b, z
WHERE ST_Intersects(b.geom, z.b1k);
2026-09-24T20:53:30.095776 image/svg+xml Matplotlib v3.11.2, https://matplotlib.org/ 0 2,500 5,000 7,500 10,000 12,500 15,000 17,500 20,000 Overture building footprints Inside a high-hazard basin Inside or within 500 m (1,640 ft) Inside or within 1 km (0.6 mi) 696 7,713 16,748
Overture building footprints inside, or within 500 m (1,640 ft) and 1 km (0.6 mi) of, the 42 Eaton basins that USGS rates high combined hazard.

696 footprints sit inside a high-hazard basin, 7,713 inside or within 500 m (1,640 ft), and 16,748 within 1 km (0.6 mi). The map below plots the 7,713 on Sentinel-2 imagery from 5 August 2026. They form a solid band along the scar’s lower edge, where the basins drain onto the valley floor. Nineteen months on, the scar is still bare.

Raster vector join result: 7,713 building footprints near 42 high-hazard debris-flow basins below the Eaton burn scar in Altadena

The 7,713 footprints (cyan) inside or within 500 m (1,640 ft) of the 42 high-hazard basins (red), on Sentinel-2 imagery from 5 August 2026.

Joining parcels to put a value on the screen

Building counts measure how many structures sit in each screen, not what they are worth. Regrid’s parcel sample in wherobots_open_data.partner_samples adds use codes, acreage and assessed values, and a join to the exposed footprints carries those onto both screens. The join runs in two steps so nothing is counted twice: first the distinct parcel and footprint pairs, keyed on Regrid’s stable ll_uuid, then one row per parcel for its use code, acreage and assessed value. A footprint that crosses a parcel line counts once in the total and once in each use it touches, so the rows sum slightly above the total.

Show the code
-- exposed holds the footprints from either screen; parcels is Regrid's sample with the same bbox filter as bld
WITH pairs AS (
  SELECT DISTINCT p.ll_uuid, p.usedesc, e.id AS bid
  FROM parcels p JOIN exposed e ON ST_Intersects(p.geometry, e.geom)
),
parcel AS (  -- one row per parcel, so acres and value are counted once
  SELECT ll_uuid, first(usedesc) AS usedesc, first(ll_gisacre) AS acres, first(parval) AS parval
  FROM parcels WHERE ll_uuid IN (SELECT ll_uuid FROM pairs) GROUP BY ll_uuid
),
by_use AS (
  SELECT usedesc, count(*) AS parcels, round(sum(acres)) AS acres, round(sum(parval) / 1e6, 1) AS assessed_musd
  FROM parcel GROUP BY usedesc
),
bld_use AS (SELECT usedesc, count(distinct bid) AS buildings FROM pairs GROUP BY usedesc)
SELECT u.usedesc, u.parcels, b.buildings, u.acres, u.assessed_musd
FROM by_use u JOIN bld_use b USING (usedesc) ORDER BY u.parcels DESC;

Inside or within 500 m (1,640 ft) of a high-hazard basin: 5,934 parcels, 7,709 footprints

Use Parcels Buildings Acres Assessed value
Residential single 5,686 7,157 2,910 $4.69 B
Residential two units 72 186 34 $51 M
Government parcel 67 175 944 $4 M
Residential three units 29 86 14 $24 M
Utility pumping plants 15 33 25 $0.1 M
Apartments, five or more 14 82 10 $13 M
Miscellaneous 10 48 3,383 $0
Motion picture, radio, TV 8 52 868 $12 M

Below the scar’s lowest point, 223 m (732 ft): 10,513 parcels, 15,227 footprints

Use Parcels Buildings Acres Assessed value
Residential single 9,032 12,791 3,013 $8.49 B
Residential two units 373 988 77 $210 M
Apartments, five or more 172 498 72 $499 M
Commercial stores 129 187 50 $222 M
Commercial office buildings 105 184 62 $511 M
Residential three units 86 233 17 $56 M
Light manufacturing 84 123 27 $119 M
Residential four units 81 210 18 $66 M

The 7,713 footprints inside or within 500 m (1,640 ft) of a high-hazard basin sit on 5,934 parcels with a combined assessed value of $4.8 billion. Single-family homes account for $4.7 billion of it. The downhill screen is larger: the 15,233 footprints below the scar’s lowest point sit on 10,513 parcels assessed at $11.1 billion, and here the mix changes. Single-family homes are $8.5 billion of that total, and the remaining $2.6 billion is the commercial and apartment stock along the valley floor, led by $511 million in office buildings, $499 million in apartments and $222 million in stores. The two screens overlap, so the totals do not add. Every value comes from the assessor’s 2023 tax roll, so it describes the neighbourhood before the fire and before any rebuilding. Proposition 13 holds California assessed values well below market value, so these totals are a floor.

RS_Value reads the elevation under a point, so a cross-section is a line of generated points.

Show the code
-- a point every 0.00025 degrees, about 28 m, along 118.125° W, and the elevation under each
WITH pts AS (
  SELECT i, ST_SetSRID(ST_Point(-118.125, 34.265 - i * 0.00025), 4326) AS pt
  FROM (SELECT explode(sequence(0, 460)) AS i)
)
SELECT p.i, ST_Y(p.pt) AS lat, RS_Value(d.rast, p.pt, 1) AS elev
FROM pts p
JOIN wherobots_open_data.aster_gdem.v3_30m d ON ST_Intersects(d.footprint, p.pt)
ORDER BY p.i;
2026-09-24T20:53:30.305805 image/svg+xml Matplotlib v3.11.2, https://matplotlib.org/ 0 1 2 3 4 5 6 7 miles 0 1,000 2,000 3,000 4,000 5,000 feet 0 2 4 6 8 10 12 distance south along 118.125° W (km) 0 200 400 600 800 1000 1200 1400 1600 elevation (m) burn scar high-hazard basins 1,112 buildings within 60 m (197 ft) of the line north south
Terrain sampled every 28 m (92 ft) along 118.125° W, from the ridge north of the scar to the valley floor. Red bands mark the burn scar and the basins USGS rates high; cyan ticks mark buildings within 60 m (197 ft) of the line.

The line climbs to 1,577 m (5,174 ft) on the ridge. Then it crosses 2.4 km (1.5 mi) of high-hazard basin and 4.5 km (2.8 mi) of scar, down to 248 m (814 ft). 1,112 footprints lie within 60 m (197 ft) of it. The first of those buildings already sit inside the lower end of the high-hazard stretch.

Scope and limits of the exposure screen

These counts are an exposure screen. Neither they nor the USGS assessment predict where a debris flow travels or what it buries. USGS notes that most come within two years of a fire. So this winter falls near the end of that window, while the burned slopes are still recovering.

Three limits carry over. The elevation model is a surface model, so rooftops and treetops count as ground. Proximity is straight-line distance rather than flow path. Finally, the perimeter is simplified to 47 vertices, which changes its area by 0.4%. Flow routing would turn distance bands into runout paths. However, WherobotsDB has no built-in flow-direction function, so that takes a custom function over the rasters. It is the natural next step.

From one burn scar to the western United States

This screen covers one burn scar: the Eaton fire above Altadena. The same queries extend to every recent burn scar in California, and across the western United States, where USGS publishes debris-flow hazard assessments after major fires. Overture buildings and roads and the ASTER elevation model in wherobots_open_data already cover that whole region. Each new fire adds a perimeter and a set of basins to the same joins, once the study-area bounding box and UTM zone are updated to match. Governor Newsom declared the emergency ahead of El Niño rain this winter. The open question for every recently burned community is how many homes and roads sit below a high-hazard basin when that rain arrives.

One engine processed buildings, roads, elevation, and USGS hazard polygons, as vector and raster, with every file left where it sits. An AI agent grounded only in text, documents, databases, and the internet holds none of those layers, so it cannot say what sits below a basin. WherobotsDB supplies that physical-world context by joining them in place.

Imagery: Copernicus Sentinel data 2026. Elevation: ASTER GDEM is a product of METI and NASA. Buildings and roads: Overture Maps Foundation. Parcels: Regrid parcel sample. Fire perimeter: CAL FIRE FRAP, CC BY. Hazard basins: USGS Landslide Hazards Program, public domain.

Get Started with Wherobots

Frequently Asked Questions

What is a raster vector join?

A raster vector join is a spatial join where one dataset is raster and the other is vector. Here the raster is the ASTER GDEM elevation model and the vectors are Overture building footprints and road segments. RS_ZonalStats reads one elevation value under each footprint, and RS_MapAlgebra with RS_Polygonize turns low-elevation pixels into a polygon that roads clip against.

What is a spatial join?

A spatial join matches records from two datasets by their location instead of a shared key column. It answers questions like which buildings fall inside a hazard basin, or which roads cross a zone of low ground. A raster vector join is the case where one side is raster.

How many buildings sit below the Eaton burn scar's high-hazard basins?

About 7,700 building footprints sit inside the 42 high-hazard basins or within 500 m (1,640 ft) of them. About 700 sit inside a basin, and about 16,700 sit within 1 km (0.6 mi).

Why filter on a bounding box before a raster vector join?

A bounding-box range check tests four numbers per row and is cheap, so the costly geometry test runs only on rows that pass. With the filter, the first building count finished in 22.8 seconds. Without it, the run was cancelled at 7% after three minutes.

Does this screen predict where a debris flow will go?

No. It is an exposure screen. Neither it nor the USGS assessment predicts where a debris flow travels or what it buries.

Flow routing would turn straight-line distance bands into runout paths, which takes a custom function over the rasters.