Sentinel-2 EOPF GeoZarr

docs/_quarto.yml

@@ -21,6 +21,8 @@ website:
         text: Web-optimized Zarr (WOZ)
       - href: woz-diagram.qmd
         text: Visual WOZ comparison
+      - href: geozarr-eopf-specification.qmd
+        text: GeoZarr EOPF specification
       - section: Presentations
         contents:
         - href: slides/2025-02.qmd
@@ -41,6 +43,8 @@ website:
             text: WebMercatorQuad overviews (Zarr V2)
           - href: examples/04_multiscales_as_WebMercatorQuad_ZarrV3.ipynb
             text: WebMercatorQuad overviews (Zarr V3)
+          - href: examples/06_multiscales_EOPF_ZarrV3.ipynb
+            text: GeoZarr for EOPF
 
 format:
 

docs/examples/06_embedded_STAC_block.ipynb

@@ -0,0 +1,142 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "201df5fa",
+   "metadata": {},
+   "source": [
+    "# Sentinel-2 RGB to GeoZarr with embedded STAC\n",
+    "\n",
+    "This example demonstrates how to:\n",
+    "\n",
+    "1. Query Earth-Search for a low-cloud Sentinel-2 L2A scene over Vienna.\n",
+    "2. Stream the 10 m red/green/blue COGs into xarray/rioxarray.\n",
+    "3. Write a consolidated, 512 × 512-chunk GeoZarr dataset.\n",
+    "4. Embed a minimal STAC Item (`proj:code`, `bbox`, `geometry`, `gsd`, etc.) into `root/.zattrs`.\n",
+    "5. Re‐consolidate metadata so that `xr.open_zarr(consolidated=True)` exposes the `stac` block.\n",
+    "6. Coarsen the RGB by 4 × 4 and display a quick preview."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "be2dcbf1",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "\"\"\"Create a STAC-aware GeoZarr RGB tile of a Sentinel-2 scene over Vienna.\"\"\"\n",
+    "\n",
+    "import shutil\n",
+    "from datetime import date\n",
+    "from pathlib import Path\n",
+    "\n",
+    "import jsonschema\n",
+    "import matplotlib.pyplot as plt\n",
+    "import pystac_client\n",
+    "import rioxarray as rxr\n",
+    "import xarray as xr\n",
+    "import zarr\n",
+    "\n",
+    "# 1. Query Earth-Search for a recent, low-cloud Sentinel-2 L2A scene \n",
+    "API = \"https://earth-search.aws.element84.com/v1\"\n",
+    "coll = \"sentinel-2-l2a\"\n",
+    "bbox = [16.20, 48.10, 16.45, 48.30]  # Vienna\n",
+    "today, last_year = date.today(), date.today().replace(year=date.today().year - 1)\n",
+    "daterange = f\"{last_year:%Y-%m-%d}/{today:%Y-%m-%d}\"\n",
+    "\n",
+    "item = next(\n",
+    "    pystac_client.Client.open(API)\n",
+    "    .search(collections=[coll], bbox=bbox, datetime=daterange,\n",
+    "            query={\"eo:cloud_cover\": {\"lt\": 5}}, limit=1)\n",
+    "    .items(),\n",
+    "    None,\n",
+    ")\n",
+    "assert item, \"No Sentinel-2 scene found\"\n",
+    "print(\"Scene:\", item.id, \"cloud\", item.properties[\"eo:cloud_cover\"])\n",
+    "\n",
+    "# 2. Stack RGB bands lazily \n",
+    "bands = [\"red\", \"green\", \"blue\"]\n",
+    "rgb = xr.concat(\n",
+    "    [rxr.open_rasterio(item.assets[b].href,\n",
+    "                       chunks={\"band\": 1, \"x\": 2048, \"y\": 2048},\n",
+    "                       masked=True).assign_coords(band=[b])\n",
+    "     for b in bands],\n",
+    "    dim=\"band\",\n",
+    ")\n",
+    "rgb.name = \"radiance\"\n",
+    "rgb = rgb.rio.write_crs(item.properties[\"proj:code\"])\n",
+    "rgb.attrs[\"transform\"] = list(rgb.rio.transform())\n",
+    "\n",
+    "# 3. Write as consolidated GeoZarr (Dataset) \n",
+    "store = Path(f\"{coll}_{'_'.join(bands)}_{item.id}.zarr\")\n",
+    "if store.exists():\n",
+    "    shutil.rmtree(store)\n",
+    "\n",
+    "(radiance_ds := rgb.drop_vars([\"x\", \"y\"]).to_dataset()) \\\n",
+    "    .chunk({\"y\": 512, \"x\": 512}) \\\n",
+    "    .to_zarr(store, mode=\"w\", consolidated=True)\n",
+    "\n",
+    "# 4. Embed a minimal STAC Item in .zattrs, then re-consolidate metadata \n",
+    "gsd = min(item.assets[b].to_dict().get(\"gsd\", 10) for b in bands)\n",
+    "mini = {\n",
+    "    \"type\": \"Item\",\n",
+    "    \"stac_version\": \"1.0.0\",\n",
+    "    \"id\": item.id,\n",
+    "    \"bbox\": item.bbox,\n",
+    "    \"geometry\": item.geometry,\n",
+    "    \"properties\": {\n",
+    "        \"datetime\": item.properties[\"datetime\"],\n",
+    "        \"proj:code\": item.properties[\"proj:code\"],\n",
+    "        \"proj:bbox\": item.bbox,\n",
+    "        \"platform\": item.properties[\"platform\"],\n",
+    "        \"instruments\": item.properties[\"instruments\"],\n",
+    "        \"eo:cloud_cover\": item.properties[\"eo:cloud_cover\"],\n",
+    "        \"gsd\": gsd,\n",
+    "    },\n",
+    "    \"assets\": {\n",
+    "        \"data\": {\"href\": store.name,\n",
+    "                 \"type\": \"application/x-zarr\",\n",
+    "                 \"roles\": [\"data\"]}\n",
+    "    },\n",
+    "    \"license\": item.properties.get(\"license\", \"proprietary\"),\n",
+    "}\n",
+    "jsonschema.validate(mini, {\"type\": \"object\", \"required\": [\"type\", \"id\", \"assets\"]})\n",
+    "root = zarr.open_group(store, mode=\"a\")\n",
+    "root.attrs[\"stac\"] = mini\n",
+    "zarr.convenience.consolidate_metadata(store)  # update .zmetadata\n",
+    "\n",
+    "# 5. Re-open and show the STAC block is present \n",
+    "ds = xr.open_zarr(store, consolidated=True)\n",
+    "print(\"Embedded STAC:\", ds.attrs[\"stac\"])\n",
+    "\n",
+    "# 6. Quick preview: 4×4 coarsened RGB \n",
+    "preview = ds.radiance.coarsen(y=4, x=4, boundary=\"trim\").mean()\n",
+    "plt.imshow(preview.transpose(\"y\", \"x\", \"band\").astype(\"uint8\").values)\n",
+    "plt.axis(\"off\")\n",
+    "plt.title(\"Coarsened Sentinel-2 RGB\")\n",
+    "plt.show()"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "mm",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.12.10"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}