1 - Getting Started¶

Welcome to GPlately!

GPlately uses object-oriented programming to make life simple. In this notebook we will explore some of the main objects you will use:

  • PlateReconstruction - reconstruct features, tessellate mid ocean ridges, subduction zones
  • Points - partition points onto plates, rotate back through time
  • Raster - read in NetCDF grids, interpolation, resampling.
  • PlotTopologies - plotting topologies e.g. ridges, trenches, subduction teeth on maps
  • DataServer - downloading plate models, features, and rasters e.g. .rot, .gpml, .shp and .nc files
In [1]:
import cartopy.crs as ccrs
import gplately
import matplotlib.pyplot as plt
import numpy as np

Tectonic plate reconstructions¶

We simply supply a rotation model, plate topologies, and static polygons to initialise a plate reconstruction model. You can download these files into your machine's cache using GPlately's DataServer object.

In [2]:
# Call GPlately's DataServer object (from the download.py module) and request data from the Müller et al. 2019 study
gdownload = gplately.download.DataServer("Muller2019")
rotation_model, topology_features, static_polygons = gdownload.get_plate_reconstruction_files()

# Tessellate the subduction zones to 0.5 degrees.
tessellation_threshold_radians = np.radians(0.05)

model = gplately.PlateReconstruction(rotation_model, topology_features, static_polygons)

Checking whether the requested files need to be updated...
Requested files are up-to-date!

Now let's find the subduction zones and mid-ocean ridges at 10 Ma.

In [3]:
time = 10

# these bundle a lot of information - check PTT docs for more info
subduction_data = model.tessellate_subduction_zones(time, ignore_warnings=True)
ridge_data = model.tessellate_mid_ocean_ridges(time, ignore_warnings=True)

Mapping¶

The PlotTopologies function injests the plate model we have defined as well as the coastlines, continents, and COB. It computes all of the plate topologies for a given reconstruction time.

This object has been designed to work specifically with cartopy. Define your figure and supply your axes to these plotting routines. Some common favourites include:

  • coastlines
  • continents
  • ridges and transforms
  • trenches
  • subduction teeth (!!)
  • netCDF grids
  • plate motion vectors

You can still supply optional keywords as you normally would.

In [4]:
# Obtain features for the PlotTopologies object with DataServer
coastlines, continents, COBs = gdownload.get_topology_geometries()

# Call the PlotTopologies object
gplot = gplately.plot.PlotTopologies(model, coastlines=coastlines, continents=continents, COBs=COBs)

# Download all Muller et al. 2019 netCDF age grids with DataServer. This is returned as a Raster object.
agegrid = gdownload.get_age_grid(time)

Checking whether the requested files need to be updated...
Requested files are up-to-date!
Checking whether the requested files need to be updated...
Requested files are up-to-date!

Set the time attribute to reconstruct all topologies to the specified time.

IMPORTANT: You must set gplot.time or provide a time at initialisation before plotting anything.

In [5]:
gplot.time = 10 # Ma

Create a map with some useful geological information

In [6]:
fig = plt.figure(figsize=(16,12))

ax1 = fig.add_subplot(111, projection=ccrs.Mollweide(190))

gplot.plot_continents(ax1, facecolor='0.8')
gplot.plot_coastlines(ax1, color='0.5')
gplot.plot_ridges_and_transforms(ax1, color='red')
gplot.plot_trenches(ax1, color='k')
gplot.plot_subduction_teeth(ax1, color='k')
im = gplot.plot_grid(ax1, agegrid.data, cmap='YlGnBu', vmin=0, vmax=200)
gplot.plot_plate_motion_vectors(ax1, spacingX=10, spacingY=10, normalise=True, zorder=10, alpha=0.5)

fig.colorbar(im, orientation='horizontal', shrink=0.4, pad=0.05, label='Age (Ma)')
Out[6]:
<matplotlib.colorbar.Colorbar at 0x17be40940>
No description has been provided for this image
In [7]:
# update the time to regenerate topologies
time = 100
gplot.time = time
agegrid = gdownload.get_age_grid(time)

fig = plt.figure(figsize=(16,12))

ax1 = fig.add_subplot(111, projection=ccrs.Mollweide(190))

gplot.plot_continents(ax1, facecolor='0.8')
gplot.plot_coastlines(ax1, color='0.5')
gplot.plot_ridges_and_transforms(ax1, color='red')
gplot.plot_trenches(ax1, color='k')
gplot.plot_subduction_teeth(ax1, color='k')
im = gplot.plot_grid(ax1, agegrid.data, cmap='YlGnBu', vmin=0, vmax=200)
gplot.plot_plate_motion_vectors(ax1, spacingX=10, spacingY=10, normalise=True, zorder=10, alpha=0.5)

fig.colorbar(im, orientation='horizontal', shrink=0.4, pad=0.05, label='Age (Ma)')

Checking whether the requested files need to be updated...
Requested files are up-to-date!
Out[7]:
<matplotlib.colorbar.Colorbar at 0x17a3f8f40>
No description has been provided for this image

Working with points¶

Now that we have defined our reconstruction object, we can reconstruct point data.

In [8]:
pt_lons = np.array([140., 150., 160.])
pt_lats = np.array([-30., -40., -50.])

gpts = gplately.Points(model, pt_lons, pt_lats)

vel_x, vel_y = gpts.plate_velocity(0)
vel_mag = np.hypot(vel_x, vel_y)

print("average point velocity (cm/yr)", vel_mag)

average point velocity (cm/yr) [6.38389518 5.71007208 5.06602802]

Plot their position from time=0 to time=20

In [9]:
rlons = np.empty((21, pt_lons.size))
rlats = np.empty((21, pt_lons.size))


for time in range(0, 21):
    rlons[time], rlats[time] = gpts.reconstruct(time, return_array=True)
In [10]:
gplot.time = 0 # present day

fig = plt.figure(figsize=(16,12))
ax1 = fig.add_subplot(111, projection=ccrs.Mercator(190))
ax1.set_extent([130,180,-60,-10])

gplot.plot_coastlines(ax1, color='0.8')

for i in range(0, len(pt_lons)):
    ax1.plot(rlons[:,i], rlats[:,i], 'o', transform=ccrs.PlateCarree())
No description has been provided for this image

Rasters¶

You can initialise a Raster object by providing a raster (either a NetCDF file path or a regular numpy array), and optionally passing a PlateReconstruction model. The data attribute stores the raster data as a 2D numpy array.

In this example, we will pass a NetCDF file path to the Raster object.

In [11]:
time = 100

# download a netcdf age grid Raster and assign a plate reconstruction `model'
graster = gdownload.get_age_grid(time)
graster.plate_reconstruction = model

# the underlying numpy (masked) array can be accessed using the data attribute
print(type(graster.data))

# alternatively initialise a Raster from a numpy array
graster = gplately.Raster(data=graster.data, # 2D numpy array
                          plate_reconstruction=model, # PlateReconstruction object
                          extent='global', # equivalent to [-180, 180, -90, 90]
                          time=100, # time in Ma
                         )

Checking whether the requested files need to be updated...
Requested files are up-to-date!
<class 'numpy.ma.core.MaskedArray'>
In [12]:
gplot.time = time

fig = plt.figure(figsize=(16,12))

ax1 = fig.add_subplot(111, projection=ccrs.Robinson())
gplot.plot_grid(ax1, graster, cmap='YlGnBu', vmin=0, vmax=100)
# Alternatively:
# graster.imshow(ax=ax1, cmap="YlGnBu", vmin=0, vmax=100)
gplot.plot_coastlines(ax1, edgecolor='k', facecolor='none')
Out[12]:
<GeoAxesSubplot:>
No description has been provided for this image

There are a bunch of routines such as,

  • filling masked (NaN) regions
  • interpolation
  • resampling
  • reconstructing rasters

In-place operations can be achieved using inplace=True which will update the internal data structures.

In [ ]: