Use PyGMT instead of Cartopy in Chain tutorial (#386)

The PyGMT maps require less code and tend to look nicer than the Cartopy
ones. It's also faster for non-rectangular maps. Add PyGMT and GMT as
docs dependencies with pinned versions to avoid breaking the docs build.
Replacing in a single tutorial to establish the infrastructure. Replace
in other places in follow-ups.

Co-authored-by: Santiago Soler <santisoler@fastmail.com>

doc/Makefile

@@ -40,7 +40,8 @@ html: api
 	@echo
 	@echo "Building HTML files."
 	@echo
-	$(SPHINXBUILD) -b html $(ALLSPHINXOPTS) $(BUILDDIR)/html
+	# Set PYGMT_USE_EXTERNAL_DISPLAY to "false" to disable external display
+	PYGMT_USE_EXTERNAL_DISPLAY="false" $(SPHINXBUILD) -b html $(ALLSPHINXOPTS) $(BUILDDIR)/html
 	@echo
 	@echo "Build finished. The HTML pages are in $(BUILDDIR)/html."
 

doc/conf.py

@@ -6,6 +6,7 @@
 #
 import datetime
 
+import pygmt.sphinx_gallery
 from sphinx_gallery.sorting import FileNameSortKey
 
 import verde
@@ -50,6 +51,7 @@
     "pooch": ("https://www.fatiando.org/pooch/latest/", None),
     "matplotlib": ("https://matplotlib.org/", None),
     "dask": ("https://docs.dask.org/en/latest/", None),
+    "pygmt": ("https://www.pygmt.org/latest/", None),
 }
 
 # Autosummary pages will be generated by sphinx-autogen instead of sphinx-build
@@ -93,6 +95,8 @@
     "doc_module": "verde",
     # Insert links to documentation of objects in the examples
     "reference_url": {"verde": None},
+    # Use the PyGMT image scraper
+    "image_scrapers": (pygmt.sphinx_gallery.PyGMTScraper(),),
 }
 
 # HTML output configuration

doc/install.rst

@@ -77,5 +77,6 @@ If you wish to **run the examples in the documentation**, you will also have to
 install:
 
 * `matplotlib <https://matplotlib.org/>`__
+* `pygmt <https://www.pygmt.org>`__ for plotting maps
 * `cartopy <https://scitools.org.uk/cartopy/>`__ for plotting maps
 * `pyproj <https://jswhit.github.io/pyproj/>`__ for cartographic projections

doc/tutorials_src/chain.py

@@ -36,9 +36,9 @@
 
 For example, let's create a pipeline to grid our sample bathymetry data.
 """
-import cartopy.crs as ccrs
 import matplotlib.pyplot as plt
 import numpy as np
+import pygmt
 import pyproj
 
 import verde as vd
@@ -51,19 +51,20 @@
 projection = pyproj.Proj(proj="merc", lat_ts=data.latitude.mean())
 proj_coords = projection(data.longitude.values, data.latitude.values)
 
-plt.figure(figsize=(7, 6))
-ax = plt.axes(projection=ccrs.Mercator())
-ax.set_title("Bathymetry from Baja California")
-plt.scatter(
-    data.longitude,
-    data.latitude,
-    c=data.bathymetry_m,
-    s=0.1,
-    transform=ccrs.PlateCarree(),
+fig = pygmt.Figure()
+fig.coast(region=region, projection="M20c", land="#666666")
+pygmt.makecpt(cmap="viridis", series=[data.bathymetry_m.min(), data.bathymetry_m.max()])
+fig.plot(
+    x=data.longitude,
+    y=data.latitude,
+    color=data.bathymetry_m,
+    cmap=True,
+    style="c0.05c",
 )
-plt.colorbar().set_label("meters")
-vd.datasets.setup_baja_bathymetry_map(ax)
-plt.show()
+fig.basemap(frame=True)
+fig.colorbar(frame='af+l"bathymetric depth [m]"')
+fig.show()
+
 
 ###############################################################################
 # We'll create a chain that applies a blocked median to the data, fits a
@@ -120,20 +121,27 @@
     dims=["latitude", "longitude"],
     data_names="bathymetry",
 )
-print(grid)
+grid = vd.distance_mask(
+    data_coordinates=(data.longitude, data.latitude),
+    maxdist=spacing * 111e3,
+    grid=grid,
+    projection=projection,
+)
+grid
 
 ###############################################################################
 # Finally, we can plot the resulting grid:
 
-plt.figure(figsize=(7, 6))
-ax = plt.axes(projection=ccrs.Mercator())
-ax.set_title("Gridded result of the chain")
-pc = grid.bathymetry.plot.pcolormesh(
-    ax=ax, transform=ccrs.PlateCarree(), vmax=0, zorder=-1, add_colorbar=False
+fig = pygmt.Figure()
+fig.coast(region=region, projection="M20c", land="#666666")
+fig.grdimage(
+    grid=grid.bathymetry,
+    cmap="viridis",
+    nan_transparent=True,
 )
-plt.colorbar(pc).set_label("meters")
-vd.datasets.setup_baja_bathymetry_map(ax)
-plt.show()
+fig.basemap(frame=True)
+fig.colorbar(frame='af+l"bathymetric depth [m]"')
+fig.show()
 
 ###############################################################################
 # Each component of the chain can be accessed separately using the
@@ -155,14 +163,24 @@
     dims=["latitude", "longitude"],
     data_names="bathymetry",
 )
-print(grid_trend)
+grid_trend = vd.distance_mask(
+    data_coordinates=(data.longitude, data.latitude),
+    maxdist=spacing * 111e3,
+    grid=grid_trend,
+    projection=projection,
+)
+grid_trend
 
-plt.figure(figsize=(7, 6))
-ax = plt.axes(projection=ccrs.Mercator())
-ax.set_title("Gridded trend")
-pc = grid_trend.bathymetry.plot.pcolormesh(
-    ax=ax, transform=ccrs.PlateCarree(), zorder=-1, add_colorbar=False
+###############################################################################
+#
+
+fig = pygmt.Figure()
+fig.coast(region=region, projection="M20c", land="#666666")
+fig.grdimage(
+    grid=grid_trend.bathymetry,
+    cmap="viridis",
+    nan_transparent=True,
 )
-plt.colorbar(pc).set_label("meters")
-vd.datasets.setup_baja_bathymetry_map(ax)
-plt.show()
+fig.basemap(frame=True)
+fig.colorbar(frame='af+l"bathymetric depth [m]"')
+fig.show()

env/requirements-docs.txt

@@ -3,6 +3,9 @@ sphinx-book-theme==0.3.*
 sphinx-gallery==0.10.*
 sphinx-copybutton==0.5.*
 sphinx-design==0.1.*
-matplotlib
-cartopy>=0.18
+matplotlib==3.5.*
+cartopy>=0.20
 pyproj
+pygmt==0.6.*
+gmt==6.3.*
+ipython

environment.yml

@@ -20,8 +20,8 @@ dependencies:
     - pykdtree
     - numba
     # Test
-    - matplotlib
-    - cartopy>=0.18
+    - matplotlib==3.5.*
+    - cartopy>=0.20
     - pytest
     - pytest-cov
     - pytest-mpl
@@ -33,6 +33,9 @@ dependencies:
     - sphinx-copybutton==0.5.*
     - sphinx-design==0.1.*
     - pyproj
+    - pygmt==0.6.*
+    - gmt==6.3.*
+    - ipython
     # Style
     - black
     - pathspec