Open-source cartography with GMT

Recently, I wanted to draw a map to illustrate a presentation, but I couldn’t find anything satisfactory online. Specifically, I wanted to highlight the continental shelves surrounding the world’s landmasses, but the pictures that I found had garish colors, or were general bathymetric charts, both of which I thought were distracting from what I wanted to show.

After looking around, I found a great resource for bathymetric data: GEBCO, or the General Bathymetric Chart of the Oceans. This is a publicly-available resource that compiles bathymetry from various sources into a single grid. The chart is available at two resolutions: 30 arc-seconds and 1 arc-minute. So I downloaded the 1-minute grid, but what next? How to plot these data?

This is where GMT, or Generic Mapping Tools comes in. It’s an open-source project with a long history – the first version was released in 1988! They’re command line tools (originally for Unix, but also available for Windows and Mac OS), so they take some getting used to. However, the flexibility that they offer for producing customized maps and map overlays is fantastic.

The best way to learn how to use it is to download the tools and go through the tutorial. Another useful online tutorial is here.

Some tips which weren’t in the tutorial:

  • Instead of downloading the repository via SVN, can also just download the .tar.gz archive from the homepage or one of its mirrors.
  • GSHHG data (the coastline data) has to be placed manually in the data share directory:
    • Download the GSSHG data – a separate .tar.gz archive called gshhg-gmt-XXX.tar.gz where XXX is the version number.
    • After installing gmt, type gmt –show-datadir to show where gmt stores its shared data.
    • Unpack the gshhg archive, and rename that folder to coast
    • Place coast folder within the shared data folder for gmt, where it can find it.
    • You’ll need these data for running the pscoast command (in tutorial section 2.4.3 onwards)
  • When making overlays, use the -K parameter to keep the Postscript file “open” for further editing, otherwise any further overlay commands will not take effect.
  • Color tables (for generating color scales in plots) are in CPT files – these can be generated automatically with makecpt, as described in the tutorial, but the format is pretty easy to figure out so you can manually edit your own scales too. Just generate a dummy scale with makecpt, and edit the depth/height ranges and color codes yourself. Color scales in CPT files are in HSV format (instead of RGB), and instead of percent values for saturation and value, they take decimal fractions.

This makes it all sound more complicated than it actually is. To draw my map, here’s what I did:

  • I downloaded the GEBCO bathymetry data (you have to agree to a end-user license etc. so I’m not directly linking here). This is a file called GRIDONE_1D.nc.
  • Draw up a custom CPT file (called cont_shelf.cpt), which lists what colors I want for each height/depth range:

# COLOR_MODEL = hsv
-10600  250-1-0.78      -200    250-1-0.78
-200    194-0.29-1      0       194-0.29-1
0       45-0.37-1       8430    45-0.37-1
B       0-0-0
F       0-0-1
N       0-0-0.5

Simply put: Everything from -10600 m to -200 m depth is dark blue, from -200 to 0 m depth is light blue, and from 0 m to 8430 height (i.e. dry land) is a sickly yellow-brown.

  • Draw the base map with gmt grdimage:
    • gmt grdimage GRIDONE_1D.nc -JKf0/9i -R-180/180/-90/90 -Ba -Ccont_shelf.cpt -V -K > worldmap.ps
    • -J sets the map projection. “Kf” is the Eckert IV equal-areas projection, centered on “0” (the Greenwich meridian), and the plot has a width of 9 inches (“9i”)
    • -R specifies the map boundaries: from -180 W to 180 E, and from -90 S to 90 N.
    • -C parameter: Use the color scale contained in the file cont_shelf.cpt
    • -V verbose output
    • -K Keep Postscript file open
  • I had a list of points that I wanted to plot. This was in a file (locations.tab) with latitude and longitude (as decimal degrees), one pair per line.
  • Overlay the points with gmt psxy
    • gmt psxy -O -R-180/180/-90/90 -JKf0/9i -Sc0.05i -Gred -: locations.tab >> worldmap.ps
    • -R and -J arguments should match the arguments given to the original plot (the gmt grdimage command before)
    • -O means “overlay” mode
    • -S describes the overlay type. “c” means draw circles, and “0.05i” means the circles are of width 0.05 inches
    • -G describes color of the overlay, in this case “red”
    • -: means the coordinates supplied are in lat-long format. The default input format is long-lat.
    • Because no -K parameter is given, the file is then finalized and further commands take no effect.
  • Voila!

worldmap3_bak (copy)

Advertisements

3 comments on “Open-source cartography with GMT

  1. […] have previously blogged about an open source cartography software called GMT (Generic Mapping Tools). This tool is great […]

  2. […] gmt, which I’ve blogged about before, it’s a simple shell script to produce a set of orthographic projections (parameters adapted […]

  3. […] plugin of the dataviz library D3.js by its author Mike Bostock. Compared to gmt, which I’ve blogged about before, there are more classes of projections available if you use the plugin, and it’s designed for […]

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s