ENH: Expose target_crs in georeference accessor

docs/history.md

@@ -2,6 +2,9 @@
 
 ## 0.12.0 (2026-04-21)
 
+## Development
+
+* ENH: Expose ``target_crs`` on ``.xradar.georeference()`` accessor for DataArray, Dataset and DataTree, enabling reprojection via e.g. ``radar.xradar.georeference(target_crs=4326)`` ({issue}`243`) by [@syedhamidali](https://github.com/syedhamidali)
 * MNT: Unpin xarray, require ``xarray >= 2026.4.0`` in ``requirements.txt``, ``environment.yml``, ``ci/unittests.yml``, and ``ci/notebooktests.yml`` by [@aladinor](https://github.com/aladinor)
 * MNT: Clarify contributor, team-member, and maintainer roles in the contributing guide, including the pathway to greater project involvement ({issue}`341`) by [@kmuehlbauer](https://github.com/kmuehlbauer), ({pull}`354`) by [@syedhamidali](https://github.com/syedhamidali)
 * FIX: ``open_nexradlevel2_datatree`` keeps sweeps with interior sweep-index gaps — derive sweep names from actual indices in ``nex.data`` instead of ``range(len(...))`` so upstream-dropped interior cuts (e.g. ``[0..9, 11]``) no longer raise ``KeyError`` ({issue}`361`, {pull}`362`) by [@aladinor](https://github.com/aladinor)

docs/notebooks/Georeference_TargetCRS.ipynb

@@ -0,0 +1,290 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "eb6eb2fb",
+   "metadata": {},
+   "source": [
+    "# Reproject Radar Coordinates via the Accessor\n",
+    "\n",
+    "This notebook demonstrates the ``target_crs`` option now exposed on\n",
+    "``.xradar.georeference()`` for ``xarray.Dataset`` and ``xarray.DataTree``.\n",
+    "\n",
+    "It lets you reproject radar ``x, y`` coordinates from the radar-native\n",
+    "Azimuthal Equidistant (AEQD) projection into any ``pyproj``-compatible CRS\n",
+    "in a single accessor call:\n",
+    "\n",
+    "```python\n",
+    "radar = radar.xradar.georeference(target_crs=4326)\n",
+    "```"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "8fe22ab0",
+   "metadata": {},
+   "source": [
+    "## Imports"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "ea52d5af",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import cmweather  # noqa\n",
+    "import matplotlib.pyplot as plt\n",
+    "from open_radar_data import DATASETS\n",
+    "\n",
+    "import xradar as xd"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "de72bd74",
+   "metadata": {},
+   "source": [
+    "## Read a sample radar file"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "ddec4345",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "filename = DATASETS.fetch(\"cfrad.20080604_002217_000_SPOL_v36_SUR.nc\")\n",
+    "radar = xd.io.open_cfradial1_datatree(filename, first_dim=\"auto\")\n",
+    "radar"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "e87b83a0",
+   "metadata": {},
+   "source": [
+    "## 1. Default georeferencing (AEQD, meters)\n",
+    "\n",
+    "Without ``target_crs``, the accessor adds ``x``, ``y``, ``z`` in the\n",
+    "radar-native AEQD projection (distances in meters from the radar site)."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "61077c76",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "radar_aeqd = radar.xradar.georeference()\n",
+    "sweep = radar_aeqd[\"sweep_0\"].to_dataset(inherit=\"all_coords\")\n",
+    "crs = sweep.xradar.get_crs()\n",
+    "print(\n",
+    "    \"projection:\",\n",
+    "    crs.coordinate_operation.method_name if crs.coordinate_operation else crs.name,\n",
+    ")\n",
+    "print(\"is_projected / is_geographic:\", crs.is_projected, \"/\", crs.is_geographic)\n",
+    "print(\"x units:\", sweep.x.attrs.get(\"units\"))\n",
+    "print(\"x range:\", float(sweep.x.min()), float(sweep.x.max()))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "2ce9c5d2",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "fig, ax = plt.subplots(figsize=(7, 6))\n",
+    "radar_aeqd[\"sweep_0\"][\"DBZ\"].plot(x=\"x\", y=\"y\", cmap=\"ChaseSpectral\", ax=ax)\n",
+    "ax.set_title(\"AEQD (default) — x, y in meters from radar\")\n",
+    "ax.set_aspect(\"equal\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "103f0bf4",
+   "metadata": {},
+   "source": [
+    "## 2. Reproject to geographic lon/lat (EPSG:4326)\n",
+    "\n",
+    "Pass ``target_crs=4326`` to get ``x`` as longitude and ``y`` as latitude.\n",
+    "Attributes are updated automatically (``standard_name``, ``units``)."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "d2ea627c",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "radar_geo = radar.xradar.georeference(target_crs=4326)\n",
+    "sweep = radar_geo[\"sweep_0\"].to_dataset(inherit=\"all_coords\")\n",
+    "crs = sweep.xradar.get_crs()\n",
+    "print(\"CRS name:\", crs.name, \"| EPSG:\", crs.to_epsg())\n",
+    "print(\"x attrs:\", dict(sweep.x.attrs))\n",
+    "print(\"y attrs:\", dict(sweep.y.attrs))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "a6ce55e7",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "fig, ax = plt.subplots(figsize=(8, 6))\n",
+    "radar_geo[\"sweep_0\"][\"DBZ\"].plot(x=\"x\", y=\"y\", cmap=\"ChaseSpectral\", ax=ax)\n",
+    "ax.set_xlabel(\"Longitude\")\n",
+    "ax.set_ylabel(\"Latitude\")\n",
+    "ax.set_title(\"EPSG:4326 — x=lon, y=lat\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "9adc842c",
+   "metadata": {},
+   "source": [
+    "## 3. Reproject to a projected CRS (UTM, Web Mercator)\n",
+    "\n",
+    "Any ``pyproj``-compatible CRS input works: EPSG code, WKT string, or a\n",
+    "``pyproj.CRS`` instance."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "b68b43c6",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Web Mercator (meters)\n",
+    "radar_wm = radar.xradar.georeference(target_crs=\"EPSG:3857\")\n",
+    "sweep_wm = radar_wm[\"sweep_0\"].to_dataset(inherit=\"all_coords\")\n",
+    "crs = sweep_wm.xradar.get_crs()\n",
+    "print(\"CRS:\", crs.name, \"| EPSG:\", crs.to_epsg())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "fa85ef99",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import pyproj\n",
+    "\n",
+    "utm_crs = pyproj.CRS(\"EPSG:32633\")  # UTM zone 33N (tweak for your radar site)\n",
+    "radar_utm = radar.xradar.georeference(target_crs=utm_crs)\n",
+    "sweep_utm = radar_utm[\"sweep_0\"].to_dataset(inherit=\"all_coords\")\n",
+    "crs = sweep_utm.xradar.get_crs()\n",
+    "print(\"CRS:\", crs.name, \"| EPSG:\", crs.to_epsg())"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "342b18be",
+   "metadata": {},
+   "source": [
+    "## 4. Plot on a cartopy map\n",
+    "\n",
+    "Once data is in a known CRS, hand it off to cartopy for map plotting."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "316b50c8",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import cartopy.crs as ccrs\n",
+    "\n",
+    "radar_geo = radar.xradar.georeference(target_crs=4326)\n",
+    "\n",
+    "fig = plt.figure(figsize=(9, 8))\n",
+    "ax = fig.add_subplot(111, projection=ccrs.PlateCarree())\n",
+    "radar_geo[\"sweep_0\"][\"DBZ\"].plot(\n",
+    "    x=\"x\",\n",
+    "    y=\"y\",\n",
+    "    cmap=\"ChaseSpectral\",\n",
+    "    transform=ccrs.PlateCarree(),\n",
+    "    ax=ax,\n",
+    "    cbar_kwargs=dict(pad=0.05, shrink=0.7),\n",
+    ")\n",
+    "ax.coastlines()\n",
+    "ax.gridlines(draw_labels=True)\n",
+    "ax.set_title(\"Georeferenced PPI on PlateCarree map\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "d49ee5c1",
+   "metadata": {},
+   "source": [
+    "## 5. Works on a single Dataset too\n",
+    "\n",
+    "The accessor is available on ``Dataset`` and ``DataArray`` as well, not just\n",
+    "``DataTree``."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "014505ea",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "ds = radar[\"sweep_0\"].to_dataset(inherit=\"all_coords\")\n",
+    "ds_geo = ds.xradar.georeference(target_crs=4326)\n",
+    "ds_geo[[\"DBZ\"]]"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "7852a0d2",
+   "metadata": {},
+   "source": [
+    "## Notes\n",
+    "\n",
+    "* ``target_crs`` accepts anything ``pyproj.CRS(...)`` accepts: int\n",
+    "  [EPSG codes](https://epsg.io/), ``\"EPSG:xxxx\"`` strings, WKT, or\n",
+    "  ``pyproj.CRS`` instances.\n",
+    "* The ``spatial_ref`` / ``crs_wkt`` coordinate is updated to reflect the\n",
+    "  target CRS — ``sweep.xradar.get_crs()`` will return it.\n",
+    "* Currently only ``x`` and ``y`` are reprojected; ``z`` stays as AEQD\n",
+    "  altitude above the radar site.\n",
+    "* If you want **separate** ``lon`` / ``lat`` / ``alt`` coordinates (without\n",
+    "  overwriting ``x, y``), see the [Assign_GeoCoords](Assign_GeoCoords.md) notebook."
+   ]
+  }
+ ],
+ "metadata": {
+  "jupytext": {
+   "main_language": "python"
+  },
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "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.13.12"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}

docs/usage.md

@@ -44,6 +44,7 @@ notebooks/UF
 :caption: Examples
 
 notebooks/Assign_GeoCoords
+notebooks/Georeference_TargetCRS
 notebooks/CfRadial1_Export
 notebooks/Read-plot-Sigmet-data-from-AWS
 notebooks/plot-ppi

tests/test_accessors.py

@@ -51,6 +51,27 @@ def test_georeference_datatree():
     )
 
 
+def test_georeference_target_crs_dataset():
+    radar = xd.model.create_sweep_dataset()
+    geo = radar.xradar.georeference(target_crs=4326)
+    crs = geo.xradar.get_crs()
+    assert crs.is_geographic
+    assert geo.x.attrs["standard_name"] == "longitude"
+    assert geo.y.attrs["standard_name"] == "latitude"
+    assert_almost_equal(geo.x.values[0, 0], radar.longitude.values, decimal=3)
+    assert_almost_equal(geo.y.values[0, 0], radar.latitude.values, decimal=3)
+
+
+def test_georeference_target_crs_datatree():
+    radar = xd.model.create_sweep_dataset()
+    tree = xr.DataTree.from_dict({"sweep_0": radar})
+    geo = tree.xradar.georeference(target_crs=4326)["sweep_0"]
+    assert geo["x"].attrs["standard_name"] == "longitude"
+    assert geo["y"].attrs["standard_name"] == "latitude"
+    assert_almost_equal(geo["x"].values[0, 0], radar.longitude.values, decimal=3)
+    assert_almost_equal(geo["y"].values[0, 0], radar.latitude.values, decimal=3)
+
+
 def test_crs_dataarray():
     """Test the add_crs and get_crs methods on a DataArray."""
     # Create a sample DataArray with radar data