notebook/folium/tri_tide_movie.ipynb

tri_tide_movie-signell.ipynb

{
  "cells": [
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "#  Plotting predicted tides\n\njupyter notebook --NotebookApp.iopub_data_rate_limit=10000000 tri_tide_movie-signell.ipynb"
    },
    {
      "metadata": {
        "trusted": true
      },
      "cell_type": "code",
      "source": "import netCDF4\n\n\n# Versions >1.2.4 are failing with http://nbviewer.jupyter.org/gist/anonymous/410bf777e19ac8ca14bb06a24b03745f\nprint(netCDF4.__version__)",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "trusted": true
      },
      "cell_type": "code",
      "source": "from netCDF4 import Dataset\n\n\nurl = 'http://geoport.whoi.edu/thredds/dodsC/usgs/vault0/models/tides/vdatum_gulf_of_maine/adcirc54_38_orig.nc'\n\nnc = Dataset(url)\n\nprint(nc.variables.keys())",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "units = {\n    'knots': 1.9438,\n    'm/s': 1.0\n}\nconsts = ['STEADY', 'M2', 'S2', 'N2', 'K1', 'O1', 'P1', 'M4', 'M6']\n\n# Vineyard sound 2.\nbbox = [-70.7234, -70.4532, 41.4258, 41.5643]\n\nhalo = 0.1\n\nax2 = [\n    bbox[0] - halo * (bbox[1] - bbox[0]),\n    bbox[1] + halo * (bbox[1] - bbox[0]),\n    bbox[2] - halo * (bbox[3] - bbox[2]),\n    bbox[3] + halo * (bbox[3] - bbox[2])\n]",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "import pytz\nfrom datetime import datetime\nfrom pandas import date_range\n\n\nfmt = '%d-%b-%Y %H:%M'\ntzinfo = pytz.utc\n\nstart = datetime.strptime('18-Sep-2015 05:00', fmt).replace(tzinfo=tzinfo)\nstop = datetime.strptime('19-Sep-2015 05:00', fmt).replace(tzinfo=tzinfo)\n\nglocals = date_range(start, stop, freq='1H').to_pydatetime()\n\nntimes = len(glocals)",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "def parse_string(name):\n    try:\n        const = ''.join(name).strip()\n    except TypeError:\n        const = b''.join(name).strip()\n        const = const.decode()\n    return const ",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "data = nc['tidenames'][:]\n\nnames = []\nfor name in data.tolist():\n    names.append(parse_string(name))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "from scipy.spatial import Delaunay\n\n\ndepth = nc['depth'][:]\nlatf = nc['lat'][:]\nlonf = nc['lon'][:]\nfrequency = nc['tidefreqs'][:]\n\n# Not sure why this is not working.\n# trif = cubes.extract_strict('Horizontal Element Incidence List').data\ntrif = Delaunay(list(zip(lonf, latf))).vertices",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "# Find indices in box.\nimport numpy as np\n\n\ninbox = np.logical_and(\n    np.logical_and(lonf >= ax2[0], lonf <= ax2[1]),\n    np.logical_and(latf >= ax2[2], latf <= ax2[3])\n)\n\nlon = lonf[inbox]\nlat = latf[inbox]",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "from scipy.io import loadmat\n\n\nmat = 't_constituents.mat'\ncon_info = loadmat(mat, squeeze_me=True)\n# ignoring shallow water and sat constants.\ncon_info = con_info['const']",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "scrolled": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "from utide import _ut_constants_fname, loadbunch\n\n\ncon_info = loadbunch(_ut_constants_fname)['const']",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "# Find the indices of the tidal constituents.\n\nk = 0\nind_nc, ind_ttide = [], []\n\nconst_name = [e.strip() for e in con_info['name'].tolist()]\n\nfor name in consts:\n    try:\n        if name == 'STEADY':\n            indx = const_name.index('Z0')\n        else:\n            indx = const_name.index(name)\n        k += 1\n        ind_ttide.append(indx)\n        ind_nc.append(names.index(name))\n    except ValueError:\n        pass  # `const` not found.",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "uamp = nc['u_amp'][0, inbox, :][:, ind_nc]\nvamp = nc['v_amp'][0, inbox, :][:, ind_nc]\n\nupha = nc['u_phase'][0, inbox, :][:, ind_nc]\nvpha = nc['v_phase'][0, inbox, :][:, ind_nc]",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "freq_nc = frequency[ind_nc]\n\nfreq_ttide = con_info['freq'][ind_ttide]\n\nt_tide_names = np.array(const_name)[ind_ttide]",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "omega_ttide = 2*np.pi * freq_ttide  # Convert from radians/s to radians/hour.\n\nomega = freq_nc * 3600\n\nrllat = 55  # Reference latitude for 3rd order satellites (degrees) (55 is fine always)\n\nvscale = 10  # smaller makes arrows bigger",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "from matplotlib.dates import date2num\n\n\n# Convert to Matlab datenum.\n# (Soon UTide will take python datetime objects.)\njd_start = date2num(start) + 366.1667",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "from utide.harmonics import FUV\n\n\n# NB: I am not a 100% sure if this is identical to what we had with t_tide.\n# ngflgs -> [NodsatLint NodsatNone GwchLint GwchNone]\nv, u, f = FUV(\n    t=np.array([jd_start]),\n    tref=np.array([0]),\n    lind=np.array([ind_ttide]),\n    lat=rllat,\n    ngflgs=[0, 0, 0, 0]\n)",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "# Convert phase in radians.\nv, u, f = map(np.squeeze, (v, u, f))\nv = v * 2 * np.pi\nu = u * 2 * np.pi\n\nthours = np.array(\n    [d.total_seconds() for d in (glocals - glocals[0])]\n) / 60 / 60.",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "trusted": true,
        "collapsed": true
      },
      "cell_type": "code",
      "source": "import json\nimport mplleaflet\nimport matplotlib.pyplot as plt\n\n\nfig = plt.figure(figsize=(10,5))\nU = (f * uamp * np.cos(v + thours[k] * omega + u - upha * np.pi/180)).sum(axis=1)\nV = (f * vamp * np.cos(v + thours[k] * omega + u - vpha * np.pi/180)).sum(axis=1)\n\nw = units['knots'] * (U + 1j * V)\n\nwf = np.NaN * np.ones_like(lonf, dtype=w.dtype)\nwf[inbox] = w\n\n# FIXME: Cannot use masked arrays and tricontour!\n# wf = ma.masked_invalid(wf)\n# cs = ax.tricontour(lonf, latf, trif, np.abs(wf).filled(fill_value=0))\n# fig.colorbar(cs)\nq = plt.quiver(lon, lat, U, V, scale=vscale, color='white')\nplt.axis(bbox)  # Vineyard sound 2.\n#q.set_title('{}'.format(glocals[k]))\nmplleaflet.display(fig=fig, tiles='esri_aerial')",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "collapsed": true,
        "trusted": true
      },
      "cell_type": "code",
      "source": "js = mplleaflet.fig_to_geojson(fig=fig)\nwith open('test.json', 'w') as outfile:\n    json.dump(js, outfile, indent=4)",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {
        "trusted": true,
        "collapsed": true
      },
      "cell_type": "code",
      "source": "%matplotlib inline\nimport matplotlib.pyplot as plt\nfrom JSAnimation import IPython_display\nfrom matplotlib.animation import FuncAnimation\n\n\ndef update_figure(k):\n    global ax, fig\n    ax.cla()\n    \n    U = (f * uamp * np.cos(v + thours[k] * omega + u - upha * np.pi/180)).sum(axis=1)\n    V = (f * vamp * np.cos(v + thours[k] * omega + u - vpha * np.pi/180)).sum(axis=1)\n    \n    w = units['knots'] * (U + 1j * V)\n    \n    wf = np.NaN * np.ones_like(lonf, dtype=w.dtype)\n    wf[inbox] = w\n\n    # FIXME: Cannot use masked arrays and tricontour!\n    # wf = ma.masked_invalid(wf)\n    # cs = ax.tricontour(lonf, latf, trif, np.abs(wf).filled(fill_value=0))\n    # fig.colorbar(cs)\n    q = ax.quiver(lon, lat, U, V, scale=vscale)\n    ax.axis(bbox)  # Vineyard sound 2.\n    ax.set_title('{}'.format(glocals[k]))\n\nfig, ax = plt.subplots(figsize=(10, 5))\nFuncAnimation(fig, update_figure, interval=100, frames=ntimes)",
      "execution_count": null,
      "outputs": []
    }
  ],
  "metadata": {
    "gist": {
      "id": "",
      "data": {
        "description": "notebook/folium/tri_tide_movie.ipynb",
        "public": true
      }
    },
    "gist_id": "ab5133b7feecba1dfacd",
    "kernelspec": {
      "name": "python3",
      "display_name": "Python [default]",
      "language": "python"
    },
    "language_info": {
      "name": "python",
      "version": "3.6.1",
      "mimetype": "text/x-python",
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "pygments_lexer": "ipython3",
      "nbconvert_exporter": "python",
      "file_extension": ".py"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 1
}