This notebook demonstrates how to use the PyGMT plotting maps in GPlately.

The PyGMT integration is still in early stage, and we are working on adding more features to it. Please open GitHub issues and let us know how we could enhance the PyGMT integration. We are grateful for your feedbacks. Thank you. https://github.com/GPlates/gplately/issues

import warnings, os
from pathlib import Path
import pygmt
from gplately.auxiliary import get_gplot, get_pygmt_basemap_figure
from gplately.plot.pygmt_plot import PygmtPlotEngine
from gplately import Raster
from plate_model_manager import PresentDayRasterManager

warnings.filterwarnings("ignore", category=UserWarning, module="gplately")
warnings.filterwarnings("ignore", category=FutureWarning)
warnings.filterwarnings("ignore", category=RuntimeWarning)

reconstruction_time = 55
data_dir = Path("./gplately-example-data")
# tell PlotTopologies object to use the PygmtPlotEngine
gplot = get_gplot(
    model="Muller2025",
    model_repo_dir="plate-model-repo",
    time=reconstruction_time,
    plot_engine=PygmtPlotEngine(),
)

# reconstruct the topography grid for the specified reconstruction time,
# and use it to create an illumination grid for shading
topo_file = PresentDayRasterManager(data_dir=str(data_dir)).get_raster("topography")
topo_grid = Raster(
    data=topo_file, plate_reconstruction=gplot.plate_reconstruction
).reconstruct(time=reconstruction_time)

# download age grid CPT file from https://raw.githubusercontent.com/GPlates/gplately/refs/heads/master/tests-dir/unittest/create-age-grids-video/agegrid.cpt
cpt_file = f"{data_dir}/agegrid.cpt"
if not os.path.isfile(cpt_file):
    import urllib.request

    urllib.request.urlretrieve(
        "https://raw.githubusercontent.com/GPlates/gplately/refs/heads/master/tests-dir/unittest/create-age-grids-video/agegrid.cpt",
        cpt_file,
    )

# you need to know how to specify projection and region in GMT way
fig = get_pygmt_basemap_figure(
    projection="N180/10c", region="d", frame=["xafg30", "yafg30"]
)

gplot.plot_grid(
    fig,
    "AgeGrids",
    cmap=cpt_file,
    nan_transparent=True,
    # shading=True,
    # shading="+a315+ne0.6",
    shading=topo_grid.to_data_array(),
)

# now you can plot some features with the PlotTopologies object
gplot.plot_coastlines(
    fig,
    edgecolor="none",
    facecolor="gray",
    linewidth=0.1,
    central_meridian=180,
    gmtlabel="Coastlines",
)

gplot.plot_all_topological_sections(
    fig,
    plot_subduction_teeth=True,
    other_kwargs={
        "color": "black",
        "linewidth": 0.5,
        "gmtlabel": "Miscellaneous",
    },
    ridge_kwargs={"color": "red", "linewidth": 0.7, "gmtlabel": "Ridges"},
    transform_kwargs={"color": "green", "linewidth": 0.7, "gmtlabel": "Transforms"},
    trench_kwargs={
        "color": "blue",
        "linewidth": 0.7,
        "gmtlabel": "Subduction Zones",
    },
)

gplot.plot_plate_motion_vectors(fig, normalise=True, color="red")

# use pygmt directly to plot title and legend
fig.text(
    text=f"{reconstruction_time} Ma (Muller2025)",
    position="TC",
    no_clip=True,
    font="12p,Helvetica,black",
    offset="j0/-0.5c",
)
with pygmt.config(FONT_ANNOT_PRIMARY=4):
    fig.legend(
        position="jBL+o-1.0/0",
        box="+gwhite+p0.25p",  # pyright: ignore[reportArgumentType]
    )

out_f = f"{data_dir}/04-PlotWithPyGMT.pdf"
fig.savefig(out_f, crop="+m0.4c")  # pyright: ignore[reportArgumentType]
print(f"The map has been saved to {out_f}.")

fig.show(width=1200, crop="+m0.4c")

The map has been saved to gplately-example-data/04-PlotWithPyGMT.pdf.