adcirc_pangeo_earthsim_demo.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Read ADCIRC output from Zarr (cloud-optimized NetCDF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "from dask.distributed import Client, progress, LocalCluster\n",
    "from dask_kubernetes import KubeCluster, make_pod_spec\n",
    "import xarray as xr\n",
    "import s3fs\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Start a dask cluster to crunch the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this is the YAML approach:\n",
    "#cluster = KubeCluster.from_yaml('/home/jovyan/custom-worker-template.yaml')\n",
    "\n",
    "\n",
    "# this is the local function approach:\n",
    "pod_spec = make_pod_spec(image='esip/pangeo-notebook:2018-07-04',\n",
    "                          memory_limit='6G', memory_request='4G',\n",
    "                          cpu_limit=1, cpu_request=1,\n",
    "                          env={},\n",
    "                          extra_container_config={\n",
    "                              #'securityContext': {'capabilities': {'add': '[SYS_ADMIN]'},\n",
    "                              #                    'privileged': 'true'\n",
    "                              #                   },\n",
    "                              #'volumeMounts': {'name': 's3', 'mountPath': '/s3'},\n",
    "                              #'resources': {'limits': {'cpu': '1.75', 'memory': '6G'},\n",
    "                              #              'requests': {'cpu': '1.75', 'memory': '6G'}\n",
    "                              #             }\n",
    "                             },\n",
    "                          extra_pod_config={}\n",
    "            \n",
    "                        )\n",
    "\n",
    "cluster = KubeCluster(pod_spec)\n",
    "\n",
    "cluster.scale(30)\n",
    "cluster"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(pod_spec)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "client = Client(cluster)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "client"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# jetstream s3\n",
    "# url='https://iu.jetstream-cloud.org:8080'\n",
    "# fs = s3fs.S3FileSystem(client_kwargs=dict(endpoint_url=url), anon=True)\n",
    "# s3map = s3fs.S3Map('rsignell/adcirc_test01', s3=fs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# AWS s3\n",
    "fs = s3fs.S3FileSystem(anon=True)\n",
    "s3map = s3fs.S3Map('rsignell/adcirc_test01', s3=fs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds = xr.open_zarr(s3map)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds['zeta']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds['zeta'].nbytes/1.e9"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ds['zeta'].max(dim='time')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "max_var = ds['zeta'].max(dim='time').persist()\n",
    "progress(max_var)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize data on mesh using Datashader"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import datashader as dshade\n",
    "import holoviews as hv\n",
    "import geoviews as gv\n",
    "import cartopy.crs as ccrs\n",
    "\n",
    "from holoviews.operation.datashader import datashade, rasterize\n",
    "from colorcet import cm_n\n",
    "from matplotlib.cm import jet\n",
    "\n",
    "datashade.precompute = True\n",
    "\n",
    "hv.extension('bokeh')\n",
    "%opts Image RGB VectorField [width=800 height=600]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "v = np.vstack((ds['x'], ds['y'], max_var)).T\n",
    "verts = pd.DataFrame(v, columns=['x','y','z'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "points = gv.operation.project_points(gv.Points(verts, vdims=['z']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tris = pd.DataFrame(ds['element'].values.astype('int')-1, columns=['v0','v1','v2'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tiles = gv.WMTS('https://maps.wikimedia.org/osm-intl/{Z}/{X}/{Y}@2x.png')\n",
    "value = 'max water level'\n",
    "label = '{} (m)'.format(value)\n",
    "trimesh = gv.TriMesh((tris, points), label=label)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "%%opts Image [colorbar=True] (cmap=jet)\n",
    "\n",
    "error_meshes = rasterize(trimesh,aggregator=dshade.mean('z'))\n",
    "tiles * error_meshes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Extract a time series at a specified lon, lat location"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# find the indices of the points in (x,y) closest to the points in (xi,yi)\n",
    "def nearxy(x,y,xi,yi):\n",
    "    ind = np.ones(len(xi),dtype=int)\n",
    "    for i in range(len(xi)):\n",
    "        dist = np.sqrt((x-xi[i])**2+(y-yi[i])**2)\n",
    "        ind[i] = dist.argmin()\n",
    "    return ind"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#just offshore of Galveston\n",
    "lat = 29.2329856\n",
    "lon = -95.1535041"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ind = nearxy(ds['x'].values,ds['y'].values,[lon], [lat])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "ds['zeta'][:,ind].plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [default]",
   "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}