6 - Working with Rasters

Any numpy array can be turned into a gplately.Raster:

raster = gplately.Raster(
    plate_reconstruction=model,
    data=array,
    extent="global",  # equivalent to (-180, 180, -90, 90)
    origin="lower",  # or set extent to (-180, 180, -90, 90)
)

In this notebook, we will:

• Download rasters from EarthByte's webDAV server
• Plot rasters
• Resize and respace rasters
• Reconstruct rasters
• Linearly interpolate point data on rasters.

Import all needed packages, and create PlateReconstruction and Plot objects for the Muller et al. (2019) plate tectonic model.

import os

import cartopy.crs as ccrs
import gplately
import matplotlib.pyplot as plt
import numpy as np

Using GPlately's download module, we can easily download a variety of plate tectonic models!

gdownload = gplately.download.DataServer("Muller2019")
rotation_model, topology_features, static_polygons = gdownload.get_plate_reconstruction_files()
coastlines, continents, COBs = gdownload.get_topology_geometries()

model = gplately.PlateReconstruction(rotation_model, topology_features, static_polygons)
gplot = gplately.PlotTopologies(model, coastlines, continents, COBs)

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!

GPlately's DataServer for rasters

Let's use GPlately's DataServer object to download Muller et al. 2019 netCDF age grids.

There is a unique age grid for each millionth year - let's access the 0 Ma age grid by passing time to get_age_grid. It is returned as a gplately Raster object which we call muller_2019_age_grid.

time = 0  # Ma
muller_2019_age_grid = gdownload.get_age_grid(time)

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

Plotting rasters

The muller_2019_age_grid is a Raster object - this object allows us to work with age grids and other rasters. Let's visualise the data with imshow.

fig = plt.figure(figsize=(16, 12))
muller_2019_age_grid.imshow(cmap="YlGnBu")
plt.show()

Let's plot this netCDF grid along with coastlines, mid-ocean ridges and subduction zones (with teeth) resolved from the Muller et al. 2019 plate model.

fig = plt.figure(figsize=(16, 12))
ax1 = fig.add_subplot(111, projection=ccrs.Mollweide(central_longitude=20))

gplot.time = time

muller_2019_age_grid.imshow(ax1, cmap='YlGnBu', vmin=0, vmax=200)
gplot.plot_coastlines(ax1, edgecolor='k', facecolor='1', alpha=0.1)
gplot.plot_trenches(ax1, color='magenta', zorder=5)
gplot.plot_ridges(ax1, color='b', zorder=5)
gplot.plot_transforms(ax1, color='lightblue', linewidth=0.75)
gplot.plot_subduction_teeth(ax1, color='magenta')
ax1.set_title("{} Ma".format(time))
plt.show()

Resizing and resampling rasters

Let's resize and resample the the present-day Muller et al. 2019 agegrid. We can do this using resize.

• resize - provide the number of points in each direction to resize the raster (e.g. 100 cols by 200 rows)
• resample - provide the grid spacing in each direction (e.g. 0.1 degrees by 0.2 degrees)

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

muller_2019_age_grid = gdownload.get_age_grid(time)

# Set grid size in x and y directions
muller_2019_age_grid.resize(20, 10, inplace=True)

# Plot resampled and resized age grid
fig = plt.figure(figsize=(16,12))
ax1 = fig.add_subplot(111, projection=ccrs.Mollweide(central_longitude = 20))
muller_2019_age_grid.imshow(ax1, cmap='YlGnBu', vmin=0, vmax=200)
gplot.plot_coastlines(ax1, edgecolor='k', facecolor='1', alpha=0.1)
gplot.plot_trenches(ax1, color='magenta', zorder=5)
gplot.plot_ridges(ax1, color='b', zorder=5)
gplot.plot_subduction_teeth(ax1, color='magenta')
ax1.set_title("{} Ma".format(time))
plt.show()

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