Skip to content

Commit f1a2c03

Browse files
alxmrsclaude
andcommitted
Docs: write for the library's reader, not the review thread
Revise the prose across the suite for the future reader/user — someone from the SQL or Xarray world — rather than for the PR conversation: - Frame features by what that audience values (lazy, idiomatic) instead of the mechanics that justified review decisions ("nothing loaded or column-selected up front", "bounds passed as query parameters, not formatted into the SQL", "no table name formatted into the SQL"). - Fix stale claims: the suite compares with xarray (not `numpy.assert_allclose`), Earth Engine uses your `gcloud` login (not `earthengine authenticate`), and the run instructions cover `run_all.sh` and the in-repo build. - Trim jargon from runtime labels and tighten a couple of code comments. Touches README ("Does it work?"), docs/geospatial.md, the suite README, and the 02–06 docstrings/labels. No code behavior changes. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01AWvrZYAT2NbuETBqNAN3o9
1 parent 7e3dbd1 commit f1a2c03

8 files changed

Lines changed: 55 additions & 71 deletions

File tree

README.md

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -152,8 +152,8 @@ against an xarray/array reference** to floating-point tolerance:
152152

153153
* **Spectral indices** (NDVI) — column arithmetic over a real Sentinel-2 scene.
154154
* **Climatology, anomalies, zonal means**`GROUP BY` and self-`JOIN` over the
155-
full 0.25° **ARCO-ERA5** archive (≈1.3M hourly steps), read lazily with
156-
partition pruning and column pushdown.
155+
full 0.25° **ARCO-ERA5** archive (≈1.3M hourly steps), read lazily so each
156+
query touches only the data it needs.
157157
* **Forecast skill** — scoring the **Pangu-Weather** and **GraphCast** ML models
158158
against ERA5 (WeatherBench 2) as a `JOIN` on `valid_time = init + lead`; it
159159
reproduces the published result that GraphCast beats Pangu at every lead.

benchmarks/geospatial/02_climatology.py

Lines changed: 8 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -28,11 +28,9 @@
2828
``da.groupby("time.hour").mean()``. ERA5 is hourly, so grouping by hour of day
2929
gives a clean 24-bin **diurnal cycle**, one sample per day in the window.
3030
31-
The table is the *whole* lazily-opened ARCO-ERA5 archive — nothing is loaded or
32-
pre-selected up front. The query reads only the ``2m_temperature`` column
33-
(projection pushdown) for only the window the parameterized ``WHERE`` asks for
34-
(partition pruning). The bounds are passed as **query parameters**, not
35-
formatted into the SQL string.
31+
The table is the *whole* ARCO-ERA5 archive, opened lazily: the query reads only
32+
``2m_temperature``, and only over the window its ``WHERE`` asks for — the rest of
33+
the archive is never touched.
3634
"""
3735

3836
from __future__ import annotations
@@ -68,10 +66,10 @@
6866

6967

7068
def main() -> None:
71-
# Open the full ARCO-ERA5 archive lazily — dask off, nothing loaded or
72-
# column-selected here. ERA5 mixes surface (time, lat, lon) and atmospheric
73-
# (… level …) variables, so register it as two tables under an ``era5``
74-
# schema; the query below touches only the surface table's 2m_temperature.
69+
# Open the full ARCO-ERA5 archive lazily — no data is read here. ERA5 mixes
70+
# surface (time, lat, lon) and atmospheric (… level …) variables, so register
71+
# it as two tables under an ``era5`` schema; the query below touches only the
72+
# surface table's 2m_temperature.
7573
try:
7674
import gcsfs # noqa: F401 — required by the gs:// protocol
7775

@@ -107,7 +105,7 @@ def main() -> None:
107105

108106
# A climatology is a gridded product: round-trip the result back to an
109107
# xarray Dataset keyed by (latitude, longitude, hour) — how it is used.
110-
with timed("SQL diurnal climatology (lazy read, pushdown + pruning)"):
108+
with timed("SQL diurnal climatology (lazy read)"):
111109
got = ctx.sql(sql, param_values=_PARAMS).to_dataset(
112110
dims=["latitude", "longitude", "hour"]
113111
)

benchmarks/geospatial/03_zonal_mean.py

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -86,8 +86,8 @@ def main() -> None:
8686
},
8787
)
8888

89-
# The window bounds are passed as query parameters, not formatted into the
90-
# SQL string; pruning still kicks in, so only one day is read.
89+
# Pass the day's bounds as query parameters; the query still reads only that
90+
# one day out of the whole archive.
9191
sql = """
9292
SELECT latitude,
9393
AVG("2m_temperature") - 273.15 AS air_mean_c
@@ -99,7 +99,7 @@ def main() -> None:
9999
show_sql(sql)
100100

101101
# Round-trip the profile back to an xarray Dataset keyed by latitude.
102-
with timed("SQL zonal mean (WHERE-pruned to one day)"):
102+
with timed("SQL zonal mean (reads one day)"):
103103
got = ctx.sql(
104104
sql, param_values={"start": _START, "end": _END}
105105
).to_dataset(dims=["latitude"])

benchmarks/geospatial/04_anomaly.py

Lines changed: 3 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -27,11 +27,9 @@
2727
FROM era5 a JOIN clim c
2828
ON (a.latitude, a.longitude, hour(a.time)) = (c.latitude, c.longitude, c.hour)
2929
30-
The table is the *whole* lazily-opened ARCO-ERA5 archive — nothing loaded or
31-
column-selected up front. Both the climatology CTE and the outer scan read only
32-
``2m_temperature`` (projection pushdown) over only the window the parameterized
33-
``WHERE`` asks for (partition pruning). Bounds are query parameters, not
34-
string-formatted into the SQL.
30+
The table is the *whole* ARCO-ERA5 archive, opened lazily: both the climatology
31+
CTE and the outer scan read only ``2m_temperature``, and only over the window the
32+
``WHERE`` asks for.
3533
"""
3634

3735
from __future__ import annotations

benchmarks/geospatial/05_forecast_skill.py

Lines changed: 4 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -35,9 +35,8 @@
3535
GROUP BY f.model, f.prediction_timedelta
3636
3737
We stack the two models along a ``model`` dimension into a single forecast
38-
table, so the query groups by a ``model`` *column* — no table name is formatted
39-
into the SQL. Nothing is loaded up front either: the forecasts and ERA5 are
40-
registered lazily, and the JOIN reads only what it needs at query time.
38+
table, so one query scores them together, grouped by the ``model`` column. The
39+
forecasts and ERA5 are opened lazily, and the JOIN reads only what it needs.
4140
4241
Datasets: WeatherBench 2 **Pangu**, **GraphCast**, and **ERA5** at a coarse
4342
64×32 grid (so the demo is small and fast), read from the public ``gs://
@@ -116,8 +115,8 @@ def main() -> None:
116115
# The two models store different pressure-level sets (Pangu 13, GraphCast
117116
# 37), so we keep the common surface field 2m_temperature and stack the
118117
# models along a `model` dimension into one forecast table. Snap the grid
119-
# onto ERA5's exact coordinates (same 64×32 grid) so the equality JOIN is
120-
# bit-safe across the Zarr stores.
118+
# onto ERA5's exact coordinates (same 64×32 grid) so the join on latitude and
119+
# longitude lines up exactly across the two Zarr stores.
121120
pangu = _open(_PANGU)[[_VAR]].sel(time=_INIT)
122121
graphcast = _open(_GRAPHCAST)[[_VAR]].sel(time=_INIT)
123122
forecasts = xr.concat([pangu, graphcast], dim="model").assign_coords(

benchmarks/geospatial/06_zonal_vector.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -127,7 +127,7 @@ def main() -> None:
127127
"""
128128
show_sql(sql)
129129

130-
with timed("SQL zonal stats (raster × vector range JOIN, WHERE-pruned)"):
130+
with timed("SQL zonal stats (raster × vector range JOIN)"):
131131
got = ctx.sql(
132132
sql, param_values={"start": _START, "end": _END}
133133
).to_dataset(dims=["region_id"])

benchmarks/geospatial/README.md

Lines changed: 19 additions & 27 deletions
Original file line numberDiff line numberDiff line change
@@ -4,9 +4,9 @@
44
are, underneath, **relational** operations — `GROUP BY`, `JOIN`, window
55
functions, and `CASE`. Each script here takes one such operation, expresses it
66
in SQL against [`xarray-sql`](../../README.md), and **proves the SQL answer
7-
matches an xarray/array reference implementation** (`numpy.assert_allclose`).
8-
Wall-clock and peak memory are reported too, but the headline is correctness +
9-
clarity of the SQL.
7+
matches a plain-xarray reference** to floating-point tolerance. Wall-clock and
8+
peak memory are reported too, but the headline is correctness + clarity of the
9+
SQL.
1010

1111
This suite is *expressibility-first*: the point is that the SQL reads like the
1212
plain-English definition of the operation, and computes the same numbers.
@@ -39,13 +39,10 @@ interpolation weights — the geometry — which SQL applies but does not comput
3939
(bands B04/B08). Requires network; skips cleanly if offline.
4040
- **02–06** — the full **[ARCO-ERA5](https://github.com/google-research/arco-era5)**
4141
archive (0.25° global, ~1.3M hourly timesteps, 273 variables) read anonymously
42-
from a public GCS bucket. Every case registers the *whole* archive **lazily**
43-
(nothing loaded or column-selected up front) and filters it with a
44-
**parameterized** `WHERE` (bounds bound as query parameters, not formatted into
45-
the SQL); projection pushdown reads only `2m_temperature` and partition pruning
46-
reads only the window asked for. All require network (`gcsfs`); skip cleanly
47-
offline. ERA5 cases take roughly one to a few minutes, dominated by the GCS
48-
read (the lazy reference re-reads the same window).
42+
from a public GCS bucket. Each case opens the *whole* archive lazily, so a query
43+
reads only the variable and the window it asks for — never the other 272
44+
variables or the rest of the timesteps. All require network (`gcsfs`) and skip
45+
cleanly offline; each takes roughly one to a few minutes, dominated by the read.
4946
- **05 forecast skill** — the **[WeatherBench 2](https://weatherbench2.readthedocs.io/)**
5047
Pangu-Weather, GraphCast, and ERA5 datasets at a coarse 64×32 grid, scoring
5148
both ML models against ERA5 ground truth. Network-backed; runs in seconds
@@ -54,29 +51,24 @@ interpolation weights — the geometry — which SQL applies but does not comput
5451
07 reprojects a UTM grid and validates the SQL transform against Earth Engine's
5552
*own* per-pixel lon/lat (`ee.Image.pixelLonLat()`) — an independent reprojection
5653
reference, not PROJ-vs-PROJ. 08 regrids real **SRTM elevation** (Sierra Nevada)
57-
and validates against xarray's bilinear `.interp()`. Both need EE access
58-
(`earthengine authenticate` + an initialized project via `EARTHENGINE_PROJECT`)
59-
and skip cleanly without it.
54+
and validates against xarray's bilinear `.interp()`. Both run against Earth
55+
Engine using your existing `gcloud` login, and skip cleanly without it.
6056

6157
## Running
6258

63-
Inside the repo, use the project environment (so `import xarray_sql` resolves to
64-
the locally built native extension):
59+
Run a single case, or the whole suite, from any directory:
6560

6661
```shell
67-
python benchmarks/geospatial/03_zonal_mean.py
62+
uv run benchmarks/geospatial/03_zonal_mean.py # one case
63+
benchmarks/geospatial/run_all.sh # all of them
6864
```
6965

70-
Each script also carries [PEP 723 / `uv` inline metadata](https://docs.astral.sh/uv/guides/scripts/),
71-
so it can be run standalone against the published `xarray-sql` wheel:
66+
Each script carries [PEP 723 / `uv` inline metadata](https://docs.astral.sh/uv/guides/scripts/)
67+
and runs against the `xarray-sql` in this checkout.
7268

73-
```shell
74-
uv run benchmarks/geospatial/03_zonal_mean.py
75-
```
76-
77-
A passing case prints a `✅ … SQL matches array reference` line; a mismatch
78-
raises `AssertionError` and exits non-zero. Cases that need an unavailable
79-
dataset/dependency print `⏭ SKIPPED` and exit 0.
69+
A passing case prints a `✅ … SQL matches xarray reference` line and the result
70+
as an xarray repr; a mismatch raises `AssertionError` and exits non-zero. Cases
71+
that need data or credentials you don't have print `⏭ SKIPPED` and exit 0.
8072

81-
Shared helpers (timing, peak memory, the `assert_allclose` wrapper, SQL echo)
82-
live in [`_harness.py`](_harness.py).
73+
Shared helpers timing, peak memory, the result check and its printout, SQL
74+
echo — live in [`_harness.py`](_harness.py).

docs/geospatial.md

Lines changed: 15 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -11,14 +11,11 @@ The array paradigm (NumPy, Xarray, Dask) is a wonderful *interface* for these
1111
operations. But it is not the only one, and for a large and growing audience —
1212
the people fluent in SQL rather than in `apply_ufunc` and rechunking — it is not
1313
the most accessible one. [`xarray-sql`](../README.md) lets you pose these
14-
questions in SQL and answers them with a real query engine (DataFusion),
15-
complete with partition pruning and projection pushdown. The datasets are
16-
registered *lazily* — nothing is read or column-selected up front; each query
17-
pulls only the variable and the partitions it needs, and value filters are
18-
passed as bound **query parameters** rather than formatted into the SQL string.
19-
And because a gridded result is still gridded data, every query here round-trips
20-
its answer straight back to an `xarray.Dataset` (via `to_dataset`) — SQL in, an
21-
array out, ready to plot or save.
14+
questions in SQL and answers them with a real query engine (DataFusion). The
15+
datasets are opened *lazily*, so a query against the whole archive reads only the
16+
variable and the slice it actually needs. And because a gridded result is still
17+
gridded data, every query here round-trips its answer straight back to an
18+
`xarray.Dataset` — SQL in, an array out, ready to plot or save.
2219

2320
This page makes the argument case by case. Every claim below is backed by a
2421
runnable script in [`benchmarks/geospatial/`](../benchmarks/geospatial/) that
@@ -131,10 +128,10 @@ JOIN era5 e
131128
GROUP BY f.model, f.prediction_timedelta
132129
```
133130

134-
Both models are stacked along a `model` dimension into one forecast table, so
135-
the query scores them together by grouping on a `model` *column* — no table name
136-
formatted into the SQL. The entire evaluation — temporal alignment across three
137-
time axes, spatial matching, and the score — is one JOIN and one aggregate.
131+
Both models are stacked along a `model` dimension into one forecast table, so a
132+
single query scores them together, grouped by the `model` column. The entire
133+
evaluation — temporal alignment across three time axes, spatial matching, and the
134+
score — is one JOIN and one aggregate.
138135
[`05_forecast_skill.py`](../benchmarks/geospatial/05_forecast_skill.py) runs it
139136
for both models, matches an xarray reference, and reproduces the published result
140137
that GraphCast edges out Pangu at every lead — the classic "error grows with
@@ -168,8 +165,8 @@ literal: the raster is the full ERA5 archive (the `WHERE` prunes it to a day),
168165
the regions are a second SQL table, and the spatial relationship is an ordinary
169166
`BETWEEN`. See [`06_zonal_vector.py`](../benchmarks/geospatial/06_zonal_vector.py)
170167
— it reports e.g. Sahara 33 °C vs Greenland −8 °C for a June day. (Rectangular
171-
regions keep the demo dependency-free; arbitrary polygons are the natural next
172-
step, via a point-in-polygon UDF — see below.)
168+
regions keep this simple; arbitrary polygons would follow the same shape, with a
169+
point-in-polygon test in the join.)
173170

174171
## 6. The hard cases: where a UDF fits, and where it doesn't
175172

@@ -192,10 +189,10 @@ this against **Earth Engine itself**: it opens a UTM grid through
192189
[Xee](https://github.com/google/Xee) carrying `ee.Image.pixelLonLat()`, so EE's
193190
own geodesy engine reports the true lon/lat of every pixel — an *independent*
194191
reprojection reference, not PROJ-vs-PROJ. The SQL UDF and EE agree to sub-metre
195-
precision. Two caveats, both documented in the script: PROJ's context is not
196-
thread-safe (so the UDF returns both coordinates from *one* call and runs on a
197-
single partition), and reprojection moves coordinates without resampling onto a
198-
grid — which is the next operation.
192+
precision. The script flags one practical gotcha (PROJ is not thread-safe, so the
193+
UDF runs serially), but the caveat that matters here is conceptual: reprojection
194+
moves the coordinates without resampling the data onto a new grid — and *that* is
195+
the next operation.
199196

200197
**Regridding is not** row-independent: each output cell is a weighted blend of
201198
several input cells. That is a *many-to-many* relationship — and a many-to-many

0 commit comments

Comments
 (0)