Geometries

import pandas as pd
import holoviews as hv
import geoviews as gv
import geoviews.feature as gf
import cartopy
import cartopy.feature as cf

from geoviews import opts
from cartopy import crs as ccrs

gv.extension('matplotlib', 'bokeh')

gv.output(dpi=120, fig='svg')

Cartopy and shapely make working with geometries and shapes very simple, and GeoViews provides convenient wrappers for the various geometry types they provide. In addition to Path and Polygons types, which draw geometries from lists of arrays or a geopandas DataFrame, GeoViews also provides the Feature and Shape types, which wrap cartopy Features and shapely geometries respectively.

Feature

The Feature Element provides a very convenient means of overlaying a set of basic geographic features on top of or behind a plot. The cartopy.feature module provides various ways of loading custom features, however geoviews provides a number of default features which we have imported as gf, amongst others this includes coastlines, country borders, and land masses. Here we demonstrate how we can plot these very easily, either in isolation or overlaid:

(gf.ocean + gf.land + gf.ocean * gf.land * gf.coastline * gf.borders).cols(3)

These default features simply wrap around cartopy Features, therefore we can easily load a custom NaturalEarthFeature such as graticules at 30 degree intervals:

graticules = cf.NaturalEarthFeature(
    category='physical',
    name='graticules_30',
    scale='110m')

(gf.ocean() * gf.land() * gv.Feature(graticules, group='Lines') * gf.borders * gf.coastline).opts(
    opts.Feature('Lines', projection=ccrs.Robinson(), facecolor='none', edgecolor='gray'))

The scale of features may be controlled using the scale plot option, the most common options being '10m', '50m' and '110m'. Cartopy will downloaded the requested resolution as needed.

gv.output(backend='bokeh')

(gf.ocean * gf.land.options(scale='110m', global_extent=True) * gv.Feature(graticules, group='Lines') + 
 gf.ocean * gf.land.options(scale='50m', global_extent=True) * gv.Feature(graticules, group='Lines'))

Zoom in using the bokeh zoom widget and you should see that the right hand panel is using a higher resolution dataset for the land feature.

Instead of displaying a Feature directly it is also possible to request the geometries inside a Feature using the Feature.geoms method, which also allows specifying a scale and a bounds to select a subregion:

gf.land.geoms('50m', bounds=(-10, 40, 10, 60))

When working interactively with higher resolution datasets it is sometimes necessary to dynamically update the geometries based on the current viewport. The resample_geometry operation is an efficient way to display only polygons that intersect with the current viewport and downsample polygons on-the-fly.

gv.operation.resample_geometry(gf.coastline.geoms('10m')).opts(width=400, height=400, color='black')

Try zooming into the plot above and you will see the coastline geometry resolve to a higher resolution dynamically (this requires a live Python kernel).

Shape

The gv.Shape object wraps around any shapely geometry, allowing finer grained control over each polygon. We can, for example, access the geometries on the LAND feature and display them individually. Here we will get the geometry corresponding to the Australian continent and display it using shapely’s inbuilt SVG repr (not yet a HoloViews plot, just a bare SVG displayed by Jupyter directly):

land_geoms = gf.land.geoms(as_element=False)
land_geoms[21]

Instead of letting shapely render it as an SVG, we can now wrap it in the gv.Shape object and let matplotlib or bokeh render it, alone or with other GeoViews or HoloViews objects:

australia = gv.Shape(land_geoms[21])
alice_springs = gv.Text(133.870,-21.5, 'Alice Springs')

australia * gv.Points([(133.870,-23.700)]).opts(color='black', width=400) * alice_springs

We can also supply a list of geometries directly to a Polygons or Path element:

gv.Polygons(land_geoms) + gv.Path(land_geoms)

This makes it possible to create choropleth maps, where each part of the geometry is assigned a value that will be used to color it. However, constructing a choropleth by combining a bunch of shapes can be a lot of effort and is error prone. For that reason, the Shape Element provides convenience methods to load geometries from a shapefile. Here we load the boundaries of UK electoral districts directly from an existing shapefile:

hv.output(backend='matplotlib')

shapefile = '../assets/boundaries/boundaries.shp'
gv.Shape.from_shapefile(shapefile, crs=ccrs.PlateCarree())