Query common feature types
This example shows how to query various types of features.
See also
Query a coastline feature
In this example we query a coastline feature.
Sample code
import pygplates
# Load the global coastline features.
coastline_features = pygplates.FeatureCollection('coastlines.gpml')
# Iterate over the coastline features.
for feature in coastline_features:
# Print the feature type (Coastline) and the name of the coastline.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the coastline.
print ' description: %s' % feature.get_description()
# Print the plate ID of the coastline.
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
# Print the length of the coastline geometry(s).
# There could be more than one geometry per feature.
for geometry in feature.get_geometries():
print ' length: %f Kms' % (geometry.get_arc_length() * pygplates.Earth.mean_radius_in_kms)
# Print the valid time period of the coastline.
print ' valid time period: %f -> %f' % feature.get_valid_time()
Details
Print the coastline
pygplates.FeatureType,
name,
description,
plate ID,
length of the geometries and
valid time period.Some coastlines have more than one geometry per feature so we need to call
pygplates.Feature.get_geometries()
instead of pygplates.Feature.get_geometry() - the latter will fail if there’s not exactly one geometry.
So we print out the length of each coastline geometry in a coastline feature.print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
print ' description: %s' % feature.get_description()
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
for geometry in feature.get_geometries():
print ' length: %f Kms' % (geometry.get_arc_length() * pygplates.Earth.mean_radius_in_kms)
print ' valid time period: %f -> %f' % feature.get_valid_time()
Note
Each geometry length is converted from
radians to Kms by multiplying with the pygplates.Earth.mean_radius_in_kms attribute in
pygplates.Earth.
Output
Coastline: Pacific
description:
plate ID: 982
length: 143.477788 Kms
valid time period: 86.000000 -> -inf
Coastline: Nazca
description:
plate ID: 911
length: 55.408784 Kms
valid time period: 5.000000 -> -inf
Coastline: Nazca
description:
plate ID: 911
length: 421.757226 Kms
valid time period: 10.900000 -> -inf
...
Coastline: Mexico
description:
plate ID: 104
length: 138.553857 Kms
length: 186.286987 Kms
length: 13.861494 Kms
valid time period: 0.000000 -> -inf
Coastline: Mexico
description:
plate ID: 104
length: 135.613787 Kms
valid time period: 0.000000 -> -inf
...
Query an isochron feature
In this example we query an isochron feature.
This is the same as the Query a coastline feature example except we also query
the gpml:conjugatePlateId supported by
isochron features.
Sample code
import pygplates
# Load the global isochron features.
isochron_features = pygplates.FeatureCollection('isochrons.gpml')
# Iterate over the isochron features.
for feature in isochron_features:
# Print the feature type (Isochron) and the name of the isochron.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the isochron.
print ' description: %s' % feature.get_description()
# Print the plate ID of the isochron.
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
# Print the conjugate plate ID of the isochron.
print ' conjugate plate ID: %d' % feature.get_conjugate_plate_id()
# Print the length of the isochron geometry(s).
# There could be more than one geometry per feature.
for geometry in feature.get_geometries():
print ' length: %f Kms' % (geometry.get_arc_length() * pygplates.Earth.mean_radius_in_kms)
# Print the valid time period of the isochron.
print ' valid time period: %f -> %f' % feature.get_valid_time()
Details
Print the isochron
pygplates.FeatureType,
name,
description,
plate ID,
conjugate plate ID,
length of the geometries and
valid time period.Some isochrons have more than one geometry per feature so we need to call
pygplates.Feature.get_geometries()
instead of pygplates.Feature.get_geometry() - the latter will fail if there’s not exactly one geometry.
So we print out the length of each isochron geometry in a isochron feature.print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
print ' description: %s' % feature.get_description()
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
print ' conjugate plate ID: %d' % feature.get_conjugate_plate_id()
for geometry in feature.get_geometries():
print ' length: %f Kms' % (geometry.get_arc_length() * pygplates.Earth.mean_radius_in_kms)
print ' valid time period: %f -> %f' % feature.get_valid_time()
Note
Each geometry length is converted from
radians to Kms by multiplying with the pygplates.Earth.mean_radius_in_kms attribute in
pygplates.Earth.
Output
Some features mixed in with the isochrons were
passive continental boundaries
that were missing a conjugate plate ID.
In this case calls to
pygplates.Feature.get_conjugate_plate_id() returns a default value of zero
as can be seen in the following output:Isochron: CARLSBERG RIDGE, INDIA-AFRICA ANOMALY 5 ISOCHRON
description:
plate ID: 501
conjugate plate ID: 701
length: 2981.517878 Kms
valid time period: 10.900000 -> -inf
Isochron: CARLSBERG RIDGE, INDIA-AFRICA ANOMALY 6 ISOCHRON
description:
plate ID: 501
conjugate plate ID: 701
length: 2438.969434 Kms
valid time period: 20.100000 -> -inf
...
PassiveContinentalBoundary: NATL-ROCKALL ISO
description:
plate ID: 102
conjugate plate ID: 0
length: 1389.956801 Kms
valid time period: 55.000000 -> -inf
Isochron: LAB SEA COB
description:
plate ID: 102
conjugate plate ID: 101
length: 676.087194 Kms
valid time period: 65.000000 -> -inf
...
Query a mid-ocean ridge feature
In this example we query a mid-ocean ridge feature.
Sample code
import pygplates
# Load the global mid-ocean ridge features.
ridge_features = pygplates.FeatureCollection('ridges.gpml')
# Iterate over the mid-ocean ridge features.
for feature in ridge_features:
# Print the feature type (MidOceanRidge) and the name of the mid-ocean ridge.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the mid-ocean ridge.
print ' description: %s' % feature.get_description()
# A mid-ocean ridge can either reconstruct by plate ID or by half stage rotation.
# The former uses the reconstruction plate ID.
# The latter uses the left/right plate IDs.
if feature.get_reconstruction_method() == 'ByPlateId':
# Print the plate ID of the mid-ocean ridge.
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
# Print the conjugate plate ID of the mid-ocean ridge (if it has one - use None to test this).
conjugate_plate_id = feature.get_conjugate_plate_id(None)
if conjugate_plate_id is not None:
print ' conjugate plate ID: %d' % conjugate_plate_id
else:
# Print the left plate ID of the mid-ocean ridge.
print ' left plate ID: %d' % feature.get_left_plate()
# Print the right plate ID of the mid-ocean ridge.
print ' right plate ID: %d' % feature.get_right_plate()
# Print the length of the mid-ocean ridge geometry(s).
# There could be more than one geometry per feature.
for geometry in feature.get_geometries():
print ' length: %f Kms' % (geometry.get_arc_length() * pygplates.Earth.mean_radius_in_kms)
# Print the valid time period of the mid-ocean ridge.
print ' valid time period: %f -> %f' % feature.get_valid_time()
Output
MidOceanRidge: IS GRN_EUR, RI Fram Strait <identity>GPlates-0ebaff40-f6d5-475b-8521-
description:
plate ID: 102
length: 301.220813 Kms
valid time period: 0.000000 -> -inf
MidOceanRidge: IS GRN_EUR, RI GRN Sea <identity>GPlates-f0c10a34-758f-48f3-b577-9062
description:
plate ID: 102
length: 119.012880 Kms
valid time period: 0.000000 -> -inf
MidOceanRidge: ISO CANADA BAS XR <identity>GPlates-c218d19f-6acb-4ffb-ae8c-3dc240d4ec
description:
plate ID: 101
length: 478.158406 Kms
valid time period: 118.000000 -> -inf
Query a subduction zone feature
In this example we query a subduction zone feature.
Sample code
import pygplates
# Load the subduction zone features.
subduction_zone_features = pygplates.FeatureCollection('subduction_zones.gpml')
# Iterate over the subduction zone features.
for feature in subduction_zone_features:
# Print the feature type (SubductionZone) and the name of the feature.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the feature.
print ' description: %s' % feature.get_description()
# Print the plate ID of the feature.
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
# Print the conjugate plate ID of the feature (if it has one - use None to test this).
conjugate_plate_id = feature.get_conjugate_plate_id(None)
if conjugate_plate_id is not None:
print ' conjugate plate ID: %d' % conjugate_plate_id
# Print the subduction polarity of the feature.
# Default to 'Unknown' if there is no polarity property.
polarity = feature.get_enumeration(
pygplates.PropertyName.gpml_subduction_polarity,
'Unknown')
print ' polarity: %s' % polarity
# Print the length of the feature geometry(s).
# There could be more than one geometry per feature.
for geometry in feature.get_geometries():
print ' length: %f Kms' % (geometry.get_arc_length() * pygplates.Earth.mean_radius_in_kms)
# Print the valid time period of the feature.
print ' valid time period: %f -> %f' % feature.get_valid_time()
Output
SubductionZone: MACQUARIE
description:
plate ID: 901
polarity: Right
length: 2190.013773 Kms
valid time period: 15.000000 -> -inf
SubductionZone: Junction East seg for closure
description:
plate ID: 801
polarity: Right
length: 3210.625757 Kms
valid time period: 76.000000 -> 73.100000
SubductionZone: Junction East seg for closure
description:
plate ID: 801
polarity: Right
length: 2941.014242 Kms
valid time period: 78.000000 -> 76.100000
SubductionZone: Junction East seg for closure
description:
plate ID: 801
polarity: Right
length: 3056.550172 Kms
valid time period: 77.100000 -> 76.100000
SubductionZone: Alaska margin subduction
description:
plate ID: 625
polarity: Left
length: 7278.344241 Kms
valid time period: 250.000000 -> 200.100000
SubductionZone: Alaska Subduction from COB GS
description:
plate ID: 182
polarity: Right
length: 4147.145987 Kms
valid time period: 250.000000 -> 128.100000
Query a virtual geomagnetic pole feature
In this example we query a virtual geomagnetic pole feature.
Sample code
import pygplates
# Load the virtual geomagnetic pole features.
vgp_features = pygplates.FeatureCollection('vvirtual_geomagnetic_poles.gpml')
# Iterate over the virtual geomagnetic pole features.
for feature in vgp_features:
# Print the feature type (VirtualGeomagneticPole) and the name of the feature.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the feature.
print ' description: %s' % feature.get_description()
# Print the plate ID of the feature.
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
# Print the average inclination of the feature (if there is one - use None to test this).
average_inclination = feature.get_double(pygplates.PropertyName.gpml_average_inclination, None)
if average_inclination is not None:
print ' average inclination: %f' % average_inclination
# Print the average declination of the feature (if there is one - use None to test this).
average_declination = feature.get_double(pygplates.PropertyName.gpml_average_declination, None)
if average_declinationn is not None:
print ' average declinationn: %f' % average_declination
# Print the pole position uncertainty (if there is one - use None to test this)
pole_position_uncertainty = feature.get_double(pygplates.PropertyName.gpml_pole_a95, None)
if pole_position_uncertainty is not None:
print ' pole position uncertainty: %f' % pole_position_uncertainty
# Print the average age (if there is one - use None to test this)
average_age = feature.get_double(pygplates.PropertyName.gpml_average_age, None)
if average_age is not None:
print ' average age: %f' % average_age
# Print the VGP pole position.
# The default geometry is the pole position so we don't have to specify a property name.
pole_lat, pole_lon = feature.get_geometry().to_lat_lon()
print ' pole lat: %f, pole lon: %f' % (pole_lat, pole_lon)
# Print the average sample site position.
# We need to specify a property name otherwise we'll get the VGP pole position.
average_sample_site_lat, average_sample_site_lon = feature.get_geometry(
pygplates.PropertyName.gpml_average_sample_site_position).to_lat_lon()
print ' average sample site lat: %f, average sample site lon: %f' % (
average_sample_site_lat, average_sample_site_lon)
Output
VirtualGeomagneticPole: RM:-10 - 10Ma N= 3 (Dp col.) Lat Range: 50.4 to 48 (Dm col.)
description:
plate ID: 302
average inclination: 186.770000
average declinationn: -65.480000
pole position uncertainty: 8.550000
average age: 0.000000
pole lat: 85.130000, pole lon: 118.180000
average sample site lat: 49.540000, average sample site lon: 7.690000
VirtualGeomagneticPole: RM: 0 - 20Ma N= 4 (Dp col.) Lat Range: 50.4 to 45 (Dm col.)
description:
plate ID: 302
average inclination: 187.010000
average declinationn: -65.390000
pole position uncertainty: 5.770000
average age: 10.000000
pole lat: 85.220000, pole lon: 105.110000
average sample site lat: 48.410000, average sample site lon: 6.760000
VirtualGeomagneticPole: RM: 10 - 30Ma N= 2 (Dp col.) Lat Range: 50.8 to 45 (Dm col.)
description:
plate ID: 302
average inclination: 190.960000
average declinationn: -63.540000
pole position uncertainty: 23.100000
average age: 20.000000
pole lat: 81.970000, pole lon: 112.190000
average sample site lat: 47.910000, average sample site lon: 5.990000
Query a motion path feature
In this example we query a motion path feature.
See also
See also
Sample code
import pygplates
# Load the motion path features.
motion_path_features = pygplates.FeatureCollection('motion_paths.gpml')
# Iterate over the motion path features.
for feature in motion_path_features:
# Print the feature type (MotionPath) and the name of the motion path.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the motion path.
print ' description: %s' % feature.get_description()
# Print the plate ID of the motion path.
print ' plate ID: %d' % feature.get_reconstruction_plate_id()
# Print the relative plate ID of the motion path.
print ' relative plate ID: %d' % feature.get_relative_plate()
# Print the times of the motion path.
print ' times: ', feature.get_times()
# Print the seed points of the motion path.
for seed_point in feature.get_geometry().get_points():
print ' seed point lat: %f, seed point lon: %f' % seed_point.to_lat_lon()
# Print the valid time period of the motion path.
print ' valid time period: %f -> %f' % feature.get_valid_time()
Output
MotionPath:
description:
plate ID: 701
relative plate ID: 201
times: [[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0]
seed point lat: -19.000000, seed point lon: 12.500000
seed point lat: -28.000000, seed point lon: 15.700000
valid time period: 90.000000 -> 0.000000
Query a flowline feature
In this example we query a flowline feature.
See also
Sample code
import pygplates
# Load the flowline features.
flowline_features = pygplates.FeatureCollection('flowlines.gpml')
# Iterate over the flowline features.
for feature in flowline_features:
# Print the feature type (Flowline) and the name of the flowline.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the flowline.
print ' description: %s' % feature.get_description()
# Print the left plate ID of the flowline.
print ' left plate ID: %d' % feature.get_left_plate()
# Print the right plate ID of the flowline.
print ' right plate ID: %d' % feature.get_right_plate()
# Print the times of the flowline.
print ' times: ', feature.get_times()
# Print the seed points of the flowline.
for seed_point in feature.get_geometry().get_points():
print ' seed point lat: %f, seed point lon: %f' % seed_point.to_lat_lon()
# Print the valid time period of the flowline.
print ' valid time period: %f -> %f' % feature.get_valid_time()
Output
Flowline:
description:
left plate ID: 201
right plate ID: 701
times: [0.0, 10.0, 20.0, 30.0, 40.0]
seed point lat: -35.547600, seed point lon: -17.873000
seed point lat: -46.208000, seed point lon: -13.623000
valid time period: 40.000000 -> 0.000000
Query a total reconstruction sequence (rotation) feature
In this example we query a total reconstruction sequence feature.
See also
Sample code
import pygplates
# Load the rotation features.
rotation_features = pygplates.FeatureCollection('rotations.rot')
# Iterate over the rotation features.
for feature in rotation_features:
fixed_plate_id, moving_plate_id, total_reconstruction_pole = feature.get_total_reconstruction_pole()
# Ignore moving plate IDs equal to 999 since these are commented lines in the PLATES4 rotation format.
if moving_plate_id == 999:
continue
# Print the feature type (TotalReconstructionSequence) and the name of the rotation feature.
print '%s: %s' % (feature.get_feature_type().get_name(), feature.get_name())
# Print the description of the rotation feature.
print ' description: %s' % feature.get_description()
# Print the moving plate ID of the rotation feature.
print ' moving plate ID: %d' % moving_plate_id
# Print the fixed plate ID of the rotation feature.
print ' fixed plate ID: %d' % fixed_plate_id
# Print the time period of the rotation feature.
# This is the times of the first and last enabled rotation time samples.
enabled_time_samples = total_reconstruction_pole.get_enabled_time_samples()
if enabled_time_samples:
print ' enabled time period: %f -> %f' % (
enabled_time_samples[0].get_time(), enabled_time_samples[-1].get_time())
# Print the rotation pole information from the enabled rotation time samples.
print ' time samples:'
for time_sample in enabled_time_samples:
# Get the finite rotation from the GpmlFiniteRotation property value instead the GpmlTimeSample.
finite_rotation = time_sample.get_value().get_finite_rotation()
# Extract the pole and angle (in degrees) from the finite rotation.
pole_lat, pole_lon, pole_angle = finite_rotation.get_lat_lon_euler_pole_and_angle_degrees()
# The time and optional description come from the GpmlTimeSample.
time = time_sample.get_time()
description = time_sample.get_description()
# Print the pole data as it would appear in a PLATES4 rotation file
# (except without the moving and fixed plate IDs).
print ' %f %f %f %f %s' % (time, pole_lat, pole_lon, pole_angle, description)
Output
TotalReconstructionSequence:
description:
moving plate ID: 1
fixed plate ID: 0
enabled time period: 0.000000 -> 600.000000
time samples:
0.000000 90.000000 0.000000 0.000000 AHS-HOT Present day Atlantic-Indian hotspots fixed to 000
200.000000 90.000000 0.000000 0.000000 AHS-HOT
600.000000 90.000000 0.000000 0.000000 AHS-HOT
TotalReconstructionSequence:
description:
moving plate ID: 2
fixed plate ID: 901
enabled time period: 0.000000 -> 200.000000
time samples:
0.000000 90.000000 0.000000 0.000000 PHS-PAC Pacific Hotspots
0.780000 49.300000 -49.500000 -1.020000 PHS-PAC WK08-A Wessel & Kroenke 2008
2.580000 53.720000 -56.880000 -2.660000 PHS-PAC WK08-A Wessel & Kroenke 2008
5.890000 59.650000 -66.050000 -5.390000 PHS-PAC WK08-A Wessel & Kroenke 2008
8.860000 62.870000 -70.870000 -8.230000 PHS-PAC WK08-A Wessel & Kroenke 2008
12.290000 65.370000 -68.680000 -10.300000 PHS-PAC WK08-A Wessel & Kroenke 2008
17.470000 68.250000 -61.530000 -15.500000 PHS-PAC WK08-A Wessel & Kroenke 2008
24.060000 68.780000 -69.830000 -20.400000 PHS-PAC WK08-A Wessel & Kroenke 2008
28.280000 67.720000 -70.800000 -23.600000 PHS-PAC WK08-A Wessel & Kroenke 2008
33.540000 66.570000 -68.730000 -27.700000 PHS-PAC WK08-A Wessel & Kroenke 2008
40.100000 65.430000 -64.250000 -31.600000 PHS-PAC WK08-A Wessel & Kroenke 2008
47.910000 63.020000 -66.680000 -34.600000 PHS-PAC WK08-A Wessel & Kroenke 2008
53.350000 60.600000 -69.670000 -36.100000 PHS-PAC WK08-A Wessel & Kroenke 2008
61.100000 56.930000 -72.930000 -38.400000 PHS-PAC WK08-A Wessel & Kroenke 2008
74.500000 50.030000 -78.350000 -44.000000 PHS-PAC WK08-A Wessel & Kroenke 2008
83.500000 47.300000 -82.100000 -48.800000 PHS-PAC WK08-A Wessel & Kroenke 2008
95.000000 46.900000 -82.680000 -54.100000 PHS-PAC WK08-A Wessel & Kroenke 2008
106.200000 51.320000 -85.120000 -60.100000 PHS-PAC WK08-A Wessel & Kroenke 2008
112.300000 52.170000 -85.800000 -62.400000 PHS-PAC WK08-A Wessel & Kroenke 2008
118.400000 52.530000 -80.330000 -66.500000 PHS-PAC WK08-A Wessel & Kroenke 2008
125.000000 54.130000 -88.180000 -69.600000 PHS-PAC WK08-A Wessel & Kroenke 2008
131.900000 56.220000 -112.250000 -78.600000 PHS-PAC WK08-A Wessel & Kroenke 2008
144.000000 54.430000 -123.570000 -84.400000 PHS-PAC WK08-A Wessel & Kroenke 2008
200.000000 54.430000 -123.570000 -84.400000 PHS-PAC WK08-A Wessel & Kroenke 2008 extended rotation SZ2015
...