hs-paleoElevation-access.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "48d76092",
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")\n",
    "\n",
    "import numpy as np\n",
    "import xarray as xr\n",
    "import pandas as pd\n",
    "\n",
    "%config InlineBackend.figure_formats = ['svg']\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import cartopy\n",
    "import cartopy.crs as ccrs\n",
    "import cartopy.feature as cfeature\n",
    "from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER\n",
    "cartopy.config['data_dir'] = os.getenv('CARTOPY_DIR', cartopy.config.get('data_dir'))\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d678cc98",
   "metadata": {},
   "source": [
    "# List of available time intervals"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1e363433",
   "metadata": {},
   "outputs": [],
   "source": [
    "bristol_model = pd.read_csv('dataset/bristol_sim_list.csv')\n",
    "bristol_model"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "19b5f4aa",
   "metadata": {},
   "source": [
    "# Visualise the paleo-elevation dataset\n",
    "\n",
    "### Dataset link\n",
    "\n",
    "Salles, T., L. Husson, M. Lorcery, B. Halder Boggiani (2022). Paleo-Physiography Elevation-only Dataset, HydroShare, http://www.hydroshare.org/resource/b3f1e3581d174bf58b00ba5672604710\n",
    "\n",
    "### THREDDS Data Service\n",
    "\n",
    "https://thredds.hydroshare.org/thredds/catalog/hydroshare/resources/b3f1e3581d174bf58b00ba5672604710/data/contents/catalog.html\n",
    "\n",
    "\n",
    "Plotting a specific time interval is done by specifying a time "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a050770f",
   "metadata": {},
   "outputs": [],
   "source": [
    "time = 145"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "abe40a38",
   "metadata": {},
   "source": [
    "We then load the netcdf dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d3ad985f",
   "metadata": {},
   "outputs": [],
   "source": [
    "url = 'http://thredds.hydroshare.org/thredds/dodsC/hydroshare/resources/b3f1e3581d174bf58b00ba5672604710/data/contents/z'+str(time)+'Ma.nc'\n",
    "elev_ds = xr.open_dataset(url)\n",
    "dsz = elev_ds.z.load()\n",
    "dsz"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "04fc213b",
   "metadata": {},
   "source": [
    "## Global map\n",
    "\n",
    "For plotting global dataset, we will coarsen the 0.05 degree resolution, by specifying a coarsen factor (e.g. 2,3,4,...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dabc136b",
   "metadata": {},
   "outputs": [],
   "source": [
    "coarsen_factor = 4\n",
    "coarsez = dsz.coarsen(longitude=coarsen_factor, \n",
    "                      boundary='trim').mean().coarsen(latitude=coarsen_factor, \n",
    "                                                      boundary='trim').mean()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f3b65850",
   "metadata": {},
   "source": [
    "Let's make the plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7f5c85f0",
   "metadata": {},
   "outputs": [],
   "source": [
    "cbar_kwargs = {'orientation':'vertical', 'shrink':0.3}\n",
    "fig = plt.figure(figsize=(10, 8))\n",
    "ax = fig.add_subplot(111, projection=ccrs.Mollweide(), aspect='auto')\n",
    "coarsez.plot.imshow(ax=ax, transform=ccrs.PlateCarree(),\n",
    "                x='longitude', y='latitude',\n",
    "                robust=True, vmin=-6000, vmax=6000,\n",
    "                cmap='RdBu_r',\n",
    "                cbar_kwargs=cbar_kwargs,\n",
    "                add_colorbar=True,\n",
    "                extend='both')\n",
    "\n",
    "ax.gridlines(color='white',linestyle='--', \n",
    "             xlocs=np.arange(-180,180,45), ylocs=np.arange(-90,90,45))\n",
    "\n",
    "cnplot = ax.contour(coarsez.longitude, coarsez.latitude, coarsez,\n",
    "                        linewidths=1.5,\n",
    "                        colors='k',\n",
    "                        levels=1, \n",
    "                        transform=ccrs.PlateCarree())\n",
    "\n",
    "ax.set_title('Elevation @ '+str(time)+'Ma', fontsize=10, fontweight=\"bold\", loc='left')\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "795ba1a5",
   "metadata": {},
   "source": [
    "## Local map\n",
    "\n",
    "For a local plot, we will use a higher resolution and slice the dataset to the region of interest:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "91c080db",
   "metadata": {},
   "outputs": [],
   "source": [
    "coarsen_factor = 1\n",
    "midz = dsz.coarsen(longitude=coarsen_factor, \n",
    "                      boundary='trim').mean().coarsen(latitude=coarsen_factor, \n",
    "                                                      boundary='trim').mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5969c5a9",
   "metadata": {},
   "outputs": [],
   "source": [
    "min_lon = -20    # lower left longitude\n",
    "min_lat = -30    # lower left latitude\n",
    "max_lon = 0      # upper right longitude\n",
    "max_lat = -10    # upper right latitude\n",
    "\n",
    "# Defining region boundaries\n",
    "lon_bnds = [min_lon, max_lon]\n",
    "lat_bnds = [min_lat, max_lat]\n",
    "\n",
    "# Performing the reduction\n",
    "dsz_clip = midz.sel(latitude=slice(*lat_bnds), longitude=slice(*lon_bnds))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ee7182a",
   "metadata": {},
   "outputs": [],
   "source": [
    "cbar_kwargs = {'orientation':'vertical', 'shrink':0.3}\n",
    "fig = plt.figure(figsize=(8, 8))\n",
    "ax = fig.add_subplot(111, projection=ccrs.PlateCarree(), aspect='auto')\n",
    "dsz_clip.plot.imshow(ax=ax, transform=ccrs.PlateCarree(),\n",
    "                x='longitude', y='latitude',\n",
    "                robust=True, vmin=-500, vmax=3000,\n",
    "                cmap='terrain',\n",
    "                cbar_kwargs=cbar_kwargs,\n",
    "                add_colorbar=True,\n",
    "                extend='both')\n",
    "\n",
    "cnplot = ax.contour(dsz_clip.longitude, dsz_clip.latitude, dsz_clip,\n",
    "                        linewidths=1.5,\n",
    "                        colors='k',\n",
    "                        levels=[0], \n",
    "                        transform=ccrs.PlateCarree())\n",
    "\n",
    "gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True,\n",
    "                  linewidth=0.1, color='k', alpha=1, \n",
    "                  linestyle='--')\n",
    "gl.top_labels = False\n",
    "gl.right_labels = False\n",
    "gl.xformatter = LONGITUDE_FORMATTER\n",
    "gl.yformatter = LATITUDE_FORMATTER\n",
    "gl.xlabel_style = {'size': 8}\n",
    "gl.ylabel_style = {'size': 8} \n",
    "\n",
    "ax.set_title('Local elevation @ '+str(time)+'Ma', fontsize=10, fontweight=\"bold\", loc='left')\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a44423cb",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2c54aa9d",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.9.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}