gridgeo tour

import gridgeo


url = 'http://crow.marine.usf.edu:8080/thredds/dodsC/FVCOM-Nowcast-Agg.nc'

grid = gridgeo.GridGeo(
    url,
    standard_name='sea_water_potential_temperature'
)

import fiona

schema = {
    'geometry': 'MultiPolygon',
    'properties': {'name': f'str:{len(grid.mesh)}'}
}

with fiona.open('grid.shp', 'w', 'ESRI Shapefile', schema) as f:
    f.write(
        {
            'geometry': grid.__geo_interface__,
            'properties': {'name': grid.mesh}
        }
    )

CPLE_NotSupported in driver ESRI Shapefile does not support creation option ENCODING

Methods

[s for s in dir(grid) if not s.startswith('_')]

['geometry',
 'mesh',
 'outline',
 'polygons',
 'save',
 'to_geojson',
 'triang',
 'x',
 'y']

grid.mesh

'ugrid'

print(f'The grid has {len(grid.geometry)} polygons, showing the first 5.')

grid.geometry[:5]

The grid has 98818 polygons, showing the first 5.

grid.outline

Displaying all the polygons as vectors can be costly and crash jupyter! Let’s make a raster representation of them using cartopy instead.

%matplotlib inline
import cartopy.crs as ccrs
import matplotlib.pyplot as plt


fig, ax = plt.subplots(
    figsize=(12, 12),
    subplot_kw={'projection': ccrs.PlateCarree()}
)

kw = dict(linestyle='-', alpha=0.25, color='darkgray')
ax.triplot(grid.triang, **kw)
ax.coastlines(resolution='10m');

to_geojson() method returns a styled geojson-like dict

See https://github.com/mapbox/simplestyle-spec/tree/master/1.1.0 for styling options.

kw = {
    'fill': '#fd7d11',
    'fill_opacity': 0.2,
    'stroke_opacity': 1,
    'float_precision': 2,
}

geojson = grid.to_geojson(**kw)
geojson['properties']

{'title': 'ugrid',
 'description': '',
 'marker-size': 'medium',
 'marker-symbol': '',
 'marker-color': '7e7e7e',
 'stroke': '555555',
 'stroke-opacity': 1,
 'stroke-width': 2,
 'fill': '#fd7d11',
 'fill-opacity': 0.6}

or just use the __geo_interface__.

grid.__geo_interface__.keys()

dict_keys(['type', 'coordinates'])

Saving the grid to as geojson file

grid.save('grid.geojson', **kw)

shapefile

grid.save('grid.shp')

or just plot using folium ;-)

import folium


x, y = grid.outline.centroid.xy

m = folium.Map(location=[y[0], x[0]])
folium.GeoJson(grid.outline.__geo_interface__).add_to(m)
m.fit_bounds(m.get_bounds())

m