{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "f15f38bb-57f7-4104-8824-9feec9deee47",
   "metadata": {},
   "source": [
    "# Create a Plan Position Indicator (PPI) Plot\n",
    "A Plan Position Indicator (PPI) plot is a common plot requested by radar scientists. Let's show how to create this plot using `xradar`"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "91562ade-05d4-4b1e-ab7b-720b04641c61",
   "metadata": {},
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ce9365db-30f0-4df1-944f-b3c739c0ef99",
   "metadata": {},
   "outputs": [],
   "source": [
    "import xarray as xr\n",
    "import xradar as xd\n",
    "import cmweather\n",
    "from open_radar_data import DATASETS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "02f07bb2-d1ad-40c2-9a32-3795de327003",
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import pyproj\n",
    "import cartopy"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "100584b5-d600-4f93-980b-6c1dde9335a8",
   "metadata": {},
   "source": [
    "## Read in some data\n",
    "\n",
    "Fetching CfRadial1 radar data file from [open-radar-data](https://github.com/openradar/open-radar-data) repository."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cf955c1f-fa3b-4545-b90a-7afcd9b11b31",
   "metadata": {},
   "outputs": [],
   "source": [
    "filename = DATASETS.fetch(\"cfrad.20080604_002217_000_SPOL_v36_SUR.nc\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b88a784c-8a18-4a75-973e-5723acdcb4b0",
   "metadata": {},
   "source": [
    "Read the data using the `cfradial1` engine"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e97d5717-e6aa-4155-abe5-0960d846ed1c",
   "metadata": {},
   "outputs": [],
   "source": [
    "radar = xd.io.open_cfradial1_datatree(filename, first_dim=\"auto\")\n",
    "display(radar)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7eb28292-d412-4b9f-b506-b983f750ae18",
   "metadata": {},
   "source": [
    "## Add georeferencing\n",
    "We can use the georeference function, or the accessor to add our georeference information!"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "998f378e-7bed-473c-bab9-fa2fffb9513c",
   "metadata": {},
   "source": [
    "### Georeference Accessor\n",
    "If you prefer the accessor (`.xradar.georefence()`), this is how you would add georeference information to your radar object."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eaa2ee4c-98f2-41d2-b16e-9d05ad8608da",
   "metadata": {},
   "outputs": [],
   "source": [
    "radar = radar.xradar.georeference()\n",
    "display(radar)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "53683020-5503-4be6-8eb5-2fc40c340ef3",
   "metadata": {},
   "source": [
    "Please observe, that the additional coordinates `x`, `y`, `z` have been added to the dataset. This will also add `spatial_ref` CRS information on the used Azimuthal Equidistant Projection."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a2579869-1e4f-48d2-93eb-d57084d78fe6",
   "metadata": {},
   "outputs": [],
   "source": [
    "radar[\"sweep_0\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cb6af121-0c94-4508-84bb-7d089e5259af",
   "metadata": {},
   "source": [
    "### Use the Function\n",
    "We can also use the function `xd.geoference.get_x_y_z_tree` function if you prefer that method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9c51a29c-77ba-4b0a-a01b-f5ae0fd18681",
   "metadata": {},
   "outputs": [],
   "source": [
    "radar = xd.georeference.get_x_y_z_tree(radar)\n",
    "display(radar[\"sweep_0\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "72231f7d-2b93-4f15-a8c7-66e4f64f7f8f",
   "metadata": {},
   "source": [
    "## Plot our Data\n",
    "\n",
    "### Plot simple PPI\n",
    "\n",
    "Now, let's create our PPI plot! We just use the newly created 2D-coordinates `x` and `y` to create a meshplot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c2f4cfaa-eba9-4b9d-8d84-7b60027ae0fb",
   "metadata": {},
   "outputs": [],
   "source": [
    "radar[\"sweep_0\"][\"DBZ\"].plot(x=\"x\", y=\"y\", cmap=\"ChaseSpectral\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "505e3ad9-ff77-4d7c-b09e-28b87c3328b5",
   "metadata": {},
   "source": [
    "### Plot PPI on map with cartopy\n",
    "\n",
    "If you have `cartopy` installed, you can easily plot on maps. We first have to extract the CRS from the dataset and to wrap it in a `cartopy.crs.Projection`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2b67992f-3ba5-4579-a72e-7303623ccad4",
   "metadata": {},
   "outputs": [],
   "source": [
    "proj_crs = xd.georeference.get_crs(radar[\"sweep_0\"].ds)\n",
    "cart_crs = cartopy.crs.Projection(proj_crs)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9a3d9c1e-0eac-4fc2-bf0e-c8feffc520c3",
   "metadata": {},
   "source": [
    "Second, we create a matplotlib GeoAxes and a nice map."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b9b7680c-9ca2-4b93-8d19-69b9bdb64ae3",
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(10, 10))\n",
    "ax = fig.add_subplot(111, projection=cartopy.crs.PlateCarree())\n",
    "radar[\"sweep_0\"][\"DBZ\"].plot(\n",
    "    x=\"x\",\n",
    "    y=\"y\",\n",
    "    cmap=\"ChaseSpectral\",\n",
    "    transform=cart_crs,\n",
    "    cbar_kwargs=dict(pad=0.075, shrink=0.75),\n",
    ")\n",
    "ax.coastlines()\n",
    "ax.gridlines(draw_labels=True)"
   ]
  }
 ],
 "metadata": {
  "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.9.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}