notebooks/geocat-comp/02-more_advanced_functionality.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "# A More Complicated Example: Calculating Climatological Averages"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "The previous example showed one of the most basic GeoCAT-comp functions. This example showcases how some functions in GeoCAT-comp have a variety of options that can be selected using keyword arguments. The particular function we'll be looking at is `climatology_average`.\n",
    "\n",
    "The datafile we're using contains surface temperature data at one location generated by numerous weather models in a model ensemble. The temporal resolution is 6 hours from January 1st 1990 to January 1st 1996. We want to take the data from each model and summarize it. We can do this in a few different ways with `climatology_average`. We can find the daily, monthly, or seasonal climatological averages.\n",
    "\n",
    "## Dependencies\n",
    "- geocat-comp\n",
    "- geocat-datafiles\n",
    "- xarray\n",
    "- nc-time-axis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "from geocat.comp import climatology_average\n",
    "\n",
    "import geocat.datafiles as gdf\n",
    "import xarray as xr\n",
    "import nc_time_axis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## Retrieving sample data from GeoCAT-datafiles\n",
    "[GeoCAT-datafiles](https://github.com/ncar/geocat-datafiles) is a software package that makes it simple to download datafiles from the pacakge's GitHub repository. The repository holds data files of different types such as NetCDF and GRIB. They are used for GeoCAT tutorials, demonstrations, and test cases. You won't need to use GeoCAT-datafiles when working with your own data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "nc_file = gdf.get('netcdf_files/atm.20C.hourly6-1990-1995-TS.nc')\n",
    "nc_file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## Parsing the data with xarray\n",
    "GeoCAT-comp is designed to work with [`xarray`](https://xarray.dev/) and [`numpy`](https://numpy.org/) data structures. Here we will be using `xarray` to import the netCDF data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "ds = xr.open_dataset(nc_file)\n",
    "ds"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "We can see here that the datafile has surface temperature data from 40 model ensemble members. Let's simplify things and just look at one model for now. We use `xarray`'s `DataArray.isel` method to pick out the data for the first member_id at index 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "ds = ds.isel(member_id=0)  # select one model from the ensemble\n",
    "temp = ds.TS  # surface temperature data\n",
    "temp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Now that we have the temperature data for the first ensemble member, let's plot it. `xarray` has built in methods to quickly plot data. The generated plot isn't publication quality, but it's a quick and dirty way to see what your dataset contains."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "temp.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "This plot shows that the temperature doesn't change much throughout the year. This makes sense since the data is for the location 0.4712N 180E, which is in the Pacific Ocean just north of the equator. There is some variation though, and it'd be interesting to know what the average temperature is in each month. We can calculate this using `climatology_average`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## Using climatology_average\n",
    "We can see in the [documentation](https://geocat-comp.readthedocs.io/en/latest/user_api/generated/geocat.comp.climatologies.climatology_average.html) that `climatology_average` calculates long term hourly, daily, monthly, or seasonal averages across all years in a given dataset. This is a good function to quickly summarize our data, so let's use it. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "<iframe src=\"https://geocat-comp.readthedocs.io/en/latest/user_api/generated/geocat.comp.climatologies.climatology_average.html\" title=\"geocat-comp climatology_average documentation\" width=100% height=600></iframe>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "### Related NCL Functions\n",
    "Notice at the bottom of the documentation page there is a box with links to related NCL and GeoCAT functions. This indicates all of the NCL functions that `climatology_average` is replicating through it's keyword arguments. We want something similar to the `clmMon` family of functions to summarize our data using monthly climate averages. We can see that our chosen GeoCAT function replicates that NCL capability.\n",
    "\n",
    "<img src=\"../images/geocat-comp/climatology_average_note.png\" width=\"100%\" alt=\"A note on the climatology_average documentation page about related NCL functions clmDayHourTLL, clmDauHourTLLL, clmDayTLL, clmDayTLLL, clmMonLLLT, clmMonLLT, clmMonTLL, clmMonTLLL, month_to_season\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Let's summarize our data using monthly climate averages. Using the keyword argument `freq`, we can specify the time span for the average. The accepted options are `hour`, `day`, `month`, and `season`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "monthly_temp = climatology_average(temp, freq='month')\n",
    "monthly_temp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Notice how much smaller the time dimension is. Before, the time dimension had a length of 8761 points, one for every 6 hour interval over 5 years. Now, the time dimension has a length of 12, one for every month of the year. Notice that the value of the time coordinate is for the middle of each month in 1993. Why 1993? It's approximatly the middle year in the dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "monthly_temp.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "When we plot this, the data is more granular and has even less variation in temperature from month to month. But we can now see which months are coldest and warmest on average for this location.\n",
    "\n",
    "## Adjusting Parameters\n",
    "Monthly climate averages are great, but what else can `climatology_average` do? We can also calculate hourly, daily, and seasonal climate averages. For this dataset, we don't have a fine enough temporal resolution to compute an hourly average, but we can compute daily averages. The only change to the function call is passing in `day` instead of `month`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "daily_temp = climatology_average(temp, freq='day')\n",
    "daily_temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "daily_temp.plot()"
   ]
  }
 ],
 "metadata": {
  "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.10.6"
  },
  "vscode": {
   "interpreter": {
    "hash": "7961b2108ee7d1ae8b8cda2085110605487fed8e1d2fc7de0f719df8d0b8bd85"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}