{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Transforming Vertical Coordinates\n", "\n", "A common need in the analysis of ocean and atmospheric data is to transform the vertical coordinate from its original coordinate (e.g. depth) to a new coordinate (e.g. density).\n", "Xgcm supports this sort of one-dimensional coordinate transform on `Axis` and `Grid` objects using the `transform` method.\n", "Two algorithms are implemented:\n", "\n", "- _Linear interpolation:_ Linear interpolation is designed to interpolate intensive quantities (e.g. temperature) from one coordinate to another. This method is suitable when the target coordinate is monotonically increasing or decreasing and the data variable is intensive. For example, you want to visualize oxygen on density surfaces from a z-coordinate ocean model.\n", " - _Logarithmic interpolation:_ Logarithmic interpolation (which is linear interpolation done after applying a logarithm to the target coordinate) is also available. This method is suitable when variation of the intensive quantity is best related to the logarithm of the target coordinate, rather than the target coordinate itself. For example, you want to analyze data from a sigma-coordinate atmospheric model on isobaric (constant pressure) surfaces.\n", "- _Conservative remapping:_ This algorithm is designed to conserve extensive quantities (e.g. transport, heat content). It requires knowledge of cell bounds in both the source and target coordinate. It also handles non-monotonic target coordinates.\n", "\n", "On this page, we explain how to use these coordinate transformation capabilities." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from xgcm import Grid\n", "import xarray as xr\n", "import numpy as np\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 1D Toy Data Example\n", "\n", "First we will create a simple, one-dimensional dataset to illustrate how the `transform` function works.\n", "This dataset contains\n", "\n", "- a coordinate called `z`, representing the original depth coordinate\n", "\n", "- a data variable called `theta`, a function of `z`, which we want as our new vertical coordinate\n", "\n", "- a data variable called `phi`, also a function of `z`, which represents the data we want to transform into this new coordinate space\n", "\n", "In an oceanic context `theta` might be density and `phi` might be oxygen.\n", "In an atmospheric context, `theta` might be potential temperature and `phi` might be potential vorticity." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
<xarray.Dataset>\n", "Dimensions: (z: 10)\n", "Coordinates:\n", " * z (z) int64 2 3 4 5 6 7 8 9 10 11\n", "Data variables:\n", " phi (z) float64 1.633 1.165 1.233 1.637 ... 1.657 1.358 1.423 1.133\n", " theta (z) float64 0.6931 1.099 1.386 1.609 ... 2.079 2.197 2.303 2.398
<xarray.DataArray 'phi' (theta: 20)>\n", "array([ nan, nan, nan, nan, nan,\n", " 1.52174096, 1.33957156, 1.16664694, 1.20395286, 1.29590695,\n", " 1.5814605 , 1.56436105, 1.30629331, 1.56431596, 1.36623482,\n", " 1.22265461, nan, nan, nan, nan])\n", "Coordinates:\n", " * theta (theta) float64 0.0 0.1579 0.3158 0.4737 ... 2.526 2.684 2.842 3.0
<xarray.Dataset>\n", "Dimensions: (z: 10, zc: 11)\n", "Coordinates:\n", " * z (z) int64 2 3 4 5 6 7 8 9 10 11\n", " * zc (zc) float64 1.5 2.5 3.5 4.5 5.5 6.5 7.5 8.5 9.5 10.5 11.5\n", "Data variables:\n", " phi (z) float64 1.633 1.165 1.233 1.637 ... 1.657 1.358 1.423 1.133\n", " theta (z) float64 0.6931 1.099 1.386 1.609 ... 2.079 2.197 2.303 2.398
<xarray.DataArray 'theta_outer' (zc: 11)>\n", "array([0.34657359, 0.89587973, 1.24245332, 1.49786614, 1.70059869,\n", " 1.86883481, 2.01267585, 2.13833306, 2.24990484, 2.35024018,\n", " 1.19894764])\n", "Coordinates:\n", " * zc (zc) float64 1.5 2.5 3.5 4.5 5.5 6.5 7.5 8.5 9.5 10.5 11.5
<xarray.DataArray 'phi' (theta_outer: 19)>\n", "array([0. , 0. , 0.37785116, 0.46936054, 0.46936054,\n", " 0.48939379, 0.53079432, 0.62434095, 0.91765916, 1.18077243,\n", " 1.46775848, 1.57323631, 1.6610291 , 2.16457924, 2.03956194,\n", " 0. , 0. , 0. , 0. ])\n", "Coordinates:\n", " * theta_outer (theta_outer) float64 0.07895 0.2368 0.3947 ... 2.763 2.921
<xarray.DataArray 'temperature' (pressure: 3)>\n", "array([271.80163709, 275.48086957, 281.01789905])\n", "Coordinates:\n", " * pressure (pressure) float64 1e+05 8.5e+04 7e+04
<xarray.Dataset>\n", "Dimensions: (y: 294, x: 362, nvertex: 4, lev: 75, axis_nbounds: 2,\n", " time: 1980, y_c: 294)\n", "Coordinates: (12/13)\n", " bounds_lat (y, x, nvertex) float64 dask.array<chunksize=(294, 362, 4), meta=np.ndarray>\n", " bounds_lon (y, x, nvertex) float64 dask.array<chunksize=(294, 362, 4), meta=np.ndarray>\n", " lat (y, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " * lev (lev) float64 0.5058 1.556 2.668 ... 5.698e+03 5.902e+03\n", " lev_bounds (lev, axis_nbounds) float64 dask.array<chunksize=(75, 2), meta=np.ndarray>\n", " lon (y, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " ... ...\n", " time_bounds (time, axis_nbounds) object dask.array<chunksize=(1980, 2), meta=np.ndarray>\n", " bounds_lat_v (y_c, x, nvertex) float64 dask.array<chunksize=(294, 362, 4), meta=np.ndarray>\n", " bounds_lon_v (y_c, x, nvertex) float64 dask.array<chunksize=(294, 362, 4), meta=np.ndarray>\n", " lat_v (y_c, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " lon_v (y_c, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " areacello (y, x) float32 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", "Dimensions without coordinates: y, x, nvertex, axis_nbounds, y_c\n", "Data variables:\n", " thetao (time, lev, y, x) float32 dask.array<chunksize=(4, 75, 294, 362), meta=np.ndarray>\n", " so (time, lev, y, x) float32 dask.array<chunksize=(5, 75, 294, 362), meta=np.ndarray>\n", " vo (time, lev, y_c, x) float32 dask.array<chunksize=(3, 75, 294, 362), meta=np.ndarray>\n", "Attributes: (12/57)\n", " CMIP6_CV_version: cv=6.2.3.0-7-g2019642\n", " Conventions: CF-1.7 CMIP-6.2\n", " EXPID: CNRM-ESM2-1_historical_r1i1p1f2\n", " activity_id: CMIP\n", " arpege_minor_version: 6.3.2\n", " branch_method: standard\n", " ... ...\n", " xios_commit: 1442-shuffle\n", " status: 2019-11-05;created;by nhn2@columbia.edu\n", " netcdf_tracking_ids: hdl:21.14100/9c34b796-c31d-4c1f-be90-21d032267f6...\n", " version_id: v20181206\n", " intake_esm_varname: None\n", " intake_esm_dataset_key: CMIP.CNRM-CERFACS.CNRM-ESM2-1.historical.r1i1p1f...
<xarray.DataArray 'so' (time: 1980, y: 294, x: 362, lev: 10)>\n", "dask.array<transpose, shape=(1980, 294, 362, 10), dtype=float32, chunksize=(5, 294, 362, 10), chunktype=numpy.ndarray>\n", "Coordinates:\n", " lat (y, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " lon (y, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " * time (time) object 1850-01-16 12:00:00 ... 2014-12-16 12:00:00\n", " areacello (y, x) float32 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " * lev (lev) int64 0 50 100 150 200 250 300 350 400 450\n", "Dimensions without coordinates: y, x
<xarray.DataArray 'so' (time: 1980, y: 294, x: 362, thetao: 38)>\n", "dask.array<transpose, shape=(1980, 294, 362, 38), dtype=float32, chunksize=(4, 294, 362, 38), chunktype=numpy.ndarray>\n", "Coordinates:\n", " lat (y, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " lon (y, x) float64 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " * time (time) object 1850-01-16 12:00:00 ... 2014-12-16 12:00:00\n", " areacello (y, x) float32 dask.array<chunksize=(294, 362), meta=np.ndarray>\n", " * thetao (thetao) int64 -2 -1 0 1 2 3 4 5 6 ... 27 28 29 30 31 32 33 34 35\n", "Dimensions without coordinates: y, x