Last active
June 15, 2017 14:52
-
-
Save ocefpaf/a20b13940ee9bcfc8dd7c5ec5790ebe1 to your computer and use it in GitHub Desktop.
notebook/folium/tri_tide_movie.ipynb
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
{ | |
"cells": [ | |
{ | |
"metadata": {}, | |
"cell_type": "markdown", | |
"source": "# Plotting predicted tides\n\n\nMust call with higher `NotebookApp.iopub_data_rate_limit` to display the animations.\n\nE.g: `jupyter 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\n# http://nbviewer.jupyter.org/gist/anonymous/410bf777e19ac8ca14bb06a24b03745f\n# on Python 3.\nprint(netCDF4.__version__)", | |
"execution_count": 1, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": "1.2.4\n" | |
} | |
] | |
}, | |
{ | |
"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": 2, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": "odict_keys(['tidenames', 'tidefreqs', 'ele', 'lon', 'lat', 'depth', 'depth-averaged', 'u_amp', 'v_amp', 'u_phase', 'v_phase'])\n" | |
} | |
] | |
}, | |
{ | |
"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": 3, | |
"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": 4, | |
"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": 5, | |
"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": 6, | |
"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\ntrif = Delaunay(list(zip(lonf, latf))).vertices", | |
"execution_count": 7, | |
"outputs": [] | |
}, | |
{ | |
"metadata": { | |
"collapsed": true, | |
"trusted": true | |
}, | |
"cell_type": "code", | |
"source": "import 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": 8, | |
"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": 9, | |
"outputs": [] | |
}, | |
{ | |
"metadata": {}, | |
"cell_type": "markdown", | |
"source": "Find the indices of the tidal constituents." | |
}, | |
{ | |
"metadata": { | |
"collapsed": true, | |
"trusted": true | |
}, | |
"cell_type": "code", | |
"source": "k = 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": 10, | |
"outputs": [] | |
}, | |
{ | |
"metadata": {}, | |
"cell_type": "markdown", | |
"source": "Download everything first to speed the slice below." | |
}, | |
{ | |
"metadata": { | |
"collapsed": true, | |
"trusted": true | |
}, | |
"cell_type": "code", | |
"source": "uamp = nc['u_amp'][0, ...]\nvamp = nc['v_amp'][0, ...]\n\nupha = nc['u_phase'][0, ...]\nvpha = nc['v_phase'][0, ...]", | |
"execution_count": 11, | |
"outputs": [] | |
}, | |
{ | |
"metadata": {}, | |
"cell_type": "markdown", | |
"source": "Slice the data for plotting." | |
}, | |
{ | |
"metadata": { | |
"collapsed": true, | |
"trusted": true | |
}, | |
"cell_type": "code", | |
"source": "uamp = uamp[inbox, :][:, ind_nc]\nvamp = vamp[inbox, :][:, ind_nc]\n\nupha = upha[inbox, :][:, ind_nc]\nvpha = vpha[inbox, :][:, ind_nc]", | |
"execution_count": 12, | |
"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": 13, | |
"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": 14, | |
"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": 15, | |
"outputs": [] | |
}, | |
{ | |
"metadata": { | |
"collapsed": true, | |
"trusted": true | |
}, | |
"cell_type": "code", | |
"source": "from utide.harmonics import FUV\n\n\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": 16, | |
"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": 17, | |
"outputs": [] | |
}, | |
{ | |
"metadata": { | |
"trusted": true | |
}, | |
"cell_type": "code", | |
"source": "import json\nimport mplleaflet\nimport matplotlib.pyplot as plt\n\n\nfig, ax = plt.subplots(figsize=(10, 5))\n\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 = ax.quiver(lon, lat, U, V, scale=vscale, color='white')\nax.axis(bbox) # Vineyard sound 2.\nmplleaflet.display(fig=fig, tiles='esri_aerial')", | |
"execution_count": 18, | |
"outputs": [ | |
{ | |
"data": { |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
I cant get it work