Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Added ipymesh example #10

Open
wants to merge 8 commits into
base: master
Choose a base branch
from
Open

Added ipymesh example #10

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
bafed07
added ipymesh example
cekees Jun 6, 2019
920b654
Update README.md
cekees Jun 6, 2019
64fcd6e
added leaflet example
cekees Jun 9, 2019
bffeb57
add ipympl requirement
cekees Jun 9, 2019
2c55d3a
updated index
cekees Jun 9, 2019
a824c94
added ipyleaflet
cekees Jun 10, 2019
fa4bf43
udpate versions in requirements
cekees Jun 10, 2019
ce79c05
removed h5py import and cleaned up other imports
cekees Jun 10, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion README.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,4 +3,4 @@

Available on [http://voila-gallery.pascalbugnion.net](http://voila-gallery.pascalbugnion.net/voila/render/index).

To play with the code, open this in [Binder](https://mybinder.org/v2/gh/pbugnion/voila-gallery/master?urlpath=voila%2Ftree).
To play with the code, open this in [Binder](https://mybinder.org/v2/gh/erdc/voila-gallery/master?urlpath=voila%2Ftree).
6 changes: 4 additions & 2 deletions index.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,9 @@
"\n",
"- [country-indicators](country-indicators/index) shows how to integrate Matplotlib plots using output widgets.\n",
"- [gaussian-density](gaussian-density/index) shows how to build a dashboard around a bqplot plot.\n",
"- [stl-display](stl-display/index) shows how to read STL files and display them with [ipyvolume](https://github.com/maartenbreddels/ipyvolume)."
"- [stl-display](stl-display/index) shows how to read STL files and display them with [ipyvolume](https://github.com/maartenbreddels/ipyvolume).\n",
"- [mesh-generation-ipymesh](mesh-generation-ipymesh/index) shows how to construct a Planar Straight Line Graph with [ipymesh](https://github.com/erdc/ipymesh) and triangulate it with [triangle](https://www.cs.cmu.edu/~quake/triangle.html).\n",
"- [mesh-generation-ipyleaflet](mesh-generation-ipyleaflet/index) shows how to construct a Planar Straight Line Graph with [ipyleaflet](https://github.com/jupyter-widgets/ipymesh) and triangulate it with [triangle](https://www.cs.cmu.edu/~quake/triangle.html), adaptively, to approximate complex bathymetry/terrain from the GEBCO global grid.\n"
]
},
{
Expand All @@ -41,7 +43,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.4"
"version": "3.7.2"
}
},
"nbformat": 4,
Expand Down
256 changes: 256 additions & 0 deletions mesh-generation-ipyleaflet/index.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,256 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In this app we:\n",
"* Create a very simple [Planar Straight Line Graph](https://en.wikipedia.org/wiki/Planar_straight-line_graph) with the [ipyleaflet](https://github.com/jupygter-widgets/ipyleaflet) widget. In fact, it's just a polygon, and the rest of the app can only handle a single polygon at a time right now--so trash your existing polygon before creating another one. \n",
"* Pull the terrain/bathymetry data for the area of interest from a global database. The notebook has code to use [GEBCO 2019](https://www.gebco.net). In fact, the app is pulling from a small grid around the initial coordinates to keep the docker container small.\n",
"* Generate a triangular mesh with the [triangle](https://www.cs.cmu.edu/~quake/triangle.html) mesh generator using the [triangle python bindings](https://rufat.be/triangle/index.html) and adaptively refine it to resolve the terrain data.\n",
"* Visualize the resulting mesh in several ways.\n",
"Find the code [here](https://github.com/pbugnion/voila-gallery/blob/master/mesh-generation-ipyleaflet/index.ipynb)."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"%matplotlib widget\n",
"import ipywidgets\n",
"from ipyleaflet import (\n",
" Map,\n",
" Polygon,\n",
" DrawControl,\n",
" WMSLayer\n",
")\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import triangle as tr\n",
"from matplotlib import tri as mpltr\n",
"from scipy import interpolate as scipy_interpolate\n",
"import ipyvolume as ipv\n",
"import ipympl\n",
"\n",
"#initial map center and zoom level\n",
"center = [34.6252978589571, -77.34580993652344]\n",
"vertices=np.array([[-77.521791, 34.582935],\n",
" [-77.380628, 34.356501],\n",
" [-77.11704 , 34.561348],\n",
" [-77.321914, 34.735944]])\n",
"zoom = 9\n",
"#use the GEBCO tile layers in the map widget--this is not necessary, just for viz\n",
"wms = WMSLayer(url='http://www.gebco.net/data_and_products/gebco_web_services/web_map_service/mapserv', \n",
" layers='GEBCO_LATEST', \n",
" )\n",
"#create a drawing control \n",
"dc = DrawControl(polyline={},#\n",
" marker={},#{'shapeOptions': {'color': '#0000FF'}},\n",
" rectangle={},#{'shapeOptions': {'color': '#0000FF'}},\n",
" circle={},#{'shapeOptions': {'color': '#0000FF'}},\n",
" circlemarker={}#,\n",
" )\n",
"\n",
"m = Map(center=center, zoom=zoom, layers=(wms,), controls=(dc,))"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"plt.ioff()\n",
"gridfig=plt.figure(num=0)\n",
"gridfig.set_label(' ')\n",
"#import h5py\n",
"#f = h5py.File('GEBCO_2019.nc','r')\n",
"f = np.load(\"terrain_slice_save.npz\")\n",
"def plotGrid(button=None):\n",
" global vertices,dc, j_min, j_max, i_min, i_max, lonr, latr, gridfig\n",
" #get the vertices of the polygon\n",
" if dc.last_draw['geometry'] != None:\n",
" vertices=np.array(dc.last_draw['geometry']['coordinates'][0][:-1])\n",
" lonr = [vertices[:,0].min(),vertices[:,0].max()]\n",
" latr = [vertices[:,1].min(),vertices[:,1].max()]\n",
" j_min = np.abs(f['lon'][:] - lonr[0]).argmin() - 1\n",
" j_max = np.abs(f['lon'][:] - lonr[1]).argmin() + 1\n",
" i_min = np.abs(f['lat'][:] - latr[0]).argmin() - 1\n",
" i_max = np.abs(f['lat'][:] - latr[1]).argmin() + 1\n",
" plt.figure(gridfig.number)\n",
" plt.ion()\n",
" plt.clf()\n",
" plt.contourf(f['lon'][j_min:j_max],\n",
" f['lat'][i_min:i_max],\n",
" f['elevation'][i_min:i_max,\n",
" j_min:j_max])\n",
" plt.colorbar()\n",
" plt.title('')\n",
"plotGrid()\n",
"terrain_button = ipywidgets.Button(\n",
" description='Get Terrain',\n",
" disabled=False,\n",
" button_style='', # 'success', 'info', 'warning', 'danger' or ''\n",
" tooltip='Pull the terrain for the area of interest'\n",
")\n",
"terrain_button.on_click(plotGrid)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"plt.ioff()\n",
"meshfig=plt.figure(num=1,frameon=True)\n",
"meshfig.set_label(' ')\n",
"surffig=ipv.figure()\n",
"def genMesh(button=None):\n",
" global vertices,meshfig\n",
" nv = vertices.shape[0]\n",
" vertexFlags=[1 for i in range(nv)]\n",
" segments=[]\n",
" segmentFlags=[]\n",
" for i in range(nv):\n",
" segments.append([i,(i+1)%nv])\n",
" segmentFlags.append(1)\n",
" nx = j_max - j_min\n",
" ny = i_max - i_min\n",
" bathy = np.zeros((nx*ny,3),'d')\n",
" for i in range(ny):\n",
" bathy[i*nx:(i+1)*nx,0] = f['lon'][j_min:j_max]\n",
" bathy[i*nx:(i+1)*nx,1] = f['lat'][i_min+i]\n",
" bathy[:,2] = np.reshape(f['elevation'][i_min:i_max,j_min:j_max], nx*ny)\n",
"\n",
" he = min((lonr[1]-lonr[0])/40.,(latr[1]-latr[0])/40.0)\n",
" bathyGridDim = (ny,nx)\n",
" x = bathy[:bathyGridDim[1],0]\n",
" y = bathy[:bathyGridDim[0]*bathyGridDim[1]:bathyGridDim[1],1]\n",
" z = bathy[:,2].reshape(bathyGridDim).transpose()\n",
" bathyInterpolant = scipy_interpolate.RectBivariateSpline(x,y,z,kx=1,ky=1)\n",
" pslg = dict(vertices=vertices,\n",
" vertex_attributes=vertexFlags,\n",
" segments=segments,\n",
" segment_markers=segmentFlags)\n",
" mesh = tr.triangulate(pslg,'pq30Da{0:f}'.format(0.5*he**2))\n",
" area_max = 0.5*he**2\n",
" mplmesh=mpltr.Triangulation(mesh['vertices'][:,0],\n",
" mesh['vertices'][:,1],\n",
" mesh['triangles'])\n",
" not_converged=True\n",
" tol=0.5\n",
" while not_converged:\n",
" mplmesh=mpltr.Triangulation(mesh['vertices'][:,0],\n",
" mesh['vertices'][:,1],\n",
" mesh['triangles'])\n",
" trifinder = mpltr.TrapezoidMapTriFinder(mplmesh)\n",
" zNodes = bathyInterpolant.ev(mesh['vertices'][:,0],\n",
" mesh['vertices'][:,1])\n",
" interpolator = mpltr.LinearTriInterpolator(mplmesh, zNodes, trifinder=trifinder)\n",
" zInterp = interpolator(bathy[:,0],bathy[:,1])\n",
" xE = bathy[np.abs(bathy[:,2] - zInterp) > tol,0]\n",
" yE = bathy[np.abs(bathy[:,2] - zInterp) > tol,1]\n",
" e = trifinder(xE,yE)\n",
" area=np.ones((mesh['triangles'].shape[0],),'d')*0.5*he**2\n",
" area_max *= 0.5\n",
" area[e[e>=0]] = area_max\n",
" mesh['triangle_max_area']=area\n",
" mesh2 = tr.triangulate(mesh,'rpq30Da'.format(0.5*he**2))\n",
" err = np.abs(bathy[:,2] - zInterp) > tol\n",
" #print(\"L_infty_error\",np.abs(bathy[:,2] - zInterp).max())\n",
" not_converged = err.any()\n",
" mesh=mesh2\n",
" fig=plt.figure(meshfig.number)\n",
" plt.ion()\n",
" plt.clf()\n",
" plt.triplot(mplmesh)\n",
" fig = ipv.figure(surffig)\n",
" if len(fig.meshes):\n",
" del fig.meshes[-1]\n",
" ipv.style.use('minimal')\n",
" ipv.plot_trisurf(x=mesh['vertices'][:,0],\n",
" y=mesh['vertices'][:,1],\n",
" z=zNodes,\n",
" triangles=mesh['triangles'],\n",
" color='grey')\n",
" #,color=color, wireframe=wireframe)\n",
" ipv.xlim(lonr[0], lonr[1])\n",
" ipv.ylim(latr[0], latr[1])\n",
" ipv.zlim(zNodes.min()*10,zNodes.max()*10)\n",
" mesh3d = ipv.gcc()\n",
"button = ipywidgets.Button(\n",
" description='Generate Mesh',\n",
" disabled=False,\n",
" button_style='', # 'success', 'info', 'warning', 'danger' or ''\n",
" tooltip='Regenerate the mesh'\n",
")\n",
"genMesh()\n",
"button.on_click(genMesh)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "4a9030355a504927b16682c93ce8b80b",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Tab(children=(Map(basemap={'url': 'https://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png', 'max_zoom': 19, 'attri…"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"from ipywidgets import Tab,VBox\n",
"app=Tab(children=[m,\n",
" VBox(children=[terrain_button, gridfig.canvas]),\n",
" VBox(children=[button, meshfig.canvas]),\n",
" surffig])\n",
"app.set_title(0,\"Map\")\n",
"app.set_title(1,\"Terrain Data\")\n",
"app.set_title(2,\"Mesh\")\n",
"app.set_title(3,\"Terrain Surface\")\n",
"display(app)\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"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.7.2"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Binary file added mesh-generation-ipyleaflet/terrain_slice_save.npz
View file
Open in desktop
Binary file not shown.
Loading