update examples

Author: michaelchin

examples/output/plot_subduction_zones.png

@@ -40,8 +40,8 @@ def plot_teeth(ax, lons, lats, polarity):
 
 
 def main(show=True):
-    time = 0
-    model_name = "Merdith2021"
+    time = 100
+    model_name = "Alfonso2024"
 
     fig, ax = get_basemap_with_coastlines(model=model_name, time=time)
 
@@ -71,7 +71,7 @@ def main(show=True):
                 ax.plot(lons, lats, transform=ccrs.PlateCarree(), color="blue")
                 plot_teeth(ax, lons, lats, polarity)
 
-    plt.title(f"{time} Ma (Merdith2021)")
+    plt.title(f"{time} Ma ({model_name})")
 
     fig.text(
         0.5,
@@ -86,7 +86,7 @@ def main(show=True):
         save_fig(__file__)
 
 
-def get_basemap_with_coastlines(model="Muller2019", crs=ccrs.Robinson(), time=140):
+def get_basemap_with_coastlines(model, crs=ccrs.Robinson(), time=140):
     coastlines_shapely = coastlines.get_paleo_coastlines(
         time=time, model=model, format="shapely"
     )
@@ -110,9 +110,9 @@ def get_basemap_with_coastlines(model="Muller2019", crs=ccrs.Robinson(), time=14
     ax.add_geometries(
         coastlines_shapely,
         crs=ccrs.PlateCarree(),
-        facecolor="grey",
+        facecolor="lightgrey",
         edgecolor="none",
-        alpha=0.5,
+        alpha=1,
     )
     return fig, ax
 

examples/output/plot_topological_plate_boundaries.png

@@ -7,7 +7,7 @@
 import numpy as np
 from shapely.geometry import MultiLineString, MultiPolygon, Polygon, shape
 
-from gwspy import PlateModel
+from gwspy import PlateModel, coastlines
 
 OUTPUT_DIR = "output"
 Path(OUTPUT_DIR).mkdir(parents=True, exist_ok=True)
@@ -16,41 +16,55 @@
 # export GWS_URL=http://localhost:18000/
 # micromamba run -n gplates-ws-example ./plot_topological_plate_boundaries.py
 
+model_name = "Alfonso2024"
+
 
 def main(show=True):
-    model = PlateModel("Muller2019")
-    time = 10
-    topology_10 = model.get_topology(time)
+    model = PlateModel(model_name)
+    time = 100
+    topology_100 = model.get_topology(time)
 
-    _lines = topology_10.get_plate_boundaries()
+    _lines = topology_100.get_plate_boundaries()
     colors = []
     feature_types = []
     polylines = []
-    color_map = {}
+    color_map = {
+        "MidOceanRidge": "red",
+        "Transform": "red",
+        "SubductionZone": "blue",
+        "Fault": "orange",
+    }
     for feature in _lines["features"]:
         s = shape(feature["geometry"])
         polylines.append(s)
         f_type = feature["properties"]["type"]
-        feature_types.append(f_type)
+
         if f_type not in color_map:
-            cc = list(np.random.choice(range(256), size=3) / 256)
-            color_map[f_type] = cc
+            cc = "grey"
+            feature_types.append("other")
+            # cc = list(np.random.choice(range(256), size=3) / 256)
+            # color_map[f_type] = cc
         else:
+            feature_types.append(f_type)
             cc = color_map[f_type]
         colors.append(cc)
         # if f_type == "gpml:TopologicalNetwork":
         #    colors.append("green")
 
     print(set(feature_types))
 
+    coastlines_shapely = coastlines.get_paleo_coastlines(
+        time=time, model=model_name, format="shapely"
+    )
+
     fig = plt.figure(figsize=(12, 6), dpi=120)
     ax = plt.axes(projection=ccrs.Robinson())
 
     ax.set_global()
     gl = ax.gridlines(
         crs=ccrs.PlateCarree(),
         draw_labels=True,
-        color="grey",
+        color="black",
         alpha=0.5,
         linestyle="--",
     )
@@ -59,17 +73,26 @@ def main(show=True):
     gl.xlabel_style = {"size": 7, "color": "gray"}
     gl.ylabel_style = {"size": 7, "color": "gray"}
 
+    ax.add_geometries(
+        coastlines_shapely,
+        crs=ccrs.PlateCarree(),
+        facecolor="lightgrey",
+        edgecolor="none",
+        alpha=0.5,
+    )
+
     # plot the lines
     for line, c, t in zip(polylines, colors, feature_types):
         if isinstance(line, MultiLineString):
             for g in line.geoms:
                 ax.plot(
-                    *g.xy, transform=ccrs.PlateCarree(), color=c, label=t, linewidth=0.7
+                    *g.xy, transform=ccrs.PlateCarree(), color=c, label=t, linewidth=1
                 )
         else:
             ax.plot(
-                *line.xy, transform=ccrs.PlateCarree(), color=c, label=t, linewidth=0.7
+                *line.xy, transform=ccrs.PlateCarree(), color=c, label=t, linewidth=1
             )
+
     # plot the legend
     handles, labels = ax.get_legend_handles_labels()
     unique = [
@@ -80,16 +103,16 @@ def main(show=True):
         title="Feature Types",
         prop={"size": 6},
         loc="lower right",
-        bbox_to_anchor=(1.14, 0.6),
+        bbox_to_anchor=(1.08, 0.8),
     )
     plt.setp(legend.get_title(), fontsize="xx-small")
 
-    plt.title(f"{time} Ma (Muller2019)")
+    plt.title(f"{time} Ma ({model_name})")
 
     fig.text(
         0.5,
         0.03,
-        "topological plate boundaries with random colours for each feature type",
+        "the topological plate boundaries are plotted as line segments",
         ha="center",
     )
 

examples/plot_subduction_zones.py

@@ -15,25 +15,31 @@
 # export GWS_URL=http://localhost:18000/
 # micromamba run -n gplates-ws-example ./plot_topological_plate_polygons.py
 
+model_name = "Alfonso2024"
+
 
 def main(show=True):
-    model = PlateModel("Muller2019")
-    time = 10
-    topology_10 = model.get_topology(time)
+    model = PlateModel(model_name)
+    time = 100
+    topology_100 = model.get_topology(time)
 
     # LOOK HERE!!!!
     # plot polygon without edges first
     # and then get the polygon edges as lines and  plot the edges
     # this is to get rid of the vertical line along the dateline
-    _polygons = topology_10.get_plate_polygons()
-    _lines = topology_10.get_plate_polygons(as_lines=True)
+    _polygons = topology_100.get_plate_polygons()
+    _lines = topology_100.get_plate_polygons(as_lines=True)
 
     # print(_polygons)
 
     polygons = []
     for feature in _polygons["features"]:
-        # if feature["properties"]["type"] != "gpml:TopologicalClosedPlateBoundary":
+        # if "Africa Basins temporary rigid block" in feature["properties"]["name"]:
+        #    print(feature["properties"]["name"])
         #    continue
+        # print(feature["properties"]["name"])
+        if feature["properties"]["type"] == "gpml:OceanicCrust":
+            continue
         s = shape(feature["geometry"])
         if isinstance(s, Polygon) or isinstance(s, MultiPolygon):
             polygons.append(s.buffer(0))
@@ -43,16 +49,18 @@ def main(show=True):
     feature_types = []
     for feature in _lines["features"]:
         s = shape(feature["geometry"])
-        lines.append(s)
         f_type = feature["properties"]["type"]
+        if f_type == "gpml:OceanicCrust":
+            continue
+
+        lines.append(s)
+
         if f_type == "gpml:TopologicalNetwork":
             colors.append("green")
         elif f_type == "gpml:TopologicalClosedPlateBoundary":
             colors.append("blue")
         elif f_type == "gpml:TopologicalSlabBoundary":
             colors.append("yellow")
-        elif f_type == "gpml:OceanicCrust":
-            colors.append("orange")
         else:
             colors.append("red")
 
@@ -110,7 +118,7 @@ def main(show=True):
         bbox_to_anchor=(1.1, 0.9),
     )
     plt.setp(legend.get_title(), fontsize="xx-small")
-    plt.title(f"{time} Ma (Muller2019)")
+    plt.title(f"{time} Ma ({model_name})")
 
     fig.text(
         0.5,