Hi Danny,
So it sounds like you want the points of the reconstructed plate (at a particular time). In this case you can use pyGPlates to resolve (reconstruct) a topological model to a reconstruction time, and then find the Farallon plate using its feature ID, and then you extract its reconstructed geometry.
So first you would need the feature ID which can be found in GPlates (right-click it and Copy
)…
Then use something like the following code with pyGPlates:
def find_resolved_plate(topological_snapshot, plate_feature_id):
"""Find a resolved topological boundary in a topological snapshot using the specified feature ID."""
resolved_plates = topological_snapshot.get_resolved_topologies(
pygplates.ResolveTopologyType.boundary)
for resolved_plate in resolved_plates:
if (resolved_plate.get_feature().get_feature_id().get_string() ==
plate_feature_id):
return resolved_plate
return None # not found
def get_resolved_plate_id_and_name_and_geometry(resolved_plate):
"""Return the plate ID, name and list of boundary (lon, lat) points of a resolved plate."""
# Get plate ID and name.
plate_feature = resolved_plate.get_feature()
plate_name = plate_feature.get_name()
plate_id = plate_feature.get_reconstruction_plate_id()
# Get reconstructed plate polygon's points as a list of (lon, lat).
resolved_plate_polygon = resolved_plate.get_resolved_boundary()
plate_boundary_list_of_lon_lat = []
for point in resolved_plate_polygon.get_exterior_ring_points():
lat, lon = point.to_lat_lon()
plate_boundary_list_of_lon_lat.append((lon, lat))
return plate_id, plate_name, plate_boundary_list_of_lon_lat
reconstruction_time = 205
farallon_plate_feature_id = 'GPlates-ec8fe0bf-2fc8-4dd6-8115-04f8d2532e35'
# Create a topological model from topological and rotation features.
topological_model = pygplates.TopologicalModel(topology_features, rotation_model)
# Create a topological snapshot at a specific time.
topological_snapshot = topological_model.topological_snapshot(reconstruction_time)
# Find a particular plate using its feature ID.
farallon_resolved_plate = find_resolved_plate(topological_snapshot, farallon_plate_feature_id)
if farallon_resolved_plate:
# Get plate ID and name and the reconstructed plate polygon's points as a list of (lon, lat).
farallon_plate_id, farallon_plate_name, farallon_plate_boundary_lon_lats = get_resolved_plate_id_and_name_and_geometry(farallon_resolved_plate)
…I haven’t tried the above code so there might be a typo in there. But that’s the general idea.
If you wanted to just get all large plates you could filter using the plate area rather than searching by feature IDs. So something like:
def find_large_resolved_plates(topological_snapshot, min_plate_area_steradians):
"""Find resolved topological boundaries in a topological snapshot exceeding a minimum area."""
large_resolved_plates = []
resolved_plates = topological_snapshot.get_resolved_topologies(
pygplates.ResolveTopologyType.boundary)
for resolved_plate in resolved_plates:
if (resolved_plate.get_resolved_boundary().get_area() >
min_plate_area_steradians):
large_resolved_plates.append(resolved_plate)
return large_resolved_plates
min_plate_area_square_kms = 20e+6
min_plate_area_steradians = min_plate_area_square_kms / pygplates.Earth.mean_radius_in_kms**2
# Create a topological model from topological and rotation features.
topological_model = pygplates.TopologicalModel(topology_features, rotation_model)
# Reconstruct topologies through a sequence of times.
for time in times:
# Create a topological snapshot at a specific time.
topological_snapshot = topological_model.topological_snapshot(time)
# Filter those rigid plates exceeding a minimum area.
large_resolved_topological_boundaries = find_large_resolved_plates(topological_snapshot, min_plate_area_steradians)
for large_resolved_topological_boundary in large_resolved_topological_boundaries:
# Get plate ID and name and the reconstructed plate polygon's points as a list of (lon, lat).
plate_id, plate_name, plate_boundary_lon_lats = get_resolved_plate_id_and_name_and_geometry(large_resolved_topological_boundary)
That’s mostly for getting present day geometry. But there are easier functions to use, documented at the top of the Feature class, that mean you don’t have to iterate through properties. For example, to get a regular present-day geometry use get_geometry()
or to get a topological geometry use get_topological_geometry()
. The latter would retrieve the intersecting lines that form a topological plate polygon (note that you don’t get an actual polygon until the topologies are reconstructed/resolved to a particular reconstruction time which is what the above code does).