fix: embed EPSG in Zarr attrs for TiTiler CRS detection

TiTiler defaults to EPSG:4326 when ds.rio.crs is None, causing validation
errors on UTM coordinates. Extract helper to embed proj:epsg in attrs so
rioxarray reconstructs CRS after Zarr I/O.

Author: wietzesuijker

src/eopf_geozarr/conversion/geozarr.py

@@ -352,6 +352,13 @@ def iterative_copy(
     return dt_result if isinstance(dt_result, xr.DataTree) else xr.DataTree(dt_result)
 
 
+def _embed_crs_in_attrs(ds: xr.Dataset) -> None:
+    """Embed EPSG code in dataset attrs for Zarr persistence (modifies in-place)."""
+    if epsg_code := ds.rio.crs.to_epsg():
+        ds.attrs["proj:epsg"] = epsg_code
+        ds.attrs["proj:code"] = f"EPSG:{epsg_code}"
+
+
 def prepare_dataset_with_crs_info(
     ds: xr.Dataset, reference_crs: str | None = None
 ) -> xr.Dataset:
@@ -380,6 +387,7 @@ def prepare_dataset_with_crs_info(
         print(f"  Adding CRS information: {reference_crs}")
         ds = ds.rio.write_crs(reference_crs)
         ds.attrs["grid_mapping"] = "spatial_ref"
+        _embed_crs_in_attrs(ds)
 
         # Ensure spatial_ref variable has proper attributes
         if "spatial_ref" in ds:
@@ -1003,9 +1011,10 @@ def create_overview_dataset_all_vars(
     # Create overview dataset
     overview_ds = xr.Dataset(overview_data_vars, coords=overview_coords)
 
-    # Set CRS using rioxarray first
+    # Set CRS using rioxarray
     overview_ds.rio.write_crs(native_crs, inplace=True)
     overview_ds.attrs["grid_mapping"] = grid_mapping_var_name
+    _embed_crs_in_attrs(overview_ds)
 
     # Add grid_mapping variable after setting CRS
     # TODO: refactor? grid mapping attributes and variables are handled

src/eopf_geozarr/tests/test_crs_metadata.py

@@ -0,0 +1,53 @@
+"""Tests for CRS metadata embedding in Zarr attributes."""
+
+import numpy as np
+import pytest
+import xarray as xr
+
+from eopf_geozarr.conversion.geozarr import prepare_dataset_with_crs_info
+
+EPSG_UTM32N = 32632
+CRS_UTM32N = f"EPSG:{EPSG_UTM32N}"
+
+
+@pytest.fixture
+def utm_dataset():
+    """Minimal UTM Zone 32N dataset."""
+    return xr.Dataset(
+        {"data": (["y", "x"], np.random.rand(10000, 11000))},
+        coords={"x": np.arange(0, 110000, 10), "y": np.arange(5300000, 5200000, -10)},
+    )
+
+
+@pytest.fixture
+def zarr_roundtrip(utm_dataset, tmp_path):
+    """CRS-prepared dataset after Zarr write/read cycle."""
+    ds = prepare_dataset_with_crs_info(utm_dataset, reference_crs=CRS_UTM32N)
+    path = tmp_path / "test.zarr"
+    ds.drop_encoding().to_zarr(path, mode="w", zarr_format=3, consolidated=True)
+    return xr.open_zarr(path, zarr_format=3, consolidated=True)
+
+
+def test_crs_attrs_embedded(utm_dataset):
+    """EPSG metadata written to dataset attrs."""
+    ds = prepare_dataset_with_crs_info(utm_dataset, reference_crs=CRS_UTM32N)
+    assert ds.attrs["proj:epsg"] == EPSG_UTM32N
+    assert ds.attrs["proj:code"] == CRS_UTM32N
+
+
+def test_crs_persists_through_zarr(zarr_roundtrip):
+    """CRS metadata survives Zarr write/read."""
+    assert zarr_roundtrip.attrs["proj:epsg"] == EPSG_UTM32N
+    assert zarr_roundtrip.attrs["proj:code"] == CRS_UTM32N
+
+
+def test_rioxarray_reads_crs(zarr_roundtrip):
+    """rioxarray reads CRS from Zarr metadata."""
+    assert zarr_roundtrip.rio.crs.to_epsg() == EPSG_UTM32N
+
+
+def test_utm_bounds_preserved(zarr_roundtrip):
+    """Bounds remain UTM meters (not lat/lon degrees)."""
+    left, bottom, right, top = zarr_roundtrip.rio.bounds()
+    assert left < 200000 and right < 200000  # easting in valid UTM range
+    assert bottom > 5000000 and top > 5000000  # northing in valid UTM range