Thank you for giving GMT/Python and GMT6 a try!
This is a Jupyter notebook. It's an interactive computing environment where you can mix text (like this), code, and figures. The notebook is organized into cells. This is a Markdown cell (click on it to see the source) and it can contain text, hyperlinks, images, and even Latex equations.
To execute the code cells, select it and type shift+enter
or click on the Run button above.
The following are examples of what you can currently do with GMT/Python. We'll also make the image for the poster background. There are some empty cells for you to experiment on the bottom of the notebook along with an example using the command line GMT6 in the new "modern mode".
The GMT modules are available as functions and classes in the gmt
Python package.
So we'll start by importing it:
import gmt
All figure generation in GMT/Python is handled by the gmt.Figure
class.
It has methods to add layers to your figure, like a basemap, coastlines, and data.
We start a new figure by creating an instance of gmt.Figure
:
fig = gmt.Figure()
We add elements to the figure using its methods. For example, lets add the coastlines of Central America to a 6 inch wide map using the Mercator projection (M
). Our figure will also have a nice frame with automatic ticks.
fig.coast(region=[-90, -70, 0, 20], projection='M6i', land='chocolate',
frame=True)
You can see a preview of the figure directly in the Jupyter notebook using fig.show()
.
fig.show()
You'll probably have noticed several things that are different from classic command-line GMT. Many of these changes reflect the new GMT modern execution mode that will be part of the future 6.0 release. A few are GMT/Python exclusive (like the long argument names).
coast
instead of pscoast
. As a general rule, all ps*
modules had their ps
removed. The exceptions are: psxy == plot
, psxyz == plot3d
, and psscale == colorbar
.region
can take lists instead of strings like 1/2/3/4
. You can still use the string form but the list form is easier in Python.-B
instead of -Baf
), use a True
value instead. An empty string would also be acceptable.We could have generated the figure above using the classic GMT argument names (but not the module names):
fig_alias = gmt.Figure()
fig_alias.coast(R='-90/-70/0/20', J='M6i', G='chocolate', B=True)
fig_alias.show()
Unlike the GMT command-line interface, no figure file was generated until you ask for one.
Use method fig.savefig
(based on the matplotlib function) to save your figure to a file:
fig.savefig('first-steps-central-america.png')
If you're running a Python script, you can tell fig.savefig
to open the figure in an external viewer:
fig.savefig('first-steps-central-america.png', show=True)
We can use gmt.Figure.plot
to plot data on our map.
First, lets create some fake data to plot using numpy:
import numpy as np
# See the random number generator so we always
# get the same numbers
np.random.seed(42)
ndata = 30
region = [150, 240, -30, 60]
lon = np.random.uniform(region[0], region[1], ndata)
lat = np.random.uniform(region[2], region[3], ndata)
magnitude = np.random.uniform(1, 9, ndata)
depth = np.random.uniform(0, 1, ndata)
There are 3 ways to pass data into Figure.plot
:
x
and y
arguments.data
argument.data
argument.Let's explore all of these options.
Now we can plot the data using Figure.plot
and passing the x and y coordinate arrays:
fig = gmt.Figure()
fig.coast(region='g', projection='G200/30/6i', frame='ag',
resolution='i', area_thresh=5000, land='white',
water='DarkTurquoise')
# Plot using circles (c) of 0.3 cm
fig.plot(x=lon, y=lat, style='c0.3c', color='black')
fig.show()
We can make the size of the markers follow a data value by passing in the argument sizes
:
fig = gmt.Figure()
fig.coast(region='g', projection='G200/30/6i', frame='ag',
resolution='i', area_thresh=5000, land='white',
water='DarkTurquoise')
fig.plot(x=lon, y=lat, sizes=0.1*magnitude, style='cc', color='black')
fig.show()
We can also map the colors of the markers to data by passing color
as an array and providing a colormap name (cmap
). We can even use the new matplotlib colormap "viridis":
fig = gmt.Figure()
fig.coast(region='g', projection='G200/30/6i', frame='ag',
resolution='i', area_thresh=5000, land='white',
water='DarkTurquoise')
fig.plot(x=lon, y=lat, style='cc', sizes=0.1*magnitude,
color=depth, cmap='viridis')
fig.show()
Sometimes, you'll have data that you loaded from a file in the form of a single numpy 2d array (like a table). You can plot these data by passing in the data
argument.
table = np.transpose([lon, lat, magnitude, depth])
We can use the columns
argument to specify the columns that correspond to x, y, color, and size, respectively. We'll use it to specify that we want to set the colors using the depths:
fig = gmt.Figure()
fig.coast(region=region, projection='M6i', frame='ag',
resolution='i', area_thresh=5000, land='grey')
fig.plot(data=table, style='c0.8c', cmap='viridis', columns='0,1,3')
fig.show()
If you have data in a file and you don't need to do any calculations on it, you can pass the file name to plot
directly. The syntax for plotting will be the same as the example for a data matrix.
First, we must save our sample data to a file:
np.savetxt('first-steps-data.txt', table)
fig = gmt.Figure()
fig.coast(region='g', projection='N-165/6i', frame='ag',
resolution='i', area_thresh=5000, land='chocolate')
fig.plot(data='first-steps-data.txt', style='c0.4c',
cmap='viridis', columns='0,1,3')
fig.show()
GMT ships some sample data that can be accessed by passing @data_file
as the data
argument. For example, we can plot the earthquake epicenters from the tut_quakes.ngdc
sample dataset:
fig_quakes = gmt.Figure()
fig_quakes.coast(
region=[130, 150, 35, 50], projection='M6i',
frame='afg', shorelines=True, land='gray',
water='lightblue')
fig_quakes.plot(data='@tut_quakes.ngdc', style='c0.5c',
color='blue', pen='faint', columns=[4, 3])
fig_quakes.show()
I'll use the USGS quake data that comes with GMT. I downloaded it and saved to a numpy-friendly format in file usgs_quakes.txt
.
lon, lat, magnitude = np.loadtxt('usgs_quakes.txt', unpack=True)
Make a global Mercator map without any borders. Plot the quakes mapping the marker size and color to the magnitude.
fig = gmt.Figure()
fig.coast(region=[-270, 90, -70, 70], projection='M10i', land='#aaaaaa',
water='white', resolution='l')
fig.plot(lon, lat, style='cc', sizes=0.005*2**magnitude,
color=magnitude/magnitude.max(), cmap='ocean')
fig.show(width=900)
If it looks OK, save it to a high resolution PNG.
fig.savefig('poster_background.png', dpi=1200)
Type anything you want in the cells below.
You can even try the GMT6 command line programs. The following are some of the plots from Paul's talk "The Generic Mapping Tools 6: Classic versus Modern Mode".
Cells that start with %%bash
run their code through the bash shell (use them to execute GMT commands). We'll use the IPython.display.Image
class to display in the notebook the PNG figures that you generate on the %%bash
cells.
from IPython.display import Image
%%bash
gmt begin Chile png
gmt pscoast -RCL+r2 -JM15c+ -BWSne -B -Gbeige -Sblue -N1,1p
gmt end
Image('Chile.png', width=120)
%%bash
gmt begin map png
gmt grdimage @earth_relief_05m -RMG+r2 -Cgeo -I+
gmt coast -Wthin -BWSne -B
gmt colorbar -DJTC -B -C+Uk
gmt end
Image('map.png', width=300)
%%bash
gmt begin islands png
gmt set MAP_ANNOT_OBLIQUE 34 FORMAT_GEO_MAP ddd:mmF
gmt subplot begin 2x2 -M0.05i -Fs3i/3i -BWSne -A+gwhite+p0.5p
gmt grdimage @earth_relief_03s -R-30/30/-30/30+uk -JA159:32W/22:03N -Ba20mf5m -Csrtm -c1,1
gmt grdimage @earth_relief_03s -R-15/15/-15/15+uk -JA109:20W/27:07S -Ba20mf5m -c1,2
gmt grdimage @earth_relief_03s -R-30/30/-30/30+uk -JA149:22W/17:43S -Ba20mf5m -c2,1
gmt grdimage @earth_relief_03s -R-10/10/-10/10+uk -JA138:39W/10:29S -Ba20mf5m -c2,2
gmt subplot end
gmt colorbar -B -Dx0/0.4i+jBC+w5i+h -Xw/2 -Yh
gmt end
Image('islands.png', width=500)
%%bash
gmt begin
gmt figure hotspots png
gmt grdimage @earth_relief_10m -JG200/30/6i -Cgeo -I+
gmt coast -W -Dc -Bafg
gmt plot @hotspots.txt -Sc0.2c -Gred
gmt end
Image('hotspots.png', width=500)
Try running different things in the cells below or edit and rerun the cells above.
%%bash
Image('XXX.png', width=500)
%%bash
Image('XXX.png', width=500)
%%bash
Image('XXX.png', width=500)