From b23e18140a3971c2af62d22ecba13a93df7ce1e8 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Thu, 8 Aug 2024 17:35:01 +1000 Subject: [PATCH 01/14] add a new folder plotting --- gplately/plotting/pygmt_plotter.py | 0 gplately/plotting/reconstructor.py | 0 2 files changed, 0 insertions(+), 0 deletions(-) create mode 100644 gplately/plotting/pygmt_plotter.py create mode 100644 gplately/plotting/reconstructor.py diff --git a/gplately/plotting/pygmt_plotter.py b/gplately/plotting/pygmt_plotter.py new file mode 100644 index 00000000..e69de29b diff --git a/gplately/plotting/reconstructor.py b/gplately/plotting/reconstructor.py new file mode 100644 index 00000000..e69de29b From 642be2d5e1867898b18012b9f11fe79781cf9052 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Thu, 17 Oct 2024 19:29:42 +1100 Subject: [PATCH 02/14] initial effort to add pygmt plot --- docker/env.yaml | 2 + gplately/mapping/plot_engine.py | 23 ++++++ gplately/mapping/pygmt_plot.py | 81 +++++++++++++++++++ gplately/plot.py | 16 +++- pyproject.toml | 1 + setup.py | 1 + tests-dir/unittest/test_plot.py | 4 + tests-dir/unittest/test_pygmt_plot.py | 111 ++++++++++++++++++++++++++ 8 files changed, 235 insertions(+), 4 deletions(-) create mode 100644 gplately/mapping/plot_engine.py create mode 100644 gplately/mapping/pygmt_plot.py create mode 100644 tests-dir/unittest/test_pygmt_plot.py diff --git a/docker/env.yaml b/docker/env.yaml index 9157740c..3e9ab660 100644 --- a/docker/env.yaml +++ b/docker/env.yaml @@ -7,3 +7,5 @@ dependencies: - gplately - jupyter - moviepy + - gmt + - pygmt diff --git a/gplately/mapping/plot_engine.py b/gplately/mapping/plot_engine.py new file mode 100644 index 00000000..77af7ecc --- /dev/null +++ b/gplately/mapping/plot_engine.py @@ -0,0 +1,23 @@ +# +# Copyright (C) 2024 The University of Sydney, Australia +# +# This program is free software; you can redistribute it and/or modify it under +# the terms of the GNU General Public License, version 2, as published by +# the Free Software Foundation. +# +# This program is distributed in the hope that it will be useful, but WITHOUT +# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +# for more details. +# +# You should have received a copy of the GNU General Public License along +# with this program; if not, write to Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. +# + +from enum import Enum + + +class PlotEngine(Enum): + CARTOPY = 1 + PYGMT = 2 diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py new file mode 100644 index 00000000..21614e6f --- /dev/null +++ b/gplately/mapping/pygmt_plot.py @@ -0,0 +1,81 @@ +# +# Copyright (C) 2024 The University of Sydney, Australia +# +# This program is free software; you can redistribute it and/or modify it under +# the terms of the GNU General Public License, version 2, as published by +# the Free Software Foundation. +# +# This program is distributed in the hope that it will be useful, but WITHOUT +# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +# for more details. +# +# You should have received a copy of the GNU General Public License along +# with this program; if not, write to Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. +# +from geopandas.geodataframe import GeoDataFrame +import pygmt + +# ----- parameters for plot +region = "d" +width = 10 +projection = "N180/" +x_offset = width + 2 + +# plate boundary stuff +plateboundary_width = "0.5p" + +age_font = "12p,Helvetica,black" + +label_font = "12p,Helvetica,black" +label_offset = "j0/-0.5c" +label_position = "TC" + + +def plot_geo_data_frame(gdf: GeoDataFrame, **kwargs): + + print(gdf) + fig = pygmt.Figure() + pygmt.config( + FONT_ANNOT=8, + FONT_LABEL=8, + FONT=8, + MAP_TICK_PEN="0.75p", + MAP_FRAME_PEN="0.75p", + MAP_TICK_LENGTH_PRIMARY="4p", + ) + + # ---- part a - transforms FROM gplot.get_transforms() + fig.basemap(region=region, projection="%s%sc" % (projection, width), frame="lrtb") + # fig.coast(shorelines=True) + fig.plot(data=gdf.geometry, pen="0.5p,blue") + """ + fig.plot(data=gdf_coastlines, fill=coastline_color, frame=["xa0", "ya0"], transparency=0) + + fig.plot(data=gdf_topo_plates.geometry, pen='%s,%s' % (plateboundary_width, plate_colour), frame="lrtb") + fig.plot(data=gdf_subduction_left, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+l+t') + fig.plot(data=gdf_subduction_right, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+r+t') + fig.plot(data=gdf_ridges_transforms, pen='%s,%s' % (plateboundary_width, ridge_colour)) + fig.plot(data=gplot.get_transforms(), pen='%s,%s' % (plateboundary_width, transform_color)) + + fig.text(text='gplot.get_transforms(): %s Ma' % age, position=label_position, no_clip=True, font=label_font, offset=label_offset) + + fig.shift_origin(xshift=x_offset) + fig.basemap(region=region, projection="%s%sc" % (projection, width), frame="lrtb") + fig.plot(data=gdf_cobs, fill=COB_color, transparency=0, ) + fig.plot(data=gdf_coastlines, fill=coastline_color, frame=["xa0", "ya0"], transparency=0) + + fig.plot(data=gdf_topo_plates.geometry, pen='%s,%s' % (plateboundary_width, plate_colour), frame="lrtb", label='other plate boundary types') + fig.plot(data=gdf_subduction_left, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+l+t', label='subduction zones') + fig.plot(data=gdf_subduction_right, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+r+t') + fig.plot(data=gdf_ridges_transforms, pen='%s,%s' % (plateboundary_width, ridge_colour), label='ridges and transforms') + + # from gpml: transforms + fig.plot(data=gdf_topo_transforms, pen='%s,%s' % (plateboundary_width, transform_color), label = 'transforms') + fig.text(text='FeatureType.gpml_transform: %s Ma' % age, position=label_position, no_clip=True, font=label_font, offset=label_offset) + + fig.legend(position='jBL+o-2.7/0', box="+gwhite+p0.5p") + """ + # fig.show(width=1200) + fig.savefig("test-pygmt-plot.pdf") diff --git a/gplately/plot.py b/gplately/plot.py index b50e1b96..ebcbee16 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -51,6 +51,8 @@ from .utils.feature_utils import shapelify_features as _shapelify_features from .utils.plot_utils import _clean_polygons, _meridian_from_ax from .utils.plot_utils import plot_subduction_teeth as _plot_subduction_teeth +from .mapping.plot_engine import PlotEngine +from .mapping.pygmt_plot import plot_geo_data_frame logger = logging.getLogger("gplately") @@ -298,6 +300,7 @@ def __init__( COBs=None, time=None, anchor_plate_id=0, + plot_engine: PlotEngine = PlotEngine.CARTOPY, ): self.plate_reconstruction = plate_reconstruction @@ -321,6 +324,7 @@ def __init__( self._topologies = None self._anchor_plate_id = self._check_anchor_plate_id(anchor_plate_id) + self._plot_engine = plot_engine # store topologies for easy access # setting time runs the update_time routine @@ -726,12 +730,16 @@ def _plot_feature(self, ax, get_feature_func, **kwargs): logger.warning("No feature found for plotting. Do nothing and return.") return ax - if hasattr(ax, "projection"): - gdf = _clean_polygons(data=gdf, projection=ax.projection) + if self._plot_engine == PlotEngine.PYGMT: + kwargs["ax"] = ax + return plot_geo_data_frame(gdf, **kwargs) else: - kwargs["transform"] = self.base_projection + if hasattr(ax, "projection"): + gdf = _clean_polygons(data=gdf, projection=ax.projection) + else: + kwargs["transform"] = self.base_projection - return gdf.plot(ax=ax, **kwargs) + return gdf.plot(ax=ax, **kwargs) @validate_reconstruction_time @append_docstring(GET_DATE_DOCSTRING.format("coastlines")) diff --git a/pyproject.toml b/pyproject.toml index 3cd2956c..c9319487 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -35,6 +35,7 @@ dependencies = [ "netcdf4", "rasterio", "geopandas", + "gmt", "stripy", "plate-model-manager>=1.2.0", "pyyaml", diff --git a/setup.py b/setup.py index cbf50e70..2eab0371 100644 --- a/setup.py +++ b/setup.py @@ -96,6 +96,7 @@ def _minimal_ext_cmd(cmd): "netcdf4", "rasterio", "geopandas", + "gmt", "stripy", "plate-model-manager", "pyyaml", diff --git a/tests-dir/unittest/test_plot.py b/tests-dir/unittest/test_plot.py index 4132b5aa..3530a284 100755 --- a/tests-dir/unittest/test_plot.py +++ b/tests-dir/unittest/test_plot.py @@ -1,4 +1,7 @@ #!/usr/bin/env python3 +# import matplotlib + +# matplotlib.use("QtAgg") import sys @@ -11,6 +14,7 @@ import gplately from gplately import PlateReconstruction, PlotTopologies + print(gplately.__file__) # test the plot function with the new PlateModel class diff --git a/tests-dir/unittest/test_pygmt_plot.py b/tests-dir/unittest/test_pygmt_plot.py new file mode 100644 index 00000000..13aeea62 --- /dev/null +++ b/tests-dir/unittest/test_pygmt_plot.py @@ -0,0 +1,111 @@ +#!/usr/bin/env python3 + +import sys + +import cartopy.crs as ccrs +import matplotlib.pyplot as plt +import numpy as np +from common import MODEL_REPO_DIR, save_fig +from plate_model_manager import PlateModel, PlateModelManager + +import gplately +from gplately import PlateReconstruction, PlotTopologies +from gplately.mapping.plot_engine import PlotEngine + +print(gplately.__file__) + +# test the plot function with the new PlateModel class + +# MODEL_NAME = "Clennett2020" +# MODEL_NAME = "Muller2019" +MODEL_NAME = "merdith2021" + + +def main(show=True): + try: + model = PlateModelManager().get_model(MODEL_NAME, data_dir=MODEL_REPO_DIR) + except: + model = PlateModel(MODEL_NAME, data_dir=MODEL_REPO_DIR, readonly=True) + + age = 55 + + test_model = PlateReconstruction( + model.get_rotation_model(), + topology_features=model.get_layer("Topologies"), + static_polygons=model.get_layer("StaticPolygons"), + ) + gplot = PlotTopologies( + test_model, + coastlines=model.get_layer("Coastlines"), + COBs=model.get_layer("COBs", return_none_if_not_exist=True), + continents=model.get_layer("ContinentalPolygons"), + time=age, + plot_engine=PlotEngine.PYGMT, + ) + + # age = 100 + # gplot.time = age + + ax = None + + all_flag = 0 + plot_flag = { + "continent_ocean_boundaries": 0, + "coastlines": 0, + "trenches": 0, + "subduction_teeth": 0, + "ridges": 0, + "all_topologies": 0, + "all_topological_sections": 0, + "plate_polygon_by_id": 0, + "unclassified_features": 0, + "slab_edges": 0, + "passive_continental_boundaries": 0, + "extended_continental_crusts": 0, + "continental_crusts": 0, + "sutures": 0, + "orogenic_belts": 0, + "transitional_crusts": 0, + "terrane_boundaries": 0, + "inferred_paleo_boundaries": 0, + "fracture_zones": 0, + "faults": 0, + "continental_rifts": 0, + "misc_boundaries": 0, + "transforms": 1, + "continents": 0, + "topological_plate_boundaries": 0, + } + + for key in plot_flag: + if key == "plate_polygon_by_id": + continue + if key == "continents": + if plot_flag["continents"]: + gplot.plot_continents(ax, color="grey", facecolor="0.8") + continue + if key == "coastlines": + if plot_flag["coastlines"]: + gplot.plot_coastlines(ax, edgecolor="blue", facecolor="0.5") + continue + if all_flag or plot_flag[key]: + getattr(gplot, f"plot_{key}")( + ax, color=list(np.random.choice(range(256), size=3) / 256) + ) + + ids = set([f.get_reconstruction_plate_id() for f in gplot.topologies]) + for id in ids: + if all_flag or plot_flag["plate_polygon_by_id"]: + gplot.plot_plate_polygon_by_id( + ax, + id, + facecolor="None", + edgecolor=list(np.random.choice(range(256), size=3) / 256), + ) + + +if __name__ == "__main__": + if len(sys.argv) == 2 and sys.argv[1] == "save": + main(show=False) + else: + main(show=True) From a2a15aae2ddf3e9e028784b5a8671d51e0c8b954 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Fri, 18 Oct 2024 19:16:16 +1100 Subject: [PATCH 03/14] add a new notebook to test pygmt --- gplately/mapping/pygmt_plot.py | 19 +---- gplately/plot.py | 3 +- tests-dir/unittest/test_pygmt.ipynb | 101 ++++++++++++++++++++++++++ tests-dir/unittest/test_pygmt_plot.py | 22 +++--- 4 files changed, 116 insertions(+), 29 deletions(-) create mode 100644 tests-dir/unittest/test_pygmt.ipynb diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py index 21614e6f..821d7b9b 100644 --- a/gplately/mapping/pygmt_plot.py +++ b/gplately/mapping/pygmt_plot.py @@ -33,22 +33,7 @@ label_position = "TC" -def plot_geo_data_frame(gdf: GeoDataFrame, **kwargs): - - print(gdf) - fig = pygmt.Figure() - pygmt.config( - FONT_ANNOT=8, - FONT_LABEL=8, - FONT=8, - MAP_TICK_PEN="0.75p", - MAP_FRAME_PEN="0.75p", - MAP_TICK_LENGTH_PRIMARY="4p", - ) - - # ---- part a - transforms FROM gplot.get_transforms() - fig.basemap(region=region, projection="%s%sc" % (projection, width), frame="lrtb") - # fig.coast(shorelines=True) +def plot_geo_data_frame(fig: pygmt.Figure, gdf: GeoDataFrame, **kwargs): fig.plot(data=gdf.geometry, pen="0.5p,blue") """ fig.plot(data=gdf_coastlines, fill=coastline_color, frame=["xa0", "ya0"], transparency=0) @@ -77,5 +62,3 @@ def plot_geo_data_frame(gdf: GeoDataFrame, **kwargs): fig.legend(position='jBL+o-2.7/0', box="+gwhite+p0.5p") """ - # fig.show(width=1200) - fig.savefig("test-pygmt-plot.pdf") diff --git a/gplately/plot.py b/gplately/plot.py index ebcbee16..1a0621e4 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -731,8 +731,7 @@ def _plot_feature(self, ax, get_feature_func, **kwargs): return ax if self._plot_engine == PlotEngine.PYGMT: - kwargs["ax"] = ax - return plot_geo_data_frame(gdf, **kwargs) + return plot_geo_data_frame(fig=ax, gdf=gdf, **kwargs) else: if hasattr(ax, "projection"): gdf = _clean_polygons(data=gdf, projection=ax.projection) diff --git a/tests-dir/unittest/test_pygmt.ipynb b/tests-dir/unittest/test_pygmt.ipynb new file mode 100644 index 00000000..78664937 --- /dev/null +++ b/tests-dir/unittest/test_pygmt.ipynb @@ -0,0 +1,101 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "id": "f5003633-363f-4d07-b857-879f9ba8726e", + "metadata": {}, + "outputs": [], + "source": [ + "import sys\n", + "\n", + "import cartopy.crs as ccrs\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "from common import MODEL_REPO_DIR, save_fig\n", + "from plate_model_manager import PlateModel, PlateModelManager\n", + "\n", + "import gplately\n", + "from gplately import PlateReconstruction, PlotTopologies\n", + "from gplately.mapping.plot_engine import PlotEngine\n", + "import pygmt\n", + "pygmt.config(\n", + " FONT_ANNOT=8,\n", + " FONT_LABEL=8,\n", + " FONT=8,\n", + " MAP_TICK_PEN=\"0.75p\",\n", + " MAP_FRAME_PEN=\"0.75p\",\n", + " MAP_TICK_LENGTH_PRIMARY=\"4p\",\n", + " )\n", + "\n", + "# MODEL_NAME = \"Clennett2020\"\n", + "# MODEL_NAME = \"Muller2019\"\n", + "MODEL_NAME = \"merdith2021\"\n", + "\n", + "try:\n", + " model = PlateModelManager().get_model(MODEL_NAME, data_dir=MODEL_REPO_DIR)\n", + "except:\n", + " model = PlateModel(MODEL_NAME, data_dir=MODEL_REPO_DIR, readonly=True)\n", + "\n", + "if model is None:\n", + " raise Exception(f\"Unable to get model ({MODEL_NAME})\")\n", + "\n", + "age = 55\n", + "\n", + "test_model = PlateReconstruction(\n", + " model.get_rotation_model(),\n", + " topology_features=model.get_layer(\"Topologies\"),\n", + " static_polygons=model.get_layer(\"StaticPolygons\"),\n", + ")\n", + "gplot = PlotTopologies(\n", + " test_model,\n", + " coastlines=model.get_layer(\"Coastlines\"),\n", + " COBs=model.get_layer(\"COBs\", return_none_if_not_exist=True),\n", + " continents=model.get_layer(\"ContinentalPolygons\"),\n", + " time=age,\n", + " plot_engine=PlotEngine.PYGMT,\n", + ")\n", + "\n", + "fig = pygmt.Figure()\n", + "fig.basemap(region=\"d\", projection=\"N180/10c\", frame=\"lrtb\")\n", + "#fig.coast(shorelines=True)\n", + "\n", + "gplot.plot_topological_plate_boundaries(fig)\n", + "gplot.plot_coastlines(fig)\n", + "\n", + "fig.show(width=1200)\n", + "# fig.savefig(\"test-pygmt-plot.pdf\")\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b25b75a3-39f1-4494-912e-7f9f502a67ce", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.9" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/tests-dir/unittest/test_pygmt_plot.py b/tests-dir/unittest/test_pygmt_plot.py index 13aeea62..19448ab6 100644 --- a/tests-dir/unittest/test_pygmt_plot.py +++ b/tests-dir/unittest/test_pygmt_plot.py @@ -27,6 +27,9 @@ def main(show=True): except: model = PlateModel(MODEL_NAME, data_dir=MODEL_REPO_DIR, readonly=True) + if model is None: + raise Exception(f"Unable to get model ({MODEL_NAME})") + age = 55 test_model = PlateReconstruction( @@ -93,15 +96,16 @@ def main(show=True): ax, color=list(np.random.choice(range(256), size=3) / 256) ) - ids = set([f.get_reconstruction_plate_id() for f in gplot.topologies]) - for id in ids: - if all_flag or plot_flag["plate_polygon_by_id"]: - gplot.plot_plate_polygon_by_id( - ax, - id, - facecolor="None", - edgecolor=list(np.random.choice(range(256), size=3) / 256), - ) + if gplot.topologies is not None: + ids = set([f.get_reconstruction_plate_id() for f in gplot.topologies]) + for id in ids: + if all_flag or plot_flag["plate_polygon_by_id"]: + gplot.plot_plate_polygon_by_id( + ax, + id, + facecolor="None", + edgecolor=list(np.random.choice(range(256), size=3) / 256), + ) if __name__ == "__main__": From e523c5a0608d206c0639ddc826b4e3912eaee54a Mon Sep 17 00:00:00 2001 From: michaelchin Date: Fri, 18 Oct 2024 19:34:06 +1100 Subject: [PATCH 04/14] add a new function get_gplot() --- gplately/utils/plot_utils.py | 50 ++++++++++++++++++++++++++++++++++++ 1 file changed, 50 insertions(+) diff --git a/gplately/utils/plot_utils.py b/gplately/utils/plot_utils.py index 80930128..5c3e1ea3 100644 --- a/gplately/utils/plot_utils.py +++ b/gplately/utils/plot_utils.py @@ -25,8 +25,13 @@ from shapely.geometry import LineString, MultiPolygon, Point, Polygon, box from shapely.geometry.base import BaseGeometry, BaseMultipartGeometry from shapely.ops import linemerge, substring +from plate_model_manager import PlateModel, PlateModelManager +from typing import Union from .io_utils import get_geometries as _get_geometries +from ..reconstruction import PlateReconstruction +from ..plot import PlotTopologies +from ..mapping.plot_engine import PlotEngine logger = logging.getLogger("gplately") @@ -632,3 +637,48 @@ def plot_subduction_teeth( else: for triangle in triangles: ax.fill(*triangle.exterior.xy, **kwargs) + + +def get_gplot( + model_name: str, model_repo_dir: str, age: Union[int, float] +) -> PlotTopologies: + try: + model = PlateModelManager().get_model(model_name, data_dir=model_repo_dir) + except: + model = PlateModel(model_name, data_dir=model_repo_dir, readonly=True) + + if model is None: + raise Exception(f"Unable to get model ({model_name})") + + topology_features = None + static_polygons = None + coastlines = None + COBs = None + continents = None + + all_layers = model.get_avail_layers() + + if "Topologies" in all_layers: + topology_features = model.get_layer("Topologies") + if "StaticPolygons" in all_layers: + static_polygons = model.get_layer("StaticPolygons") + if "Coastlines" in all_layers: + coastlines = model.get_layer("Coastlines") + if "COBs" in all_layers: + COBs = model.get_layer("COBs") + if "ContinentalPolygons" in all_layers: + continents = (model.get_layer("ContinentalPolygons"),) + + m = PlateReconstruction( + model.get_rotation_model(), + topology_features=topology_features, + static_polygons=static_polygons, + ) + return PlotTopologies( + m, + coastlines=coastlines, + COBs=COBs, + continents=continents, + time=age, + plot_engine=PlotEngine.PYGMT, + ) From 006773c349af4af89f7c5219947e5543ca6e6402 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Fri, 18 Oct 2024 19:41:47 +1100 Subject: [PATCH 05/14] some unfinished code --- gplately/__init__.py | 2 ++ gplately/utils/plot_utils.py | 1 + tests-dir/unittest/test_pygmt.ipynb | 6 ++++-- 3 files changed, 7 insertions(+), 2 deletions(-) diff --git a/gplately/__init__.py b/gplately/__init__.py index 9b5b399a..99ddabda 100644 --- a/gplately/__init__.py +++ b/gplately/__init__.py @@ -248,6 +248,7 @@ from .tools import EARTH_RADIUS from .utils import io_utils from .utils.io_utils import get_geometries, get_valid_geometries +from .utils.plot_utils import get_gplot __pdoc__ = { "data": False, @@ -267,6 +268,7 @@ "data", "download", "geometry", + "get_gplot", "gpml", "grids", "oceans", diff --git a/gplately/utils/plot_utils.py b/gplately/utils/plot_utils.py index 5c3e1ea3..1619708d 100644 --- a/gplately/utils/plot_utils.py +++ b/gplately/utils/plot_utils.py @@ -642,6 +642,7 @@ def plot_subduction_teeth( def get_gplot( model_name: str, model_repo_dir: str, age: Union[int, float] ) -> PlotTopologies: + """convenient function to get gplot object""" try: model = PlateModelManager().get_model(model_name, data_dir=model_repo_dir) except: diff --git a/tests-dir/unittest/test_pygmt.ipynb b/tests-dir/unittest/test_pygmt.ipynb index 78664937..6e565252 100644 --- a/tests-dir/unittest/test_pygmt.ipynb +++ b/tests-dir/unittest/test_pygmt.ipynb @@ -15,9 +15,9 @@ "from common import MODEL_REPO_DIR, save_fig\n", "from plate_model_manager import PlateModel, PlateModelManager\n", "\n", - "import gplately\n", "from gplately import PlateReconstruction, PlotTopologies\n", "from gplately.mapping.plot_engine import PlotEngine\n", + "from gplately import get_gplot\n", "import pygmt\n", "pygmt.config(\n", " FONT_ANNOT=8,\n", @@ -32,6 +32,7 @@ "# MODEL_NAME = \"Muller2019\"\n", "MODEL_NAME = \"merdith2021\"\n", "\n", + "\"\"\"\n", "try:\n", " model = PlateModelManager().get_model(MODEL_NAME, data_dir=MODEL_REPO_DIR)\n", "except:\n", @@ -55,7 +56,8 @@ " time=age,\n", " plot_engine=PlotEngine.PYGMT,\n", ")\n", - "\n", + "\"\"\"\n", + "gplot = get_gplot(MODEL_NAME,MODEL_REPO_DIR,55)\n", "fig = pygmt.Figure()\n", "fig.basemap(region=\"d\", projection=\"N180/10c\", frame=\"lrtb\")\n", "#fig.coast(shorelines=True)\n", From 0bbcd213880d2d1e5152b35ce34ecaf9cdce3cd8 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Sat, 19 Oct 2024 15:06:53 +1100 Subject: [PATCH 06/14] finish the unfinished code yesterday --- gplately/__init__.py | 2 -- gplately/auxiliary.py | 56 +++++++++++++++++++++++++++++ gplately/plot.py | 6 ++-- gplately/utils/plot_utils.py | 51 -------------------------- pyproject.toml | 1 - tests-dir/unittest/test_pygmt.ipynb | 39 ++------------------ 6 files changed, 62 insertions(+), 93 deletions(-) create mode 100644 gplately/auxiliary.py diff --git a/gplately/__init__.py b/gplately/__init__.py index 99ddabda..9b5b399a 100644 --- a/gplately/__init__.py +++ b/gplately/__init__.py @@ -248,7 +248,6 @@ from .tools import EARTH_RADIUS from .utils import io_utils from .utils.io_utils import get_geometries, get_valid_geometries -from .utils.plot_utils import get_gplot __pdoc__ = { "data": False, @@ -268,7 +267,6 @@ "data", "download", "geometry", - "get_gplot", "gpml", "grids", "oceans", diff --git a/gplately/auxiliary.py b/gplately/auxiliary.py new file mode 100644 index 00000000..27c89983 --- /dev/null +++ b/gplately/auxiliary.py @@ -0,0 +1,56 @@ +from typing import Union + +from plate_model_manager import PlateModel, PlateModelManager + +from .mapping.plot_engine import PlotEngine +from .plot import PlotTopologies +from .reconstruction import PlateReconstruction + + +def get_gplot( + model_name: str, + model_repo_dir: str, + age: Union[int, float], + plot_engine: PlotEngine = PlotEngine.CARTOPY, +) -> PlotTopologies: + """auxiliary function to get gplot object""" + try: + model = PlateModelManager().get_model(model_name, data_dir=model_repo_dir) + except: + model = PlateModel(model_name, data_dir=model_repo_dir, readonly=True) + + if model is None: + raise Exception(f"Unable to get model ({model_name})") + + topology_features = None + static_polygons = None + coastlines = None + COBs = None + continents = None + + all_layers = model.get_avail_layers() + + if "Topologies" in all_layers: + topology_features = model.get_layer("Topologies") + if "StaticPolygons" in all_layers: + static_polygons = model.get_layer("StaticPolygons") + if "Coastlines" in all_layers: + coastlines = model.get_layer("Coastlines") + if "COBs" in all_layers: + COBs = model.get_layer("COBs") + if "ContinentalPolygons" in all_layers: + continents = model.get_layer("ContinentalPolygons") + + m = PlateReconstruction( + model.get_rotation_model(), + topology_features=topology_features, + static_polygons=static_polygons, + ) + return PlotTopologies( + m, + coastlines=coastlines, + COBs=COBs, + continents=continents, + time=age, + plot_engine=plot_engine, + ) diff --git a/gplately/plot.py b/gplately/plot.py index 1a0621e4..37992832 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -45,14 +45,14 @@ validate_topology_availability, ) from .gpml import _load_FeatureCollection +from .mapping.plot_engine import PlotEngine +from .mapping.pygmt_plot import plot_geo_data_frame from .pygplates import FeatureCollection as _FeatureCollection from .reconstruction import PlateReconstruction as _PlateReconstruction from .tools import EARTH_RADIUS from .utils.feature_utils import shapelify_features as _shapelify_features from .utils.plot_utils import _clean_polygons, _meridian_from_ax from .utils.plot_utils import plot_subduction_teeth as _plot_subduction_teeth -from .mapping.plot_engine import PlotEngine -from .mapping.pygmt_plot import plot_geo_data_frame logger = logging.getLogger("gplately") @@ -306,7 +306,7 @@ def __init__( if self.plate_reconstruction.topology_features is None: self.plate_reconstruction.topology_features = [] - logger.warn("Plate model does not have topology features.") + logger.warning("Plate model does not have topology features.") self.base_projection = ccrs.PlateCarree() diff --git a/gplately/utils/plot_utils.py b/gplately/utils/plot_utils.py index 1619708d..80930128 100644 --- a/gplately/utils/plot_utils.py +++ b/gplately/utils/plot_utils.py @@ -25,13 +25,8 @@ from shapely.geometry import LineString, MultiPolygon, Point, Polygon, box from shapely.geometry.base import BaseGeometry, BaseMultipartGeometry from shapely.ops import linemerge, substring -from plate_model_manager import PlateModel, PlateModelManager -from typing import Union from .io_utils import get_geometries as _get_geometries -from ..reconstruction import PlateReconstruction -from ..plot import PlotTopologies -from ..mapping.plot_engine import PlotEngine logger = logging.getLogger("gplately") @@ -637,49 +632,3 @@ def plot_subduction_teeth( else: for triangle in triangles: ax.fill(*triangle.exterior.xy, **kwargs) - - -def get_gplot( - model_name: str, model_repo_dir: str, age: Union[int, float] -) -> PlotTopologies: - """convenient function to get gplot object""" - try: - model = PlateModelManager().get_model(model_name, data_dir=model_repo_dir) - except: - model = PlateModel(model_name, data_dir=model_repo_dir, readonly=True) - - if model is None: - raise Exception(f"Unable to get model ({model_name})") - - topology_features = None - static_polygons = None - coastlines = None - COBs = None - continents = None - - all_layers = model.get_avail_layers() - - if "Topologies" in all_layers: - topology_features = model.get_layer("Topologies") - if "StaticPolygons" in all_layers: - static_polygons = model.get_layer("StaticPolygons") - if "Coastlines" in all_layers: - coastlines = model.get_layer("Coastlines") - if "COBs" in all_layers: - COBs = model.get_layer("COBs") - if "ContinentalPolygons" in all_layers: - continents = (model.get_layer("ContinentalPolygons"),) - - m = PlateReconstruction( - model.get_rotation_model(), - topology_features=topology_features, - static_polygons=static_polygons, - ) - return PlotTopologies( - m, - coastlines=coastlines, - COBs=COBs, - continents=continents, - time=age, - plot_engine=PlotEngine.PYGMT, - ) diff --git a/pyproject.toml b/pyproject.toml index 5f233c3e..7ca1623b 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -35,7 +35,6 @@ dependencies = [ "netcdf4", "rasterio", "geopandas", - "gmt", "stripy", "plate-model-manager>=1.2.1", "pyyaml", diff --git a/tests-dir/unittest/test_pygmt.ipynb b/tests-dir/unittest/test_pygmt.ipynb index 6e565252..6f9be13a 100644 --- a/tests-dir/unittest/test_pygmt.ipynb +++ b/tests-dir/unittest/test_pygmt.ipynb @@ -9,15 +9,11 @@ "source": [ "import sys\n", "\n", - "import cartopy.crs as ccrs\n", - "import matplotlib.pyplot as plt\n", - "import numpy as np\n", - "from common import MODEL_REPO_DIR, save_fig\n", "from plate_model_manager import PlateModel, PlateModelManager\n", "\n", "from gplately import PlateReconstruction, PlotTopologies\n", "from gplately.mapping.plot_engine import PlotEngine\n", - "from gplately import get_gplot\n", + "from gplately.auxiliary import get_gplot\n", "import pygmt\n", "pygmt.config(\n", " FONT_ANNOT=8,\n", @@ -28,36 +24,7 @@ " MAP_TICK_LENGTH_PRIMARY=\"4p\",\n", " )\n", "\n", - "# MODEL_NAME = \"Clennett2020\"\n", - "# MODEL_NAME = \"Muller2019\"\n", - "MODEL_NAME = \"merdith2021\"\n", - "\n", - "\"\"\"\n", - "try:\n", - " model = PlateModelManager().get_model(MODEL_NAME, data_dir=MODEL_REPO_DIR)\n", - "except:\n", - " model = PlateModel(MODEL_NAME, data_dir=MODEL_REPO_DIR, readonly=True)\n", - "\n", - "if model is None:\n", - " raise Exception(f\"Unable to get model ({MODEL_NAME})\")\n", - "\n", - "age = 55\n", - "\n", - "test_model = PlateReconstruction(\n", - " model.get_rotation_model(),\n", - " topology_features=model.get_layer(\"Topologies\"),\n", - " static_polygons=model.get_layer(\"StaticPolygons\"),\n", - ")\n", - "gplot = PlotTopologies(\n", - " test_model,\n", - " coastlines=model.get_layer(\"Coastlines\"),\n", - " COBs=model.get_layer(\"COBs\", return_none_if_not_exist=True),\n", - " continents=model.get_layer(\"ContinentalPolygons\"),\n", - " time=age,\n", - " plot_engine=PlotEngine.PYGMT,\n", - ")\n", - "\"\"\"\n", - "gplot = get_gplot(MODEL_NAME,MODEL_REPO_DIR,55)\n", + "gplot = get_gplot(\"merdith2021\", \"plate-model-repo\", age=55, plot_engine=PlotEngine.PYGMT)\n", "fig = pygmt.Figure()\n", "fig.basemap(region=\"d\", projection=\"N180/10c\", frame=\"lrtb\")\n", "#fig.coast(shorelines=True)\n", @@ -95,7 +62,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.11.9" + "version": "3.11.6" } }, "nbformat": 4, From 8a6134484d01455098474b8de0d5d81fb018e313 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Sat, 19 Oct 2024 15:11:22 +1100 Subject: [PATCH 07/14] remove old placeholder files --- gplately/plotting/pygmt_plotter.py | 0 gplately/plotting/reconstructor.py | 0 2 files changed, 0 insertions(+), 0 deletions(-) delete mode 100644 gplately/plotting/pygmt_plotter.py delete mode 100644 gplately/plotting/reconstructor.py diff --git a/gplately/plotting/pygmt_plotter.py b/gplately/plotting/pygmt_plotter.py deleted file mode 100644 index e69de29b..00000000 diff --git a/gplately/plotting/reconstructor.py b/gplately/plotting/reconstructor.py deleted file mode 100644 index e69de29b..00000000 From 5aea1542e66daaa5bc6de25a544244c8f878f982 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Sat, 19 Oct 2024 15:12:57 +1100 Subject: [PATCH 08/14] rollback change in setup.py --- setup.py | 1 - 1 file changed, 1 deletion(-) diff --git a/setup.py b/setup.py index 05e81f84..b3e18d42 100644 --- a/setup.py +++ b/setup.py @@ -96,7 +96,6 @@ def _minimal_ext_cmd(cmd): "netcdf4", "rasterio", "geopandas", - "gmt", "stripy", "plate-model-manager>=1.2.1", "pyyaml", From 8753b71bb1d6dffb622b010a9ca70afc96dec2bc Mon Sep 17 00:00:00 2001 From: michaelchin Date: Thu, 16 Jan 2025 13:08:23 +1100 Subject: [PATCH 09/14] add parameters for pygmt plot --- gplately/mapping/pygmt_plot.py | 53 ++++++++++++++++++++++++++++- gplately/plot.py | 5 +++ tests-dir/unittest/test_pygmt.ipynb | 32 +++-------------- 3 files changed, 62 insertions(+), 28 deletions(-) diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py index 821d7b9b..c6cb7628 100644 --- a/gplately/mapping/pygmt_plot.py +++ b/gplately/mapping/pygmt_plot.py @@ -17,6 +17,16 @@ from geopandas.geodataframe import GeoDataFrame import pygmt +pygmt.config( + FONT_ANNOT=8, + FONT_LABEL=8, + FONT=8, + MAP_TICK_PEN="0.75p", + MAP_FRAME_PEN="0.75p", + MAP_TICK_LENGTH_PRIMARY="4p", +) + + # ----- parameters for plot region = "d" width = 10 @@ -33,8 +43,49 @@ label_position = "TC" +def get_pygmt_basemap_figure(projection="N180/10c", region="d"): + fig = pygmt.Figure() + fig.basemap(region=region, projection=projection, frame="lrtb") + return fig + + def plot_geo_data_frame(fig: pygmt.Figure, gdf: GeoDataFrame, **kwargs): - fig.plot(data=gdf.geometry, pen="0.5p,blue") + line_width = "0.1p" + line_color = "blue" + + if "edgecolor" in kwargs.keys(): + if isinstance(kwargs["edgecolor"], str): + line_color = kwargs["edgecolor"] + else: + raise Exception( + "The edgecolor parameter is not string. Currently, the pygmt plot engine only supports colour name." + ) + + if "linewidth" in kwargs.keys(): + line_width = f"{kwargs['linewidth']}p" + + fill = None + if "facecolor" in kwargs.keys() and kwargs["facecolor"].lower() != "none": + fill = f"{kwargs['facecolor']}" + + if line_color.lower() == "none": + line_width = "0" + line_color = fill + + if "fill" in kwargs.keys(): + fill = kwargs["fill"] + + if "pen" in kwargs.keys(): + pen = kwargs["pen"] + else: + pen = f"{line_width},{line_color}" + + style = None + if "style" in kwargs.keys(): + style = kwargs["style"] + + fig.plot(data=gdf.geometry, pen=pen, fill=fill, style=style) + """ fig.plot(data=gdf_coastlines, fill=coastline_color, frame=["xa0", "ya0"], transparency=0) diff --git a/gplately/plot.py b/gplately/plot.py index 37992832..5d079e2c 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -726,6 +726,11 @@ def _plot_feature(self, ax, get_feature_func, **kwargs): tessellate_degrees=tessellate_degrees, ) + if not isinstance(gdf, gpd.GeoDataFrame): + raise Exception( + f"Expecting a GeoDataFrame object, but the gdf is {type(gdf)}" + ) + if len(gdf) == 0: logger.warning("No feature found for plotting. Do nothing and return.") return ax diff --git a/tests-dir/unittest/test_pygmt.ipynb b/tests-dir/unittest/test_pygmt.ipynb index 6f9be13a..474f39cd 100644 --- a/tests-dir/unittest/test_pygmt.ipynb +++ b/tests-dir/unittest/test_pygmt.ipynb @@ -7,43 +7,21 @@ "metadata": {}, "outputs": [], "source": [ - "import sys\n", - "\n", - "from plate_model_manager import PlateModel, PlateModelManager\n", - "\n", - "from gplately import PlateReconstruction, PlotTopologies\n", "from gplately.mapping.plot_engine import PlotEngine\n", + "from gplately.mapping.pygmt_plot import get_pygmt_basemap_figure\n", "from gplately.auxiliary import get_gplot\n", - "import pygmt\n", - "pygmt.config(\n", - " FONT_ANNOT=8,\n", - " FONT_LABEL=8,\n", - " FONT=8,\n", - " MAP_TICK_PEN=\"0.75p\",\n", - " MAP_FRAME_PEN=\"0.75p\",\n", - " MAP_TICK_LENGTH_PRIMARY=\"4p\",\n", - " )\n", "\n", "gplot = get_gplot(\"merdith2021\", \"plate-model-repo\", age=55, plot_engine=PlotEngine.PYGMT)\n", - "fig = pygmt.Figure()\n", - "fig.basemap(region=\"d\", projection=\"N180/10c\", frame=\"lrtb\")\n", + "fig = get_pygmt_basemap_figure(projection=\"N180/10c\", region=\"d\")\n", "#fig.coast(shorelines=True)\n", "\n", - "gplot.plot_topological_plate_boundaries(fig)\n", - "gplot.plot_coastlines(fig)\n", + "gplot.plot_topological_plate_boundaries(fig, edgecolor=\"blue\", linewidth=1.0, central_meridian=180, style='f0.2/0.08+l+t')\n", + "gplot.plot_coastlines(fig, edgecolor=\"none\", facecolor='gray', linewidth=0.1, central_meridian=180)\n", "\n", "fig.show(width=1200)\n", "# fig.savefig(\"test-pygmt-plot.pdf\")\n", "\n" ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "b25b75a3-39f1-4494-912e-7f9f502a67ce", - "metadata": {}, - "outputs": [], - "source": [] } ], "metadata": { @@ -62,7 +40,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.11.6" + "version": "3.11.9" } }, "nbformat": 4, From cc8694cc505a919150913fcd6587b8c02124622d Mon Sep 17 00:00:00 2001 From: michaelchin Date: Fri, 17 Jan 2025 13:40:40 +1100 Subject: [PATCH 10/14] plot subduction with pygmt --- gplately/auxiliary.py | 4 ++-- gplately/mapping/plot_engine.py | 2 +- gplately/mapping/pygmt_plot.py | 29 +++++++++++++++++++++++---- gplately/plot.py | 22 ++++++++++++++------ tests-dir/unittest/test_pygmt.ipynb | 11 +++++++--- tests-dir/unittest/test_pygmt_plot.py | 4 ++-- 6 files changed, 54 insertions(+), 18 deletions(-) diff --git a/gplately/auxiliary.py b/gplately/auxiliary.py index 27c89983..43cccdc8 100644 --- a/gplately/auxiliary.py +++ b/gplately/auxiliary.py @@ -2,7 +2,7 @@ from plate_model_manager import PlateModel, PlateModelManager -from .mapping.plot_engine import PlotEngine +from .mapping.plot_engine import PlotEngineType from .plot import PlotTopologies from .reconstruction import PlateReconstruction @@ -11,7 +11,7 @@ def get_gplot( model_name: str, model_repo_dir: str, age: Union[int, float], - plot_engine: PlotEngine = PlotEngine.CARTOPY, + plot_engine: PlotEngineType = PlotEngineType.CARTOPY, ) -> PlotTopologies: """auxiliary function to get gplot object""" try: diff --git a/gplately/mapping/plot_engine.py b/gplately/mapping/plot_engine.py index 77af7ecc..d0949397 100644 --- a/gplately/mapping/plot_engine.py +++ b/gplately/mapping/plot_engine.py @@ -18,6 +18,6 @@ from enum import Enum -class PlotEngine(Enum): +class PlotEngineType(Enum): CARTOPY = 1 PYGMT = 2 diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py index c6cb7628..e48d1ee1 100644 --- a/gplately/mapping/pygmt_plot.py +++ b/gplately/mapping/pygmt_plot.py @@ -26,7 +26,6 @@ MAP_TICK_LENGTH_PRIMARY="4p", ) - # ----- parameters for plot region = "d" width = 10 @@ -35,9 +34,7 @@ # plate boundary stuff plateboundary_width = "0.5p" - age_font = "12p,Helvetica,black" - label_font = "12p,Helvetica,black" label_offset = "j0/-0.5c" label_position = "TC" @@ -49,6 +46,24 @@ def get_pygmt_basemap_figure(projection="N180/10c", region="d"): return fig +def plot_subduction_zones( + fig: pygmt.Figure, + gdf_subduction_left: GeoDataFrame, + gdf_subduction_right: GeoDataFrame, + color="blue", + **kwargs, +): + fig.plot( + data=gdf_subduction_left, pen=f"0.5p,{color}", fill=color, style="f0.2/0.08+l+t" + ) + fig.plot( + data=gdf_subduction_right, + pen=f"0.5p,{color}", + fill=color, + style="f0.2/0.08+r+t", + ) + + def plot_geo_data_frame(fig: pygmt.Figure, gdf: GeoDataFrame, **kwargs): line_width = "0.1p" line_color = "blue" @@ -84,7 +99,13 @@ def plot_geo_data_frame(fig: pygmt.Figure, gdf: GeoDataFrame, **kwargs): if "style" in kwargs.keys(): style = kwargs["style"] - fig.plot(data=gdf.geometry, pen=pen, fill=fill, style=style) + label = None + if "gmtlabel" in kwargs.keys(): + label = kwargs["gmtlabel"] + + fig.plot( + data=gdf.geometry, pen=pen, fill=fill, style=style, transparency=0, label=label + ) """ fig.plot(data=gdf_coastlines, fill=coastline_color, frame=["xa0", "ya0"], transparency=0) diff --git a/gplately/plot.py b/gplately/plot.py index 5d079e2c..6fb3e84c 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -45,8 +45,11 @@ validate_topology_availability, ) from .gpml import _load_FeatureCollection -from .mapping.plot_engine import PlotEngine -from .mapping.pygmt_plot import plot_geo_data_frame +from .mapping.plot_engine import PlotEngineType +from .mapping.pygmt_plot import ( + plot_geo_data_frame, + plot_subduction_zones as pygmt_plot_subduction_zones, +) from .pygplates import FeatureCollection as _FeatureCollection from .reconstruction import PlateReconstruction as _PlateReconstruction from .tools import EARTH_RADIUS @@ -300,7 +303,7 @@ def __init__( COBs=None, time=None, anchor_plate_id=0, - plot_engine: PlotEngine = PlotEngine.CARTOPY, + plot_engine: PlotEngineType = PlotEngineType.CARTOPY, ): self.plate_reconstruction = plate_reconstruction @@ -735,7 +738,7 @@ def _plot_feature(self, ax, get_feature_func, **kwargs): logger.warning("No feature found for plotting. Do nothing and return.") return ax - if self._plot_engine == PlotEngine.PYGMT: + if self._plot_engine == PlotEngineType.PYGMT: return plot_geo_data_frame(fig=ax, gdf=gdf, **kwargs) else: if hasattr(ax, "projection"): @@ -1083,10 +1086,10 @@ def get_subduction_direction(self): trench_right_features = shapelify_feature_lines(self.trench_right) gdf_left = gpd.GeoDataFrame( - {"geometry": trench_left_features}, geometry="geometry" + {"geometry": trench_left_features}, geometry="geometry", crs="EPSG:4326" ) gdf_right = gpd.GeoDataFrame( - {"geometry": trench_right_features}, geometry="geometry" + {"geometry": trench_right_features}, geometry="geometry", crs="EPSG:4326" ) return gdf_left, gdf_right @@ -1128,6 +1131,13 @@ def plot_subduction_teeth( See `Matplotlib` keyword arguments [here](https://matplotlib.org/stable/api/_as_gen/matplotlib.pyplot.plot.html). """ + if self._plot_engine == PlotEngineType.PYGMT: + gdf_subduction_left, gdf_subduction_right = self.get_subduction_direction() + pygmt_plot_subduction_zones( + ax, gdf_subduction_left, gdf_subduction_right, color=color, **kwargs + ) + return + if not self.plate_reconstruction.topology_features: logger.warn( "Plate model does not have topology features. Unable to plot_subduction_teeth." diff --git a/tests-dir/unittest/test_pygmt.ipynb b/tests-dir/unittest/test_pygmt.ipynb index 474f39cd..85fbb333 100644 --- a/tests-dir/unittest/test_pygmt.ipynb +++ b/tests-dir/unittest/test_pygmt.ipynb @@ -7,16 +7,21 @@ "metadata": {}, "outputs": [], "source": [ - "from gplately.mapping.plot_engine import PlotEngine\n", + "from gplately.mapping.plot_engine import PlotEngineType\n", "from gplately.mapping.pygmt_plot import get_pygmt_basemap_figure\n", "from gplately.auxiliary import get_gplot\n", "\n", - "gplot = get_gplot(\"merdith2021\", \"plate-model-repo\", age=55, plot_engine=PlotEngine.PYGMT)\n", + "gplot = get_gplot(\"merdith2021\", \"plate-model-repo\", age=55, plot_engine=PlotEngineType.PYGMT)\n", "fig = get_pygmt_basemap_figure(projection=\"N180/10c\", region=\"d\")\n", "#fig.coast(shorelines=True)\n", "\n", - "gplot.plot_topological_plate_boundaries(fig, edgecolor=\"blue\", linewidth=1.0, central_meridian=180, style='f0.2/0.08+l+t')\n", + "gplot.plot_topological_plate_boundaries(fig, edgecolor=\"black\", linewidth=0.25, central_meridian=180)\n", "gplot.plot_coastlines(fig, edgecolor=\"none\", facecolor='gray', linewidth=0.1, central_meridian=180)\n", + "gplot.plot_ridges(fig, pen='0.5p,red', gmtlabel=\"ridges\")\n", + "gplot.plot_subduction_teeth(fig, color=\"blue\")\n", + "\n", + "fig.text(text='what is this', position=\"TC\", no_clip=True, font=\"12p,Helvetica,black\", offset=\"j0/-0.5c\")\n", + "fig.legend(position='jBL+o-2.7/0', box=\"+gwhite+p0.5p\")\n", "\n", "fig.show(width=1200)\n", "# fig.savefig(\"test-pygmt-plot.pdf\")\n", diff --git a/tests-dir/unittest/test_pygmt_plot.py b/tests-dir/unittest/test_pygmt_plot.py index 19448ab6..4d1e9dd0 100644 --- a/tests-dir/unittest/test_pygmt_plot.py +++ b/tests-dir/unittest/test_pygmt_plot.py @@ -10,7 +10,7 @@ import gplately from gplately import PlateReconstruction, PlotTopologies -from gplately.mapping.plot_engine import PlotEngine +from gplately.mapping.plot_engine import PlotEngineType print(gplately.__file__) @@ -43,7 +43,7 @@ def main(show=True): COBs=model.get_layer("COBs", return_none_if_not_exist=True), continents=model.get_layer("ContinentalPolygons"), time=age, - plot_engine=PlotEngine.PYGMT, + plot_engine=PlotEngineType.PYGMT, ) # age = 100 From 0acf995001a64364d871b0a0c9ebc424d5b38398 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Sat, 18 Jan 2025 08:07:59 +1100 Subject: [PATCH 11/14] add PlotEngine class --- gplately/mapping/__init__.py | 1 + gplately/mapping/cartopy_plot.py | 0 gplately/mapping/plot_engine.py | 23 +++++++++++++++++++++++ gplately/mapping/pygmt_plot.py | 24 ++++++++++++++++++++++++ 4 files changed, 48 insertions(+) create mode 100644 gplately/mapping/__init__.py create mode 100644 gplately/mapping/cartopy_plot.py diff --git a/gplately/mapping/__init__.py b/gplately/mapping/__init__.py new file mode 100644 index 00000000..2bf22597 --- /dev/null +++ b/gplately/mapping/__init__.py @@ -0,0 +1 @@ +# This submodule contains code to plot maps. diff --git a/gplately/mapping/cartopy_plot.py b/gplately/mapping/cartopy_plot.py new file mode 100644 index 00000000..e69de29b diff --git a/gplately/mapping/plot_engine.py b/gplately/mapping/plot_engine.py index d0949397..ac0d529a 100644 --- a/gplately/mapping/plot_engine.py +++ b/gplately/mapping/plot_engine.py @@ -16,8 +16,31 @@ # from enum import Enum +from abc import ABC, abstractmethod + +from geopandas.geodataframe import GeoDataFrame class PlotEngineType(Enum): CARTOPY = 1 PYGMT = 2 + + +class PlotEngine(ABC): + @abstractmethod + def plot_geo_data_frame(self, gdf: GeoDataFrame, **kwargs): + pass # This is an abstract method, no implementation here. + + @abstractmethod + def plot_pygplates_features(self, features, **kwargs): + pass # This is an abstract method, no implementation here. + + @abstractmethod + def plot_subduction_zones( + self, + gdf_subduction_left: GeoDataFrame, + gdf_subduction_right: GeoDataFrame, + color="blue", + **kwargs, + ): + pass # This is an abstract method, no implementation here. diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py index e48d1ee1..efaacf46 100644 --- a/gplately/mapping/pygmt_plot.py +++ b/gplately/mapping/pygmt_plot.py @@ -16,6 +16,7 @@ # from geopandas.geodataframe import GeoDataFrame import pygmt +from plot_engine import PlotEngine pygmt.config( FONT_ANNOT=8, @@ -40,6 +41,29 @@ label_position = "TC" +class PygmtPlotEngine(PlotEngine): + def __init__(self, projection="N180/10c", region="d"): + self.fig = pygmt.Figure() + self.fig.basemap(region=region, projection=projection, frame="lrtb") + + def plot_geo_data_frame(self, gdf: GeoDataFrame, **kwargs): + plot_geo_data_frame(self.fig, gdf, **kwargs) + + def plot_pygplates_features(self, features, **kwargs): + pass + + def plot_subduction_zones( + self, + gdf_subduction_left: GeoDataFrame, + gdf_subduction_right: GeoDataFrame, + color="blue", + **kwargs, + ): + plot_subduction_zones( + self.fig, gdf_subduction_left, gdf_subduction_right, color=color, **kwargs + ) + + def get_pygmt_basemap_figure(projection="N180/10c", region="d"): fig = pygmt.Figure() fig.basemap(region=region, projection=projection, frame="lrtb") From e93d3531b853a83645a42b2e2ff538f63467888c Mon Sep 17 00:00:00 2001 From: michaelchin Date: Mon, 20 Jan 2025 08:29:06 +1100 Subject: [PATCH 12/14] add PygmtPlotEngine and CartopyPlotEngine --- gplately/auxiliary.py | 5 +- gplately/decorators.py | 2 +- gplately/mapping/cartopy_plot.py | 87 ++++++++++++++++ gplately/mapping/plot_engine.py | 8 +- gplately/mapping/pygmt_plot.py | 16 +-- gplately/plot.py | 150 ++++++++-------------------- tests-dir/unittest/test_plot.py | 29 +++--- tests-dir/unittest/test_pygmt.ipynb | 15 ++- 8 files changed, 177 insertions(+), 135 deletions(-) diff --git a/gplately/auxiliary.py b/gplately/auxiliary.py index 43cccdc8..cf8f918c 100644 --- a/gplately/auxiliary.py +++ b/gplately/auxiliary.py @@ -2,7 +2,8 @@ from plate_model_manager import PlateModel, PlateModelManager -from .mapping.plot_engine import PlotEngineType +from .mapping.plot_engine import PlotEngine +from .mapping.cartopy_plot import CartopyPlotEngine from .plot import PlotTopologies from .reconstruction import PlateReconstruction @@ -11,7 +12,7 @@ def get_gplot( model_name: str, model_repo_dir: str, age: Union[int, float], - plot_engine: PlotEngineType = PlotEngineType.CARTOPY, + plot_engine: PlotEngine = CartopyPlotEngine(), ) -> PlotTopologies: """auxiliary function to get gplot object""" try: diff --git a/gplately/decorators.py b/gplately/decorators.py index 932eafcb..f21e9cc5 100644 --- a/gplately/decorators.py +++ b/gplately/decorators.py @@ -60,7 +60,7 @@ def inner(func_pointer): @wraps(func_pointer) def wrapper(self, ax, **kwargs): if not self.plate_reconstruction.topology_features: - logger.warn( + logger.warning( f"Plate model does not have topology features. Unable to plot {feature_name}." ) return ax diff --git a/gplately/mapping/cartopy_plot.py b/gplately/mapping/cartopy_plot.py index e69de29b..2f58f621 100644 --- a/gplately/mapping/cartopy_plot.py +++ b/gplately/mapping/cartopy_plot.py @@ -0,0 +1,87 @@ +import logging, math +from geopandas.geodataframe import GeoDataFrame +import cartopy.crs as ccrs +from .plot_engine import PlotEngine +from ..utils.plot_utils import _clean_polygons, plot_subduction_teeth +from ..tools import EARTH_RADIUS + +logger = logging.getLogger("gplately") + +DEFAULT_CARTOPY_PROJECTION = ccrs.PlateCarree() + + +class CartopyPlotEngine(PlotEngine): + def __init__(self): + pass + + def plot_geo_data_frame(self, ax_or_fig, gdf: GeoDataFrame, **kwargs): + if hasattr(ax_or_fig, "projection"): + gdf = _clean_polygons(data=gdf, projection=ax_or_fig.projection) + else: + kwargs["transform"] = DEFAULT_CARTOPY_PROJECTION + + return gdf.plot(ax=ax_or_fig, **kwargs) + + def plot_pygplates_features(self, ax_or_fig, features, **kwargs): + pass + + def plot_subduction_zones( + self, + ax_or_fig, + gdf_subduction_left: GeoDataFrame, + gdf_subduction_right: GeoDataFrame, + color="blue", + **kwargs, + ): + if "transform" in kwargs.keys(): + logger.warning( + "'transform' keyword argument is ignored by CartopyPlotEngine." + ) + kwargs.pop("transform") + + spacing = kwargs.pop("spacing") + size = kwargs.pop("size") + aspect = kwargs.pop("aspect") + + try: + projection = ax_or_fig.projection + except AttributeError: + logger.warning( + "The ax.projection does not exist. You must set projection to plot Cartopy maps, such as ax = plt.subplot(211, projection=cartopy.crs.PlateCarree())" + ) + projection = None + + if isinstance(projection, ccrs.PlateCarree): + spacing = math.degrees(spacing) + else: + spacing = spacing * EARTH_RADIUS * 1e3 + + if aspect is None: + aspect = 2.0 / 3.0 + if size is None: + size = spacing * 0.5 + + height = size * aspect + + plot_subduction_teeth( + gdf_subduction_left, + size, + "l", + height, + spacing, + projection=projection, + ax=ax_or_fig, + color=color, + **kwargs, + ) + plot_subduction_teeth( + gdf_subduction_right, + size, + "r", + height, + spacing, + projection=projection, + ax=ax_or_fig, + color=color, + **kwargs, + ) diff --git a/gplately/mapping/plot_engine.py b/gplately/mapping/plot_engine.py index ac0d529a..ed65782f 100644 --- a/gplately/mapping/plot_engine.py +++ b/gplately/mapping/plot_engine.py @@ -28,19 +28,23 @@ class PlotEngineType(Enum): class PlotEngine(ABC): @abstractmethod - def plot_geo_data_frame(self, gdf: GeoDataFrame, **kwargs): + def plot_geo_data_frame(self, ax_or_fig, gdf: GeoDataFrame, **kwargs): + """Plot GeoPandas GeoDataFrame object""" pass # This is an abstract method, no implementation here. @abstractmethod - def plot_pygplates_features(self, features, **kwargs): + def plot_pygplates_features(self, ax_or_fig, features, **kwargs): + """Plot one or more pygplates feature(s)""" pass # This is an abstract method, no implementation here. @abstractmethod def plot_subduction_zones( self, + ax_or_fig, gdf_subduction_left: GeoDataFrame, gdf_subduction_right: GeoDataFrame, color="blue", **kwargs, ): + """Plot subduction zones with teeth""" pass # This is an abstract method, no implementation here. diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py index efaacf46..4376c931 100644 --- a/gplately/mapping/pygmt_plot.py +++ b/gplately/mapping/pygmt_plot.py @@ -16,7 +16,7 @@ # from geopandas.geodataframe import GeoDataFrame import pygmt -from plot_engine import PlotEngine +from .plot_engine import PlotEngine pygmt.config( FONT_ANNOT=8, @@ -42,25 +42,25 @@ class PygmtPlotEngine(PlotEngine): - def __init__(self, projection="N180/10c", region="d"): - self.fig = pygmt.Figure() - self.fig.basemap(region=region, projection=projection, frame="lrtb") + def __init__(self): + pass - def plot_geo_data_frame(self, gdf: GeoDataFrame, **kwargs): - plot_geo_data_frame(self.fig, gdf, **kwargs) + def plot_geo_data_frame(self, ax_or_fig, gdf: GeoDataFrame, **kwargs): + plot_geo_data_frame(ax_or_fig, gdf, **kwargs) - def plot_pygplates_features(self, features, **kwargs): + def plot_pygplates_features(self, ax_or_fig, features, **kwargs): pass def plot_subduction_zones( self, + ax_or_fig, gdf_subduction_left: GeoDataFrame, gdf_subduction_right: GeoDataFrame, color="blue", **kwargs, ): plot_subduction_zones( - self.fig, gdf_subduction_left, gdf_subduction_right, color=color, **kwargs + ax_or_fig, gdf_subduction_left, gdf_subduction_right, color=color, **kwargs ) diff --git a/gplately/plot.py b/gplately/plot.py index 6fb3e84c..2d0475de 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -45,17 +45,22 @@ validate_topology_availability, ) from .gpml import _load_FeatureCollection -from .mapping.plot_engine import PlotEngineType +from .mapping.plot_engine import PlotEngine from .mapping.pygmt_plot import ( - plot_geo_data_frame, + PygmtPlotEngine, plot_subduction_zones as pygmt_plot_subduction_zones, ) +from .mapping.cartopy_plot import CartopyPlotEngine, DEFAULT_CARTOPY_PROJECTION from .pygplates import FeatureCollection as _FeatureCollection from .reconstruction import PlateReconstruction as _PlateReconstruction from .tools import EARTH_RADIUS from .utils.feature_utils import shapelify_features as _shapelify_features -from .utils.plot_utils import _clean_polygons, _meridian_from_ax -from .utils.plot_utils import plot_subduction_teeth as _plot_subduction_teeth +from .utils.plot_utils import ( + _clean_polygons, + _meridian_from_ax, + plot_subduction_teeth as _plot_subduction_teeth, +) + logger = logging.getLogger("gplately") @@ -303,15 +308,16 @@ def __init__( COBs=None, time=None, anchor_plate_id=0, - plot_engine: PlotEngineType = PlotEngineType.CARTOPY, + plot_engine: PlotEngine = CartopyPlotEngine(), ): + self._plot_engine = plot_engine self.plate_reconstruction = plate_reconstruction if self.plate_reconstruction.topology_features is None: self.plate_reconstruction.topology_features = [] logger.warning("Plate model does not have topology features.") - self.base_projection = ccrs.PlateCarree() + self.base_projection = DEFAULT_CARTOPY_PROJECTION # store these for when time is updated # make sure these are initialised as FeatureCollection objects @@ -383,7 +389,7 @@ def __setstate__(self, state): self._COBs.add(_FeatureCollection(feature)) self._anchor_plate_id = state["plate_id"] - self.base_projection = ccrs.PlateCarree() + self.base_projection = DEFAULT_CARTOPY_PROJECTION self._time = None @property @@ -693,7 +699,7 @@ def get_feature( tessellate_degrees=tessellate_degrees, ) - return gpd.GeoDataFrame({"geometry": shp}, geometry="geometry") + return gpd.GeoDataFrame({"geometry": shp}, geometry="geometry") # type: ignore @append_docstring(PLOT_DOCSTRING.format("feature")) def plot_feature(self, ax, feature, feature_name="", color="black", **kwargs): @@ -710,7 +716,7 @@ def plot_feature(self, ax, feature, feature_name="", color="black", **kwargs): kwargs["facecolor"] = "none" return self._plot_feature(ax, partial(self.get_feature, feature), **kwargs) - def _plot_feature(self, ax, get_feature_func, **kwargs): + def _plot_feature(self, ax, get_feature_func, **kwargs) -> None: if "transform" in kwargs.keys(): warnings.warn( "'transform' keyword argument is ignored by PlotTopologies", @@ -738,15 +744,7 @@ def _plot_feature(self, ax, get_feature_func, **kwargs): logger.warning("No feature found for plotting. Do nothing and return.") return ax - if self._plot_engine == PlotEngineType.PYGMT: - return plot_geo_data_frame(fig=ax, gdf=gdf, **kwargs) - else: - if hasattr(ax, "projection"): - gdf = _clean_polygons(data=gdf, projection=ax.projection) - else: - kwargs["transform"] = self.base_projection - - return gdf.plot(ax=ax, **kwargs) + self._plot_engine.plot_geo_data_frame(ax, gdf, **kwargs) @validate_reconstruction_time @append_docstring(GET_DATE_DOCSTRING.format("coastlines")) @@ -1040,7 +1038,8 @@ def plot_subduction_teeth_deprecated( return ax.add_geometries(teeth, crs=self.base_projection, color=color, **kwargs) - def get_subduction_direction(self): + @validate_reconstruction_time + def get_subduction_direction(self, central_meridian=0.0, tessellate_degrees=None): """Create a geopandas.GeoDataFrame object containing geometries of trench directions. Notes @@ -1072,31 +1071,36 @@ def get_subduction_direction(self): If the optional `time` parameter has not been passed to `PlotTopologies`. This is needed to construct `trench_left` or `trench_right` geometries to the requested `time` and thus populate the GeoDataFrame. """ - if self._time is None: - raise ValueError( - "No miscellaneous topologies have been resolved. Set `PlotTopologies.time` to construct them." - ) - if self.trench_left is None or self.trench_right is None: - raise ValueError( - "No trench_left or trench_right topologies passed to PlotTopologies." + raise Exception( + "No subduction zone/trench data is found. Make sure the plate model has topology feature." ) - trench_left_features = shapelify_feature_lines(self.trench_left) - trench_right_features = shapelify_feature_lines(self.trench_right) + trench_left_features = shapelify_feature_lines( + self.trench_left, + tessellate_degrees=tessellate_degrees, + central_meridian=central_meridian, + ) + trench_right_features = shapelify_feature_lines( + self.trench_right, + tessellate_degrees=tessellate_degrees, + central_meridian=central_meridian, + ) gdf_left = gpd.GeoDataFrame( {"geometry": trench_left_features}, geometry="geometry", crs="EPSG:4326" - ) + ) # type: ignore gdf_right = gpd.GeoDataFrame( {"geometry": trench_right_features}, geometry="geometry", crs="EPSG:4326" - ) + ) # type: ignore return gdf_left, gdf_right + @validate_reconstruction_time + @validate_topology_availability("Subduction Zones") def plot_subduction_teeth( self, ax, spacing=0.07, size=None, aspect=None, color="black", **kwargs - ): + ) -> None: """Plot subduction teeth onto a standard map Projection. Notes @@ -1131,84 +1135,18 @@ def plot_subduction_teeth( See `Matplotlib` keyword arguments [here](https://matplotlib.org/stable/api/_as_gen/matplotlib.pyplot.plot.html). """ - if self._plot_engine == PlotEngineType.PYGMT: - gdf_subduction_left, gdf_subduction_right = self.get_subduction_direction() - pygmt_plot_subduction_zones( - ax, gdf_subduction_left, gdf_subduction_right, color=color, **kwargs - ) - return - - if not self.plate_reconstruction.topology_features: - logger.warn( - "Plate model does not have topology features. Unable to plot_subduction_teeth." - ) - return ax - if self._time is None: - raise ValueError( - "No topologies have been resolved. Set `PlotTopologies.time` to construct them." - ) - if "transform" in kwargs.keys(): - warnings.warn( - "'transform' keyword argument is ignored by PlotTopologies", - UserWarning, - ) - kwargs.pop("transform") + kwargs["spacing"] = spacing + kwargs["size"] = size + kwargs["aspect"] = aspect central_meridian = _meridian_from_ax(ax) tessellate_degrees = np.rad2deg(spacing) - - try: - projection = ax.projection - except AttributeError: - print( - "The ax.projection does not exist. You must set projection to plot Cartopy maps, such as ax = plt.subplot(211, projection=cartopy.crs.PlateCarree())" - ) - projection = None - - if isinstance(projection, ccrs.PlateCarree): - spacing = math.degrees(spacing) - else: - spacing = spacing * EARTH_RADIUS * 1e3 - - if aspect is None: - aspect = 2.0 / 3.0 - if size is None: - size = spacing * 0.5 - - height = size * aspect - - trench_left_features = shapelify_feature_lines( - self.trench_left, - tessellate_degrees=tessellate_degrees, - central_meridian=central_meridian, - ) - trench_right_features = shapelify_feature_lines( - self.trench_right, - tessellate_degrees=tessellate_degrees, - central_meridian=central_meridian, + gdf_subduction_left, gdf_subduction_right = self.get_subduction_direction( + tessellate_degrees=tessellate_degrees, central_meridian=central_meridian ) - plot_subduction_teeth( - trench_left_features, - size, - "l", - height, - spacing, - projection=projection, - ax=ax, - color=color, - **kwargs, - ) - plot_subduction_teeth( - trench_right_features, - size, - "r", - height, - spacing, - projection=projection, - ax=ax, - color=color, - **kwargs, + self._plot_engine.plot_subduction_zones( + ax, gdf_subduction_left, gdf_subduction_right, color=color, **kwargs ) def plot_plate_polygon_by_id(self, ax, plate_id, color="black", **kwargs): @@ -1471,7 +1409,7 @@ def compute_radius(ortho, radius_degrees): # adding a patch patch = ax.add_patch( - mpatches.Circle(xy=[lon, lat], radius=r_ortho, transform=proj1, **kwargs) + mpatches.Circle(xy=(lon, lat), radius=r_ortho, transform=proj1, **kwargs) ) return patch @@ -1858,7 +1796,7 @@ def get_all_topologies( "feature_name": feature_names, }, geometry="geometry", - ) + ) # type: ignore return gdf @validate_topology_availability("all topologies") diff --git a/tests-dir/unittest/test_plot.py b/tests-dir/unittest/test_plot.py index 3530a284..f9d85b63 100755 --- a/tests-dir/unittest/test_plot.py +++ b/tests-dir/unittest/test_plot.py @@ -30,12 +30,16 @@ def main(show=True): except: model = PlateModel(MODEL_NAME, data_dir=MODEL_REPO_DIR, readonly=True) + if not model: + return + age = 55 test_model = PlateReconstruction( model.get_rotation_model(), topology_features=model.get_layer("Topologies"), static_polygons=model.get_layer("StaticPolygons"), + plate_model_name=MODEL_NAME, ) gplot = PlotTopologies( test_model, @@ -51,7 +55,7 @@ def main(show=True): fig = plt.figure(figsize=(10, 5), dpi=96) ax = fig.add_subplot(111, projection=ccrs.Robinson(central_longitude=180)) - all_flag = 0 + all_flag = 1 plot_flag = { "continent_ocean_boundaries": 0, "coastlines": 0, @@ -96,17 +100,18 @@ def main(show=True): ax, color=list(np.random.choice(range(256), size=3) / 256) ) - ax.set_global() - - ids = set([f.get_reconstruction_plate_id() for f in gplot.topologies]) - for id in ids: - if all_flag or plot_flag["plate_polygon_by_id"]: - gplot.plot_plate_polygon_by_id( - ax, - id, - facecolor="None", - edgecolor=list(np.random.choice(range(256), size=3) / 256), - ) + ax.set_global() # type: ignore + + if gplot.topologies: + ids = set([f.get_reconstruction_plate_id() for f in gplot.topologies]) + for id in ids: + if all_flag or plot_flag["plate_polygon_by_id"]: + gplot.plot_plate_polygon_by_id( + ax, + id, + facecolor="None", + edgecolor=list(np.random.choice(range(256), size=3) / 256), + ) plt.title(f"{age} Ma") if show: diff --git a/tests-dir/unittest/test_pygmt.ipynb b/tests-dir/unittest/test_pygmt.ipynb index 85fbb333..5e059a0a 100644 --- a/tests-dir/unittest/test_pygmt.ipynb +++ b/tests-dir/unittest/test_pygmt.ipynb @@ -7,11 +7,10 @@ "metadata": {}, "outputs": [], "source": [ - "from gplately.mapping.plot_engine import PlotEngineType\n", - "from gplately.mapping.pygmt_plot import get_pygmt_basemap_figure\n", + "from gplately.mapping.pygmt_plot import PygmtPlotEngine, get_pygmt_basemap_figure\n", "from gplately.auxiliary import get_gplot\n", "\n", - "gplot = get_gplot(\"merdith2021\", \"plate-model-repo\", age=55, plot_engine=PlotEngineType.PYGMT)\n", + "gplot = get_gplot(\"merdith2021\", \"plate-model-repo\", age=55, plot_engine=PygmtPlotEngine())\n", "fig = get_pygmt_basemap_figure(projection=\"N180/10c\", region=\"d\")\n", "#fig.coast(shorelines=True)\n", "\n", @@ -20,13 +19,21 @@ "gplot.plot_ridges(fig, pen='0.5p,red', gmtlabel=\"ridges\")\n", "gplot.plot_subduction_teeth(fig, color=\"blue\")\n", "\n", - "fig.text(text='what is this', position=\"TC\", no_clip=True, font=\"12p,Helvetica,black\", offset=\"j0/-0.5c\")\n", + "fig.text(text='55Ma (Merdith2021)', position=\"TC\", no_clip=True, font=\"12p,Helvetica,black\", offset=\"j0/-0.5c\")\n", "fig.legend(position='jBL+o-2.7/0', box=\"+gwhite+p0.5p\")\n", "\n", "fig.show(width=1200)\n", "# fig.savefig(\"test-pygmt-plot.pdf\")\n", "\n" ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "64cb9ddb-359d-4c67-90ae-05ba8f9e35c5", + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": { From d6d1de77d004d181c76246d71982aec3f9df53fb Mon Sep 17 00:00:00 2001 From: michaelchin Date: Mon, 20 Jan 2025 09:39:00 +1100 Subject: [PATCH 13/14] clean up code, remove temporary comments --- gplately/mapping/pygmt_plot.py | 93 +++++++---------------------- gplately/plot.py | 2 +- tests-dir/unittest/test_pygmt.ipynb | 5 +- 3 files changed, 24 insertions(+), 76 deletions(-) diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py index 4376c931..e7fdb545 100644 --- a/gplately/mapping/pygmt_plot.py +++ b/gplately/mapping/pygmt_plot.py @@ -27,18 +27,7 @@ MAP_TICK_LENGTH_PRIMARY="4p", ) -# ----- parameters for plot -region = "d" -width = 10 -projection = "N180/" -x_offset = width + 2 - -# plate boundary stuff -plateboundary_width = "0.5p" -age_font = "12p,Helvetica,black" -label_font = "12p,Helvetica,black" -label_offset = "j0/-0.5c" -label_position = "TC" +# NW's example is at https://gist.github.com/nickywright/f53018a8eda29223cca6f39ab2cfa25d class PygmtPlotEngine(PlotEngine): @@ -77,8 +66,14 @@ def plot_subduction_zones( color="blue", **kwargs, ): + label = kwargs.pop("gmtlabel", None) + fig.plot( - data=gdf_subduction_left, pen=f"0.5p,{color}", fill=color, style="f0.2/0.08+l+t" + data=gdf_subduction_left, + pen=f"0.5p,{color}", + fill=color, + style="f0.2/0.08+l+t", + label=label, ) fig.plot( data=gdf_subduction_right, @@ -89,72 +84,24 @@ def plot_subduction_zones( def plot_geo_data_frame(fig: pygmt.Figure, gdf: GeoDataFrame, **kwargs): - line_width = "0.1p" - line_color = "blue" - if "edgecolor" in kwargs.keys(): - if isinstance(kwargs["edgecolor"], str): - line_color = kwargs["edgecolor"] - else: - raise Exception( - "The edgecolor parameter is not string. Currently, the pygmt plot engine only supports colour name." - ) + line_color = kwargs.pop("edgecolor", "blue") + line_width = f"{kwargs.pop('linewidth',0.1)}p" - if "linewidth" in kwargs.keys(): - line_width = f"{kwargs['linewidth']}p" - - fill = None - if "facecolor" in kwargs.keys() and kwargs["facecolor"].lower() != "none": - fill = f"{kwargs['facecolor']}" + fill = kwargs.pop("facecolor", None) + if fill and fill.lower() == "none": + fill = None + fill = kwargs.pop("fill", fill) # the "fill" parameter override the "facecolor" if line_color.lower() == "none": - line_width = "0" - line_color = fill - - if "fill" in kwargs.keys(): - fill = kwargs["fill"] - - if "pen" in kwargs.keys(): - pen = kwargs["pen"] + # line_width = "0" + # line_color = fill + pen = None else: - pen = f"{line_width},{line_color}" - - style = None - if "style" in kwargs.keys(): - style = kwargs["style"] - - label = None - if "gmtlabel" in kwargs.keys(): - label = kwargs["gmtlabel"] + pen = kwargs.pop("pen", f"{line_width},{line_color}") + style = kwargs.pop("style", None) + label = kwargs.pop("gmtlabel", None) fig.plot( data=gdf.geometry, pen=pen, fill=fill, style=style, transparency=0, label=label ) - - """ - fig.plot(data=gdf_coastlines, fill=coastline_color, frame=["xa0", "ya0"], transparency=0) - - fig.plot(data=gdf_topo_plates.geometry, pen='%s,%s' % (plateboundary_width, plate_colour), frame="lrtb") - fig.plot(data=gdf_subduction_left, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+l+t') - fig.plot(data=gdf_subduction_right, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+r+t') - fig.plot(data=gdf_ridges_transforms, pen='%s,%s' % (plateboundary_width, ridge_colour)) - fig.plot(data=gplot.get_transforms(), pen='%s,%s' % (plateboundary_width, transform_color)) - - fig.text(text='gplot.get_transforms(): %s Ma' % age, position=label_position, no_clip=True, font=label_font, offset=label_offset) - - fig.shift_origin(xshift=x_offset) - fig.basemap(region=region, projection="%s%sc" % (projection, width), frame="lrtb") - fig.plot(data=gdf_cobs, fill=COB_color, transparency=0, ) - fig.plot(data=gdf_coastlines, fill=coastline_color, frame=["xa0", "ya0"], transparency=0) - - fig.plot(data=gdf_topo_plates.geometry, pen='%s,%s' % (plateboundary_width, plate_colour), frame="lrtb", label='other plate boundary types') - fig.plot(data=gdf_subduction_left, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+l+t', label='subduction zones') - fig.plot(data=gdf_subduction_right, pen='%s,%s' % (plateboundary_width, subduction_zone_colour), fill=subduction_zone_colour, style='f0.2/0.08+r+t') - fig.plot(data=gdf_ridges_transforms, pen='%s,%s' % (plateboundary_width, ridge_colour), label='ridges and transforms') - - # from gpml: transforms - fig.plot(data=gdf_topo_transforms, pen='%s,%s' % (plateboundary_width, transform_color), label = 'transforms') - fig.text(text='FeatureType.gpml_transform: %s Ma' % age, position=label_position, no_clip=True, font=label_font, offset=label_offset) - - fig.legend(position='jBL+o-2.7/0', box="+gwhite+p0.5p") - """ diff --git a/gplately/plot.py b/gplately/plot.py index 2d0475de..240f5327 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -699,7 +699,7 @@ def get_feature( tessellate_degrees=tessellate_degrees, ) - return gpd.GeoDataFrame({"geometry": shp}, geometry="geometry") # type: ignore + return gpd.GeoDataFrame({"geometry": shp}, geometry="geometry", crs="EPSG:4326") # type: ignore @append_docstring(PLOT_DOCSTRING.format("feature")) def plot_feature(self, ax, feature, feature_name="", color="black", **kwargs): diff --git a/tests-dir/unittest/test_pygmt.ipynb b/tests-dir/unittest/test_pygmt.ipynb index 5e059a0a..085a0ca2 100644 --- a/tests-dir/unittest/test_pygmt.ipynb +++ b/tests-dir/unittest/test_pygmt.ipynb @@ -14,10 +14,11 @@ "fig = get_pygmt_basemap_figure(projection=\"N180/10c\", region=\"d\")\n", "#fig.coast(shorelines=True)\n", "\n", - "gplot.plot_topological_plate_boundaries(fig, edgecolor=\"black\", linewidth=0.25, central_meridian=180)\n", + "gplot.plot_topological_plate_boundaries(fig, edgecolor=\"black\", linewidth=0.25, central_meridian=180, gmtlabel=\"plate boundaries\")\n", "gplot.plot_coastlines(fig, edgecolor=\"none\", facecolor='gray', linewidth=0.1, central_meridian=180)\n", "gplot.plot_ridges(fig, pen='0.5p,red', gmtlabel=\"ridges\")\n", - "gplot.plot_subduction_teeth(fig, color=\"blue\")\n", + "gplot.plot_transforms(fig, pen='0.5p,red', gmtlabel=\"transforms\")\n", + "gplot.plot_subduction_teeth(fig, color=\"blue\", gmtlabel=\"subduction zones\")\n", "\n", "fig.text(text='55Ma (Merdith2021)', position=\"TC\", no_clip=True, font=\"12p,Helvetica,black\", offset=\"j0/-0.5c\")\n", "fig.legend(position='jBL+o-2.7/0', box=\"+gwhite+p0.5p\")\n", From 9092c9d1db46d69c3ee54a2f2290bad6062b4f41 Mon Sep 17 00:00:00 2001 From: michaelchin Date: Mon, 20 Jan 2025 12:53:44 +1100 Subject: [PATCH 14/14] add comments for PlotEngine classes --- gplately/grids.py | 2 +- gplately/mapping/cartopy_plot.py | 25 ++++++ gplately/mapping/plot_engine.py | 2 +- gplately/mapping/pygmt_plot.py | 130 +++++++++++++++++++------------ gplately/plot.py | 77 +++++++++--------- 5 files changed, 148 insertions(+), 88 deletions(-) diff --git a/gplately/grids.py b/gplately/grids.py index a2ea4595..4237fdec 100644 --- a/gplately/grids.py +++ b/gplately/grids.py @@ -167,7 +167,7 @@ def read_netcdf_grid( x_dimension_name: str = "", y_dimension_name: str = "", data_variable_name: str = "", -): +) -> tuple[np.ma.MaskedArray, np.ma.MaskedArray, np.ma.MaskedArray] | np.ma.MaskedArray: """Read a `netCDF` (.nc) grid from a given `filename` and return its data as a `MaskedArray`. diff --git a/gplately/mapping/cartopy_plot.py b/gplately/mapping/cartopy_plot.py index 2f58f621..f26c04bd 100644 --- a/gplately/mapping/cartopy_plot.py +++ b/gplately/mapping/cartopy_plot.py @@ -15,6 +15,16 @@ def __init__(self): pass def plot_geo_data_frame(self, ax_or_fig, gdf: GeoDataFrame, **kwargs): + """Plot GeoDataFrame object with Cartopy + + Parameters + ---------- + ax_or_fig : cartopy.mpl.geoaxes.GeoAxes + Cartopy GeoAxes instance + gdf : GeoDataFrame + GeoPandas GeoDataFrame object + + """ if hasattr(ax_or_fig, "projection"): gdf = _clean_polygons(data=gdf, projection=ax_or_fig.projection) else: @@ -23,6 +33,7 @@ def plot_geo_data_frame(self, ax_or_fig, gdf: GeoDataFrame, **kwargs): return gdf.plot(ax=ax_or_fig, **kwargs) def plot_pygplates_features(self, ax_or_fig, features, **kwargs): + """TODO""" pass def plot_subduction_zones( @@ -33,6 +44,20 @@ def plot_subduction_zones( color="blue", **kwargs, ): + """Plot subduction zones with "teeth" using pygmt + + Parameters + ---------- + ax_or_fig : cartopy.mpl.geoaxes.GeoAxes + Cartopy GeoAxes instance + gdf_subduction_left : GeoDataFrame + subduction zone with "left" polarity + gdf_subduction_right : GeoDataFrame + subduction zone with "right" polarity + color : str + The colour used to fill the "teeth". + + """ if "transform" in kwargs.keys(): logger.warning( "'transform' keyword argument is ignored by CartopyPlotEngine." diff --git a/gplately/mapping/plot_engine.py b/gplately/mapping/plot_engine.py index ed65782f..491c05f3 100644 --- a/gplately/mapping/plot_engine.py +++ b/gplately/mapping/plot_engine.py @@ -46,5 +46,5 @@ def plot_subduction_zones( color="blue", **kwargs, ): - """Plot subduction zones with teeth""" + """Plot subduction zones with "teeth" """ pass # This is an abstract method, no implementation here. diff --git a/gplately/mapping/pygmt_plot.py b/gplately/mapping/pygmt_plot.py index e7fdb545..20d15bfc 100644 --- a/gplately/mapping/pygmt_plot.py +++ b/gplately/mapping/pygmt_plot.py @@ -35,9 +35,57 @@ def __init__(self): pass def plot_geo_data_frame(self, ax_or_fig, gdf: GeoDataFrame, **kwargs): - plot_geo_data_frame(ax_or_fig, gdf, **kwargs) + """Plot GeoDataFrame object with pygmt + + Parameters + ---------- + ax_or_fig : pygmt.Figure() + pygmt Figure object + gdf : GeoDataFrame + GeoPandas GeoDataFrame object + edgecolor : str + For polygons, it is the border colour. For polylines, it is the line colour. + Currently, only colour names are tested and officially supported, for example, "red", "blue", etc. + facecolor : str + The colour used to fill the polygon. + fill : str + GMT "fill" parameter + pen : str + GMT "pen" parameter + style : str + GMT "style" parameter + gmtlabel : str + GMT "label" parameter + + """ + line_color = kwargs.pop("edgecolor", "blue") + line_width = f"{kwargs.pop('linewidth',0.1)}p" + + fill = kwargs.pop("facecolor", None) + if fill and fill.lower() == "none": + fill = None + fill = kwargs.pop("fill", fill) # the "fill" parameter override the "facecolor" + + if line_color.lower() == "none": + # line_width = "0" + # line_color = fill + pen = None + else: + pen = kwargs.pop("pen", f"{line_width},{line_color}") + style = kwargs.pop("style", None) + label = kwargs.pop("gmtlabel", None) + + ax_or_fig.plot( + data=gdf.geometry, + pen=pen, + fill=fill, + style=style, + transparency=0, + label=label, + ) def plot_pygplates_features(self, ax_or_fig, features, **kwargs): + """TODO""" pass def plot_subduction_zones( @@ -48,60 +96,40 @@ def plot_subduction_zones( color="blue", **kwargs, ): - plot_subduction_zones( - ax_or_fig, gdf_subduction_left, gdf_subduction_right, color=color, **kwargs + """Plot subduction zones with "teeth" using pygmt + + Parameters + ---------- + ax_or_fig : pygmt.Figure() + pygmt Figure object + gdf_subduction_left : GeoDataFrame + subduction zone with "left" polarity + gdf_subduction_right : GeoDataFrame + subduction zone with "right" polarity + color : str + The colour used to fill the "teeth". + gmtlabel : str + GMT "label" parameter + """ + label = kwargs.pop("gmtlabel", None) + + ax_or_fig.plot( + data=gdf_subduction_left, + pen=f"0.5p,{color}", + fill=color, + style="f0.2/0.08+l+t", + label=label, + ) + ax_or_fig.plot( + data=gdf_subduction_right, + pen=f"0.5p,{color}", + fill=color, + style="f0.2/0.08+r+t", ) def get_pygmt_basemap_figure(projection="N180/10c", region="d"): + """return a pygmt.Figure() object""" fig = pygmt.Figure() fig.basemap(region=region, projection=projection, frame="lrtb") return fig - - -def plot_subduction_zones( - fig: pygmt.Figure, - gdf_subduction_left: GeoDataFrame, - gdf_subduction_right: GeoDataFrame, - color="blue", - **kwargs, -): - label = kwargs.pop("gmtlabel", None) - - fig.plot( - data=gdf_subduction_left, - pen=f"0.5p,{color}", - fill=color, - style="f0.2/0.08+l+t", - label=label, - ) - fig.plot( - data=gdf_subduction_right, - pen=f"0.5p,{color}", - fill=color, - style="f0.2/0.08+r+t", - ) - - -def plot_geo_data_frame(fig: pygmt.Figure, gdf: GeoDataFrame, **kwargs): - - line_color = kwargs.pop("edgecolor", "blue") - line_width = f"{kwargs.pop('linewidth',0.1)}p" - - fill = kwargs.pop("facecolor", None) - if fill and fill.lower() == "none": - fill = None - fill = kwargs.pop("fill", fill) # the "fill" parameter override the "facecolor" - - if line_color.lower() == "none": - # line_width = "0" - # line_color = fill - pen = None - else: - pen = kwargs.pop("pen", f"{line_width},{line_color}") - style = kwargs.pop("style", None) - label = kwargs.pop("gmtlabel", None) - - fig.plot( - data=gdf.geometry, pen=pen, fill=fill, style=style, transparency=0, label=label - ) diff --git a/gplately/plot.py b/gplately/plot.py index 240f5327..ef114c87 100644 --- a/gplately/plot.py +++ b/gplately/plot.py @@ -46,17 +46,13 @@ ) from .gpml import _load_FeatureCollection from .mapping.plot_engine import PlotEngine -from .mapping.pygmt_plot import ( - PygmtPlotEngine, - plot_subduction_zones as pygmt_plot_subduction_zones, -) + from .mapping.cartopy_plot import CartopyPlotEngine, DEFAULT_CARTOPY_PROJECTION from .pygplates import FeatureCollection as _FeatureCollection from .reconstruction import PlateReconstruction as _PlateReconstruction from .tools import EARTH_RADIUS from .utils.feature_utils import shapelify_features as _shapelify_features from .utils.plot_utils import ( - _clean_polygons, _meridian_from_ax, plot_subduction_teeth as _plot_subduction_teeth, ) @@ -481,7 +477,7 @@ def update_time(self, time): anchor_plate_id=self.anchor_plate_id, # use ResolveTopologyType.boundary parameter to resolve rigid plate boundary only # because the Mid-ocean ridges(and transforms) should not contain lines from topological networks - resolve_topology_types=pygplates.ResolveTopologyType.boundary, + resolve_topology_types=pygplates.ResolveTopologyType.boundary, # type: ignore ) self._topologies = None @@ -502,64 +498,64 @@ def update_time(self, time): self.unclassified_features = [] for topol in self.other: - if topol.get_feature_type() == pygplates.FeatureType.gpml_continental_rift: + if topol.get_feature_type() == pygplates.FeatureType.gpml_continental_rift: # type: ignore self.continental_rifts.append(topol) - elif topol.get_feature_type() == pygplates.FeatureType.gpml_fault: + elif topol.get_feature_type() == pygplates.FeatureType.gpml_fault: # type: ignore self.faults.append(topol) - elif topol.get_feature_type() == pygplates.FeatureType.gpml_fracture_zone: + elif topol.get_feature_type() == pygplates.FeatureType.gpml_fracture_zone: # type: ignore self.fracture_zones.append(topol) elif ( topol.get_feature_type() - == pygplates.FeatureType.gpml_inferred_paleo_boundary + == pygplates.FeatureType.gpml_inferred_paleo_boundary # type: ignore ): self.inferred_paleo_boundaries.append(topol) elif ( - topol.get_feature_type() == pygplates.FeatureType.gpml_terrane_boundary + topol.get_feature_type() == pygplates.FeatureType.gpml_terrane_boundary # type: ignore ): self.terrane_boundaries.append(topol) elif ( topol.get_feature_type() - == pygplates.FeatureType.gpml_transitional_crust + == pygplates.FeatureType.gpml_transitional_crust # type: ignore ): self.transitional_crusts.append(topol) - elif topol.get_feature_type() == pygplates.FeatureType.gpml_orogenic_belt: + elif topol.get_feature_type() == pygplates.FeatureType.gpml_orogenic_belt: # type: ignore self.orogenic_belts.append(topol) - elif topol.get_feature_type() == pygplates.FeatureType.gpml_suture: + elif topol.get_feature_type() == pygplates.FeatureType.gpml_suture: # type: ignore self.sutures.append(topol) elif ( - topol.get_feature_type() == pygplates.FeatureType.gpml_continental_crust + topol.get_feature_type() == pygplates.FeatureType.gpml_continental_crust # type: ignore ): self.continental_crusts.append(topol) elif ( topol.get_feature_type() - == pygplates.FeatureType.gpml_extended_continental_crust + == pygplates.FeatureType.gpml_extended_continental_crust # type: ignore ): self.extended_continental_crusts.append(topol) elif ( topol.get_feature_type() - == pygplates.FeatureType.gpml_passive_continental_boundary + == pygplates.FeatureType.gpml_passive_continental_boundary # type: ignore ): self.passive_continental_boundaries.append(topol) - elif topol.get_feature_type() == pygplates.FeatureType.gpml_slab_edge: + elif topol.get_feature_type() == pygplates.FeatureType.gpml_slab_edge: # type: ignore self.slab_edges.append(topol) - elif topol.get_feature_type() == pygplates.FeatureType.gpml_transform: + elif topol.get_feature_type() == pygplates.FeatureType.gpml_transform: # type: ignore self.transforms.append(topol) elif ( topol.get_feature_type() - == pygplates.FeatureType.gpml_unclassified_feature + == pygplates.FeatureType.gpml_unclassified_feature # type: ignore ): self.unclassified_features.append(topol) @@ -618,7 +614,7 @@ def _tessellate_triangles( triangle_pointsX = [] triangle_pointsY = [] - date_line_wrapper = pygplates.DateLineWrapper() + date_line_wrapper = pygplates.DateLineWrapper() # type: ignore for feature in features: cum_distance = 0.0 @@ -1170,27 +1166,36 @@ def plot_plate_polygon_by_id(self, ax, plate_id, color="black", **kwargs): [here](https://matplotlib.org/3.5.1/api/_as_gen/matplotlib.axes.Axes.imshow.html). """ - return self.plot_feature( + features = [] + if self.topologies: + features = ( + [ + feature + for feature in self.topologies + if feature.get_reconstruction_plate_id() == plate_id + ], + ) + self.plot_feature( ax, - [ - feature - for feature in self.topologies - if feature.get_reconstruction_plate_id() == plate_id - ], + features, color=color, **kwargs, ) - # the old function name(plot_plate_id) is bad. we should change the name - # for backward compatibility, we have to allow users to use the old name def plot_plate_id(self, *args, **kwargs): - logger.warn( - "The class method plot_plate_id will be deprecated in the future soon. Use plot_plate_polygon_by_id instead." + """TODO: remove this function + + The function name plot_plate_id() is bad and should be changed. + The new name is plot_plate_polygon_by_id(). + For backward compatibility, we allow users to use the old name in their legcy code for now. + No new code should call this function. + + """ + logger.warning( + "The class method plot_plate_id is deprecated and will be removed in the future soon. Use plot_plate_polygon_by_id instead." ) return self.plot_plate_polygon_by_id(*args, **kwargs) - plot_plate_id.__doc__ = plot_plate_polygon_by_id.__doc__ - def plot_grid(self, ax, grid, extent=[-180, 180, -90, 90], **kwargs): """Plot a `MaskedArray` raster or grid onto a standard map Projection. @@ -1796,6 +1801,7 @@ def get_all_topologies( "feature_name": feature_names, }, geometry="geometry", + crs="EPSG:4326", ) # type: ignore return gdf @@ -1868,7 +1874,8 @@ def get_all_topological_sections( "feature_name": feature_names, }, geometry="geometry", - ) + crs="EPSG:4326", + ) # type: ignore return gdf @validate_topology_availability("all topological sections") @@ -1888,7 +1895,7 @@ def _resolve_both_boundaries_and_networks(self): if self._topologies is None: resolved_topologies = [] shared_boundary_sections = [] - pygplates.resolve_topologies( + pygplates.resolve_topologies( # type: ignore self.plate_reconstruction.topology_features, self.plate_reconstruction.rotation_model, resolved_topologies,