From 171da5248800e9c2e2dbc0b25d8128ab10d87744 Mon Sep 17 00:00:00 2001 From: Brendan Collins Date: Fri, 26 Jun 2026 20:01:03 -0700 Subject: [PATCH 1/3] Add global-projection templates to from_template (#3549) Add web_mercator (EPSG:3857), equal_earth (EPSG:8857), and wgs84/latlon aliases for the EPSG:4326 world grid. All resolve to real EPSG codes so attrs['crs'] stays an int. Reuses the existing backend dispatch, so all four backends work. --- xrspatial/_template_data.py | 23 ++++++++ xrspatial/templates.py | 35 +++++++----- xrspatial/tests/test_templates.py | 89 ++++++++++++++++++++++++++++--- 3 files changed, 128 insertions(+), 19 deletions(-) diff --git a/xrspatial/_template_data.py b/xrspatial/_template_data.py index bbd2ab4c5..50cc5ef31 100644 --- a/xrspatial/_template_data.py +++ b/xrspatial/_template_data.py @@ -52,6 +52,29 @@ default_resolution=0.5, label='World (WGS84)', lonlat=(-180.0, -85.0, 180.0, 85.0), area_epsg=8857, shape_epsg=3395), + # Global Web Mercator (EPSG:3857). Native bounds are the canonical square + # extent +/-20037508 m, which is +/-180 lon and +/-85.051129 lat (the + # latitude where the Mercator y matches the x half-extent). lonlat stops at + # +/-85.051129 because Mercator y diverges to infinity at the poles. + 'web_mercator': dict( + bounds=(-20037508, -20037508, 20037508, 20037508), crs=3857, + default_resolution=50000, label='World (Web Mercator)', + lonlat=(-180.0, -85.051129, 180.0, 85.051129), + area_epsg=8857, shape_epsg=3395), + # Global equal-area grid (EPSG:8857, Equal Earth). Native bounds are the + # +/-90 world projected into Equal Earth. + 'equal_earth': dict( + bounds=(-17243959, -8392928, 17243959, 8392928), crs=8857, + default_resolution=50000, label='World (Equal Earth)', + lonlat=(-180.0, -90.0, 180.0, 90.0), + area_epsg=8857, shape_epsg=3395), +} + +# Alternate spellings that resolve to a curated region (single source of truth). +# 'wgs84' / 'latlon' are friendly names for the EPSG:4326 'world' grid. +_REGION_ALIASES = { + 'wgs84': 'world', + 'latlon': 'world', } # Equal-area fallback when a template has no curated ``area_epsg`` diff --git a/xrspatial/templates.py b/xrspatial/templates.py index e24371eea..ace337ad8 100644 --- a/xrspatial/templates.py +++ b/xrspatial/templates.py @@ -14,7 +14,8 @@ from xrspatial._template_data import (_CITIES, _CITY_DEFAULT_RESOLUTION, _COUNTRY_BBOXES, _COUNTRY_DEFAULT_RESOLUTION, _EQUAL_AREA_FALLBACK_EPSG, - _REGIONS, _UPS_NORTH_EPSG, _UPS_SOUTH_EPSG) + _REGION_ALIASES, _REGIONS, _UPS_NORTH_EPSG, + _UPS_SOUTH_EPSG) from xrspatial.reproject._crs_utils import _require_pyproj, _resolve_crs from xrspatial.reproject._grid import _edge_samples, _make_output_coords, _transform_boundary @@ -39,11 +40,12 @@ def _resolve(name): f"got {type(name).__name__}" ) - region = _REGIONS.get(name.lower()) + name_l = name.lower() + region = _REGIONS.get(_REGION_ALIASES.get(name_l, name_l)) if region is not None: return dict( bounds=region["bounds"], crs=region["crs"], - default_resolution=region["default_resolution"], key=name.lower(), + default_resolution=region["default_resolution"], key=name_l, lonlat=region["lonlat"], area_epsg=region.get("area_epsg"), shape_epsg=region.get("shape_epsg"), ) @@ -65,7 +67,7 @@ def _resolve(name): lonlat=bbox, area_epsg=None, shape_epsg=None, ) - regions = ", ".join(sorted(_REGIONS)) + regions = ", ".join(sorted(set(_REGIONS) | set(_REGION_ALIASES))) raise ValueError( f"Unknown template {name!r}. Available named regions: {regions}. " f"{len(_CITIES)} world cities are also supported (lowercase name, e.g. " @@ -235,7 +237,7 @@ def list_templates(kind: Optional[str] = None) -> Union[Dict[str, List[str]], Li True """ groups = { - "regions": sorted(_REGIONS), + "regions": sorted(set(_REGIONS) | set(_REGION_ALIASES)), "cities": sorted(_CITIES), "countries": sorted(_COUNTRY_BBOXES), } @@ -270,13 +272,16 @@ def from_template(name: str, ---------- name : str A curated region name (case-insensitive), e.g. ``'conus'``, ``'nyc'``, - ``'europe'``, ``'world'``; a world-city name (case-insensitive), e.g. - ``'london'``, ``'tokyo'``, ``'sao_paulo'``; or an ISO-3166 / GADM alpha-3 - country code, e.g. ``'USA'``, ``'FRA'``, ``'JPN'``. Curated regions and - cities come back in a projected CRS (cities in their UTM zone); country - codes come back in EPSG:4326. Where two cities share a name the larger - keeps the bare name and the others take a ``_`` suffix - (e.g. ``'hyderabad'`` vs ``'hyderabad_pk'``). + ``'europe'``, ``'world'``; a global-projection name, e.g. + ``'web_mercator'`` (EPSG:3857), ``'wgs84'`` / ``'latlon'`` (EPSG:4326, + the same grid as ``'world'``), or ``'equal_earth'`` (EPSG:8857); a + world-city name (case-insensitive), e.g. ``'london'``, ``'tokyo'``, + ``'sao_paulo'``; or an ISO-3166 / GADM alpha-3 country code, e.g. + ``'USA'``, ``'FRA'``, ``'JPN'``. Curated regions and cities come back in + a projected CRS (cities in their UTM zone); country codes come back in + EPSG:4326. Where two cities share a name the larger keeps the bare name + and the others take a ``_`` suffix (e.g. ``'hyderabad'`` vs + ``'hyderabad_pk'``). resolution : float or tuple of float, optional Cell size in the template's CRS units (metres for projected regions, degrees for country codes). A scalar gives square cells; a @@ -330,6 +335,12 @@ def from_template(name: str, >>> agg = from_template("FRA") # France bbox in EPSG:4326 >>> agg.attrs["crs"] 4326 + >>> from_template("web_mercator").attrs["crs"] # global EPSG:3857 + 3857 + >>> from_template("latlon").attrs["crs"] # alias for wgs84 + 4326 + >>> from_template("equal_earth").attrs["crs"] # global equal-area + 8857 >>> from_template("london").attrs["crs"] # greater London, UTM 30N 32630 >>> from_template("FRA", preserve="shape").attrs["crs"] # UTM 30N diff --git a/xrspatial/tests/test_templates.py b/xrspatial/tests/test_templates.py index 9650d31f4..f09ee4163 100644 --- a/xrspatial/tests/test_templates.py +++ b/xrspatial/tests/test_templates.py @@ -6,7 +6,8 @@ import xarray as xr from xrspatial import from_template, list_templates, slope -from xrspatial._template_data import _CITIES, _CITY_DEFAULT_RESOLUTION, _COUNTRY_BBOXES, _REGIONS +from xrspatial._template_data import (_CITIES, _CITY_DEFAULT_RESOLUTION, _COUNTRY_BBOXES, + _REGION_ALIASES, _REGIONS) from xrspatial.tests.general_checks import cuda_and_cupy_available, dask_array_available @@ -144,6 +145,62 @@ def test_world_grid(): assert agg.attrs["crs"] == 4326 +# --------------------------------------------------------------------------- +# global-projection templates +# --------------------------------------------------------------------------- + +@pytest.mark.parametrize( + "name,crs", + [("web_mercator", 3857), ("wgs84", 4326), ("latlon", 4326), + ("equal_earth", 8857)], +) +def test_global_projection_contract(name, crs): + agg = from_template(name) + assert agg.attrs["crs"] == crs + assert agg.dims == ("y", "x") + assert agg.shape[0] > 1 and agg.shape[1] > 1 + assert np.isnan(agg.values).all() + assert agg.dtype == np.float32 + assert agg.name == name + # north-up (descending y), ascending x + assert agg.y.values[0] > agg.y.values[-1] + assert agg.x.values[0] < agg.x.values[-1] + + +@pytest.mark.parametrize("alias", ["wgs84", "latlon"]) +def test_wgs84_latlon_alias_world(alias): + # the aliases resolve to the EPSG:4326 'world' grid (same coords and attrs); + # only the DataArray name reflects the spelling that was asked for + a = from_template(alias) + world = from_template("world") + np.testing.assert_array_equal(a.x.values, world.x.values) + np.testing.assert_array_equal(a.y.values, world.y.values) + assert a.attrs == world.attrs + assert a.name == alias + + +@pytest.mark.parametrize("name", ["web_mercator", "equal_earth", "latlon"]) +def test_global_projection_case_insensitive(name): + a = from_template(name) + b = from_template(name.upper()) + np.testing.assert_array_equal(a.x.values, b.x.values) + assert a.attrs == b.attrs + + +def test_web_mercator_metre_coords_within_bounds(): + agg = from_template("web_mercator") + assert agg.x.attrs["units"] == "m" + assert agg.x.attrs["standard_name"] == "projection_x_coordinate" + left, bottom, right, top = _REGIONS["web_mercator"]["bounds"] + assert left <= agg.x.values.min() and agg.x.values.max() <= right + assert bottom <= agg.y.values.min() and agg.y.values.max() <= top + + +def test_global_resolution_honored_exactly(): + agg = from_template("web_mercator", resolution=100000) + assert agg.attrs["res"] == (100000.0, 100000.0) + + @pytest.mark.parametrize("bad", ["does-not-exist", "ZZZ"]) def test_unknown_name_raises(bad): with pytest.raises(ValueError, match="Unknown template"): @@ -159,18 +216,21 @@ def test_unknown_name_points_to_list_templates(): def test_list_templates_grouped(): names = list_templates() assert set(names) == {"regions", "cities", "countries"} - # each group lists exactly its registry keys, sorted - assert names["regions"] == sorted(_REGIONS) + # each group lists exactly its registry keys, sorted; the regions group also + # advertises the alias spellings (wgs84, latlon) so they are discoverable + assert names["regions"] == sorted(set(_REGIONS) | set(_REGION_ALIASES)) assert names["cities"] == sorted(_CITIES) assert names["countries"] == sorted(_COUNTRY_BBOXES) @pytest.mark.parametrize( - "kind,registry", - [("regions", _REGIONS), ("cities", _CITIES), ("countries", _COUNTRY_BBOXES)], + "kind,expected", + [("regions", sorted(set(_REGIONS) | set(_REGION_ALIASES))), + ("cities", sorted(_CITIES)), + ("countries", sorted(_COUNTRY_BBOXES))], ) -def test_list_templates_kind_filter(kind, registry): - assert list_templates(kind) == sorted(registry) +def test_list_templates_kind_filter(kind, expected): + assert list_templates(kind) == expected def test_list_templates_bad_kind_raises(): @@ -528,6 +588,21 @@ def test_grid_mapping_omitted_for_equal_earth(): assert "Equal Earth" in agg.attrs["crs_wkt"] +@pytest.mark.parametrize( + "name,crs,wkt_marker", + [("web_mercator", 3857, "Pseudo-Mercator"), + ("equal_earth", 8857, "Equal Earth")], +) +def test_global_projection_grid_mapping(name, crs, wkt_marker): + # neither Pseudo-Mercator nor Equal Earth has a CF grid_mapping_name, so the + # key is left off and crs_wkt stands alone (carrying the human CRS name). + agg = from_template(name) + assert agg.attrs["crs"] == crs + assert "grid_mapping_name" not in agg.attrs + assert wkt_marker in agg.attrs["crs_wkt"] + assert agg.x.attrs["units"] == "m" + + def test_cf_attrs_omitted_without_pyproj(monkeypatch): # Without pyproj the default (non-reproject) path stays dependency-free: # the CF grid-mapping keys are left off rather than raising. From 174b1343495a6696185a53e532660aa6ae72cc05 Mon Sep 17 00:00:00 2001 From: Brendan Collins Date: Fri, 26 Jun 2026 20:01:54 -0700 Subject: [PATCH 2/3] Document global-projection templates (#3549) --- README.md | 2 +- docs/source/reference/templates.rst | 4 +++- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/README.md b/README.md index 3765d2eda..d5581e851 100644 --- a/README.md +++ b/README.md @@ -274,7 +274,7 @@ Built-in Numba JIT and CUDA projection kernels bypass pyproj for per-pixel coord | Name | Description | Source | NumPy xr.DataArray | Dask xr.DataArray | CuPy GPU xr.DataArray | Dask GPU xr.DataArray | |:----------:|:------------|:------:|:----------------------:|:--------------------:|:-------------------:|:------:| -| [from_template](xrspatial/templates.py) | Empty study-area grid for a named region (CONUS, NYC, ...), a world city (London, Tokyo, ... in its UTM zone) or country code; `preserve='area'/'shape'` picks an EPSG projection by property; `list_templates()` lists every accepted name | Custom | ✅ | 🔼 | 🔼 | 🔼 | +| [from_template](xrspatial/templates.py) | Empty study-area grid for a named region (CONUS, NYC, ...), a world city (London, Tokyo, ... in its UTM zone), a country code, or a whole-world projection (web_mercator, wgs84/latlon, equal_earth); `preserve='area'/'shape'` picks an EPSG projection by property; `list_templates()` lists every accepted name | Custom | ✅ | 🔼 | 🔼 | 🔼 | ----------- diff --git a/docs/source/reference/templates.rst b/docs/source/reference/templates.rst index b07341698..7bb20e218 100644 --- a/docs/source/reference/templates.rst +++ b/docs/source/reference/templates.rst @@ -9,7 +9,9 @@ a region name, a world-city name, or a country code into a NaN-filled :class:`xarray.DataArray` that follows the xarray-spatial array contract, so it feeds straight into the rest of the library. Cities (national capitals, major regional metros, and recognizable US secondary cities) come back as a metro -bounding box in their UTM zone. +bounding box in their UTM zone. Whole-world canvases are available in a few +projections too: ``'web_mercator'`` (EPSG:3857), ``'wgs84'`` / ``'latlon'`` +(EPSG:4326), and ``'equal_earth'`` (EPSG:8857). Call :func:`~xrspatial.templates.list_templates` to discover every name ``from_template`` accepts (curated regions, world cities, and country codes). From b1719cbf3188a424c11c6cbbbee68ffe3b5baacd Mon Sep 17 00:00:00 2001 From: Brendan Collins Date: Fri, 26 Jun 2026 20:04:10 -0700 Subject: [PATCH 3/3] Address review: test preserve on global templates, exact merc latitude (#3549) --- xrspatial/_template_data.py | 6 +++--- xrspatial/tests/test_templates.py | 11 +++++++++++ 2 files changed, 14 insertions(+), 3 deletions(-) diff --git a/xrspatial/_template_data.py b/xrspatial/_template_data.py index 50cc5ef31..ccf011d3b 100644 --- a/xrspatial/_template_data.py +++ b/xrspatial/_template_data.py @@ -53,13 +53,13 @@ lonlat=(-180.0, -85.0, 180.0, 85.0), area_epsg=8857, shape_epsg=3395), # Global Web Mercator (EPSG:3857). Native bounds are the canonical square - # extent +/-20037508 m, which is +/-180 lon and +/-85.051129 lat (the + # extent +/-20037508 m, which is +/-180 lon and +/-85.0511287798 lat (the # latitude where the Mercator y matches the x half-extent). lonlat stops at - # +/-85.051129 because Mercator y diverges to infinity at the poles. + # that latitude because Mercator y diverges to infinity at the poles. 'web_mercator': dict( bounds=(-20037508, -20037508, 20037508, 20037508), crs=3857, default_resolution=50000, label='World (Web Mercator)', - lonlat=(-180.0, -85.051129, 180.0, 85.051129), + lonlat=(-180.0, -85.0511287798, 180.0, 85.0511287798), area_epsg=8857, shape_epsg=3395), # Global equal-area grid (EPSG:8857, Equal Earth). Native bounds are the # +/-90 world projected into Equal Earth. diff --git a/xrspatial/tests/test_templates.py b/xrspatial/tests/test_templates.py index f09ee4163..8177cc3ae 100644 --- a/xrspatial/tests/test_templates.py +++ b/xrspatial/tests/test_templates.py @@ -578,6 +578,17 @@ def test_cf_grid_mapping_preserve_path(): assert "Conus Albers" in agg.attrs["crs_wkt"] +@pytest.mark.parametrize("name", ["web_mercator", "equal_earth"]) +def test_global_preserve_picks_world_projection(name): + # the global templates carry the same area/shape EPSG hints as 'world', so + # preserve='shape' lands on World Mercator (EPSG:3395) instead of a stray + # UTM zone for the (0, 0) centroid, and preserve='area' lands on Equal Earth + agg_shape = from_template(name, preserve="shape") + assert agg_shape.attrs["crs"] == 3395 + assert _proj(3395) == "merc" + assert from_template(name, preserve="area").attrs["crs"] == 8857 + + def test_grid_mapping_omitted_for_equal_earth(): # Equal Earth (the preserve='area' fallback for the world bbox) has no # CF grid mapping, so grid_mapping_name is left off and crs_wkt stands