

In [ ]:

    
%matplotlib inline
%config InlineBackend.figure_format = "retina"

from matplotlib import rcParams

rcParams["savefig.dpi"] = 200
# rcParams["font.size"] = 8

import warnings

warnings.filterwarnings("ignore")

# turn off pandas html repr:
# does not gracefully survive the ipynb -> rst -> html conversion
import pandas as pd
pd.set_option("display.notebook_repr_html", False)

Working with geopandas (shapefiles)

regionmask includes support for regions defined as geopandas GeoDataFrame. These are often shapefiles, which can be opened in the formats .zip or .shp with geopandas.read_file(url_or_path).

There are two possibilities:

Directly create a mask from a geopandas GeoDataFrame or GeoSeries using mask_geopandas or mask_3D_geopandas.
Convert a GeoDataFrame to a Regions object (regionmask's internal data container) using from_geopandas.

As always, start with the imports:



In [ ]:

    
import cartopy.crs as ccrs
import geopandas as gp
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import pooch

import regionmask

regionmask.__version__

Opening an example shapefile

The U.S. Geological Survey (USGS) offers a shapefile containing the outlines of continens [1]. We use the library pooch to locally cache the file:



In [ ]:

    
file = pooch.retrieve(
    "https://pubs.usgs.gov/of/2006/1187/basemaps/continents/continents.zip", None
)

continents = gp.read_file("zip://" + file)

display(continents)

Create a mask from a GeoDataFrame

mask_geopandas and mask_3D_geopandas allow to directly create a mask from a GeoDataFrame or GeoSeries:



In [ ]:

    
lon = np.arange(-180, 180)
lat = np.arange(-90, 90)

mask = regionmask.mask_geopandas(continents, lon, lat)

Let's plot the new mask:



In [ ]:

    
f, ax = plt.subplots(subplot_kw=dict(projection=ccrs.PlateCarree()))
mask.plot(
    ax=ax, transform=ccrs.PlateCarree(), add_colorbar=False,
)

ax.coastlines(color="0.1");

Similarly a 3D boolean mask can be created from a GeoDataFrame:



In [ ]:

    
mask_3D = regionmask.mask_3D_geopandas(continents, lon, lat)

and plotted:



In [ ]:

    
from matplotlib import colors as mplc

cmap1 = mplc.ListedColormap(["none", "#9ecae1"])

f, ax = plt.subplots(subplot_kw=dict(projection=ccrs.PlateCarree()))

mask_3D.sel(region=0).plot(
    ax=ax, transform=ccrs.PlateCarree(), add_colorbar=False, cmap=cmap1,
)

ax.coastlines(color="0.1");

2. Convert GeoDataFrame to a Regions object

Creating a Regions object with regionmask.from_geopandas requires a GeoDataFrame:



In [ ]:

    
continents_regions = regionmask.from_geopandas(continents)
continents_regions

This creates default names ("Region0", ..., "RegionN") and abbreviations ("r0", ..., "rN").

However, it is often advantageous to use columns of the GeoDataFrame as names and abbrevs. If no column with abbreviations is available, you can use abbrevs='_from_name', which creates unique abbreviations using the names column.



In [ ]:

    
continents_regions = regionmask.from_geopandas(
    continents, names="CONTINENT", abbrevs="_from_name", name="continent"
)
continents_regions

As usual the newly created Regions object can be plotted on a world map:



In [ ]:

    
continents_regions.plot(label="name", coastlines=False);

And to create mask a mask for arbitrary latitude/ longitude grids:



In [ ]:

    
lon = np.arange(0, 360)
lat = np.arange(-90, 90)

mask = continents_regions.mask(lon, lat)

which can then be plotted



In [ ]:

    
f, ax = plt.subplots(subplot_kw=dict(projection=ccrs.PlateCarree()))
h = mask.plot(
    ax=ax,
    transform=ccrs.PlateCarree(),
    cmap="Reds",
    add_colorbar=False,
    levels=np.arange(-0.5, 8),
)

cbar = plt.colorbar(h, shrink=0.625, pad=0.025, aspect=12)
cbar.set_ticks(np.arange(8))
cbar.set_ticklabels(continents_regions.names)

ax.coastlines(color="0.2")

continents_regions.plot_regions(add_label=False);

References

[1] Environmental Systems Research , Inc. (ESRI), 20020401, World Continents: ESRI Data & Maps 2002, Environmental Systems Research Institute, Inc. (ESRI), Redlands, California, USA.