GMT Study Notes

Simple tricks to make scripting in GMT manageable

Posted by heapify on October 15, 2020

Timeline

  • 2017-02-01: created and shared on OverLeaf
  • 2017-02-20: peer reviewed by Sabber Ahamed and Chunyu Liu
  • 2020-10-15: hosted on GitHub and rendered as Jekyll pages

Preface

These techniques aim to improve your working efficiency of producing figures with Generic Mapping Tools. Although the document uses GMT 4 for demonstrations, it is not difficult for readers to follow the ideas and rewrite them in GMT 5 or higher. The author reserves the copyright of any figures in the examples (if there are any). Full scripts for producing these figures will not be released to the public.

execution_example

Preparatory knowledge

A title goes here

This section is still under rapid development…

Start working

Disposal of temporary files

In a timely manner, clean up *.intns and *.grd files after calling grdimage and grdcontour and *.cpt files after calling psxy and psscale, such as

1
2
3
echo "Recycling wastes..."
rm *.cpt          # used by psxy, psxyz and psscale
rm *.intns *.grd  # used by grdimage and grdcontour

These lines of code clean up used files in the working directory to maintain a tidy environment. There is no consequence if you do not do it. However, the following is more than just a suggestion.

At either the beginning or the end of your GMT script, use

1
2
3
echo "Cleaning GMT configurations..."  # modifications for gmt 5 needed
[[ -e .gmtcommands4 ]] && rm .gmtcommands4
[[ -e .gmtdefaults4 ]] && rm .gmtdefaults4

to remove temporary files. These files are generated once gmtset is invoked to modify page parameters. Omitting this step could be dangerous, for that the same script may produce different plots in different environments. Once it happens, it is usually challenging to troubleshoot and locate the issue.

Using these bash commands at the beginning of a script cleans up the configurations file left by other scripts in the same directory to protect your plots from being affected. Using them at the end of the script cleans up the configuration files generated by the current script to not affect other plots. Executing GMT scripts is a repetitive work that usually requires many trial-and-errors. Using these commands at the beginning and the end of a script virtually plays the same role.

Moreover, it is never recommended to share the same configuration files with multiple scripts.

Portrait Mode in a line of code

If you prefer to use the Portrait Mode to create a plot, use

1
gmtset PAGE_ORIENTATION portrait

to specify the page orientation. Make sure these lines appear after cleaning the GMT configurations. They are equivalent to using a -P argument for each plotting command. Using the Portrait Mode is much more convenient than using the -P argument many times.

Safe usage of -R and -J arguments

Declaring variables in a bash script is not only to improve the readability and portability. The significance of using variables is far more than these!

In general, study area and projection method just need to be provided once in a GMT script. According to the GMT syntax requirements, these parameters could be omitted after the first invoke, as demonstrated below:

1
2
3
4
5
6
7
8
9
10
echo "Plotting shorelines..."
pscoast -JM6i -R124/132/33/39 -Ir -N1 -Di \
    -Ba2f1:."Busan Earthquakes Epicenter and Station Map": \
    -Wdarkgrey -Xc -Yc -K > map_names.ps
    # arguments for -R and -J are specified

echo "Projecting grids..."
grdimage pattopobath.grd -R -J \
    -Cpattopo.cpt -Ipattopo.intns -O -K >> map_names.ps
    # arguments for -R and -J are omitted

However, this turns out to be a lazy omission and can cause issues when the map maker needs to move the first plotting command downward but forgets to modify the -R and -J arguments. With no parameters given for the first invoke, GMT will look for the parameters of these arguments from their last execution, to be more specific, the previous time you run a GMT script or any GMT calls involving -R and -J.

1
2
3
4
5
6
7
8
9
10
echo "Projecting grids..."
grdimage pattopobath.grd -R -J \
    -Cpattopo.cpt -Ipattopo.intns -Xc -Yc -K > map_names.ps
    # wrong, -R and -J will take the
    # parameters from their last call

echo "Plotting shorelines..."
pscoast -JM6i -R124/132/33/39 -Ir -N1 -Di \
    -Ba2f1:."Busan Earthquakes Epicenter and Station Map": \
    -Wdarkgrey -O -K >> map_names.ps

For safety, it is strongly recommended to save these parameters in variables and provide them for each -R and -J invoke. With these good habits developed, the arrangement of code sections might become much more straightforward.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
# declare some variables
J=M6i
R=124/132/33/39
ps=map_names.ps

echo "Projecting grids..."
grdimage pattopobath.grd -J$J -R$R \
    -Cpattopo.cpt -Ipattopo.intns -Xc -Yc -K > $ps
    # arguments for -R and -J are specified for each invoke

echo "Plotting shorelines..."
pscoast -R$R -J$J -Ir -N1 -Di \
    -Ba2f1:."Busan Earthquakes Epicenter and Station Map": \
    -Wdarkgrey -O -K >> $ps
    # safe for switching lines

Be aware that you still have to be careful about the usages of -O, -K, > and >>. Without upgrading to the latest GMT 6, there are thus far no easy ways to get rid of these troubles.

Advanced techniques

Avoid hard-coded filenames

In a GMT script, the output filename can be either specified explicitly or through a predefined variable implicitly to guide each plotting command. Doing this explicitly requires that the filename must be consistent throughout the script. Inconsistent filenames always result in broken PostScript documents because of making a mistake. In addition to saving the filename in a variable as one of the feasible solutions, we are looking for a method that outputs the figure into a file with the same name as your script so that you will no longer be bothered.

To do this elegantly, consider introducing a little bit usage of bash variables. Suppose a script is named example.gmt, and you would like to guide each GMT command to example.ps. A preliminary attempt could be

1
> $0.ps  # >> $0.ps for writing to open files

as GMT standard output. This will write the figure into example.gmt.ps. The \$0 variable gives you the full filename, and you have to keep the original file extension in the output filename.

To get rid of the file extension, an improvement is to declare several variables

1
2
3
# define filename
file=`echo $0 | sed -E 's/^(.*)\..*$/\1/g'`
filename="$file.ps"

somewhere at the beginning of the script; then, use

1
> $filename  # >> $filename for writing to open files

to guide each plotting command. This will write the figure into example.ps instead of example.gmt.ps.

For filenames containing multiple dots, only the last extension is treated as the “real” file extension, as illustrated in the flow chart below:

1
2
3
## e.g., $0 = ex.ample.gmt
## --> $file = ex.ample
## --> $filename = ex.ample.ps

Automatic document reloading

Most GMT users in Mac OS X find it annoying to repeatedly close and open a PostScript file in order to preview updated content. What is even more bothering, double-clicking in Finder a PostScript file currently being shown does not update the plot in the current window. The Preview app will always activate a new window.

As an explanation of the terminology, applications that possess an automatic document reloading feature perform the following actions:

  • The figure previewing program is capable of detecting changes; by the time you modify something and re-execute your GMT script, the figure in the previewing window is reloaded.
  • Further, the figure will be reloaded in the current window, not a newly activated window; more concisely, reloading the figure does not create new windows.

For Microsoft Windows, a famous light-weighted PostScript viewer called Sumatra PDF supports automatic document reloading - the ideal choice for GMT and LaTeX users who frequently need to refresh their documents. For Linux, there is a PostScript viewer with similar features named Zathura. The only difficulty thus far is to introduce automatic document reloading for Mac OS.

Looking for an alternative to Mac’s Preview is not an easy job. Several options have been tried, but they are all problematic (either slow or lack of software maintenance), including the GMT officially recommended Ghostview. Without using third-party software, is there a method to automatically reload a document when it is modified?

The built-in figure previewing program in Mac OS, Preview, locally support automatic document reloading. The only reason it does not update PostScript files is the conversion from PS to PDF before opening - what is being displayed is a converted PDF file somewhere in a temporary directory - the same reason for that repeatedly open the same figure will activate new windows. While GMT does not support different output formats other than PostScript, it provides some conversion tools for you to obtain other formats. Having these working mechanisms in mind, we have come up with an idea of realizing automatic document reloading with Linux Ghost Script converter or the conversion tools coming with GMT.

Take the PNG format as an example. Add these lines of codes at the end of a script after all the plotting and cleaning-up. The definitions of $file and $filename here remain the same as in the last section.

1
2
3
4
5
6
echo "Reloading the document..."
gs -dNOPAUSE -dBATCH -sDEVICE=pngalpha \
    -sOutputFile=$file.png $filename
open $file.png

echo "Done."

The usage of open here seems unnecessary, but it brings the Preview window in front of the terminal for Preview to reload the figure (the Preview program only reloads a figure when its window is activated). Because the conversion is not done in the Preview program, as long as the UNIX command open is called, the figure is updated in the current window instead of in a newly-created window.

Alternatively, use the conversion tool provided by GMT

1
2
3
4
5
echo "Reloading the document..."
ps2raster $filename -A -P -Tg
open $file.png

echo "Done."

to convert the file to PNG using a tight BoundingBox and rotating it back to normal orientation in case it was in Landscape mode. In GMT 5, ps2raster can also be replaced by psconvert for the same purpose.

Two ways of creating contour plots

grdcontour and pscontour are both GMT commands for plotting contours. Some differences are briefly summarized below:

  • pscontour interpolates data through an optimal triangulation method; grdcontour obtains contours from the surface created with surface or xyz2grd.
  • pscontour creates contours directly from the original data; while to use grdcontour, you have to plot a surface with grdimage as prior work.
  • The contour lines created with grdimage are smoother than those created with pscontour related to their interpolation methods. pscontour does not support extrapolation; the extrapolation of grdimage is default on and there is no native option to switch it off. This is also related to interpolation methods.

To fill the spaces with colors in between contour lines, a -C option appended to pscontour specifying a *.cpt file is needed. There are multiple ways of creating color palette tables, manually or automatically. When the range of data is clear, makecpt is a way to do it manually; xyz2grd with grd2cpt creates a color palette table automatically according to the maximum and minimum values in the data. Furthermore, creating a text file and entering a table from scratch is even possible, as an example[1] shows:

1
2
3
4
echo "Getting color palette table file..."
echo -8000 255 15 255 -7000 255 15 255  > twallnew.cpt
echo -7000 255 15 255 -6000 255 75 255 >> twallnew.cpt
echo -6000 255 75 255 -5000 255 75 225 >> twallnew.cpt

[1] from Jyr-Ching HU, Department of Geosciences, National Taiwan University

Using the RGB system, the format of the cpt-file is:

1
2
3
z0 Rmin Gmin Bmin z1 Rmax Gmax Bmax [A]
...
zn-2 Rmin Gmin Bmin zn-1 1 Rmax Gmax Bmax [A]

The colors may be specified either in the RGB (Red, Green, Blue) system or in the HSV system (Hue, Saturation, Value), and the parameter COLOR MODEL in the .gmtdefaults file must be set accordingly.

More differences between gridding methods are demonstrated in the official documentation GMT historical collection.


Appendix

A function for reporting status to the user:

1
2
3
4
5
6
7
function status {
  case $1 in
    0 ) printf "\e[0;32mOK\e[0m\n" ;; # 0 = OK
    1 ) printf "\bskip\n" ;;          # 1 = skip
    * ) printf "%-42s" "$*" ;;        # otherwise, print message
  esac
}

Here is how to use:

1
2
3
4
5
6
7
8
9
10
11
12
13
status "Cleaning temporary files..."
[[ -e .gmtcommands4 ]] && rm .gmtcommands4
[[ -e .gmtdefaults4 ]] && rm .gmtdefaults4
[[ -e  gmt.history  ]] && rm  gmt.history
status 0

status "Initializing parameters..."  # all fixed for gmt 5
gmt set FONT_ANNOT_PRIMARY 10p       # legend and axes
gmt set MAP_FRAME_TYPE fancy
gmt set FONT_LABEL 12p  # color scale title
gmt set FONT_TITLE 14p  # super plot title
gmt set PS_PAGE_ORIENTATION portrait
status 0

In the Terminal, you will see

1
2
Cleaning temporary files...               OK
Initializing parameters...                OK

The status 0 statements print “OK” texts in green, showing the completion of each code section.

List of GMT 4 page settings that are deprecated in GMT 5:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
ANNOT_FONT_PRIMARY          INTERPOLANT
ANNOT_FONT_SECONDARY        LABEL_FONT
ANNOT_FONT_SIZE_PRIMARY     LABEL_OFFSET
ANNOT_FONT_SIZE_SECONDARY   LINE_STEP
ANNOT_MIN_ANGLE             MAP_SCALE_FACTOR
ANNOT_OFFSET_PRIMARY        MEASURE_UNIT
ANNOT_OFFSET_SECONDARY      NAN_RECORDS
BASEMAP_AXES                OBLIQUE_ANNOTATION
BASEMAP_FRAME_RGB           OUTPUT_CLOCK_FORMAT
BASEMAP_TYPE                OUTPUT_DATE_FORMAT
CHAR_ENCODING               OUTPUT_DEGREE_FORMAT
D_FORMAT                    PAGE_COLOR
DEGREE_SYMBOL               PAGE_ORIENTATION
ELLIPSOID                   PAPER_MEDIA
FIELD_DELIMITER             PLOT_CLOCK_FORMAT
FRAME_PEN                   PLOT_DATE_FORMAT
FRAME_WIDTH                 PLOT_DEGREE_FORMAT
GLOBAL_Y_SCALE              POLAR_CAP
GLOBAL_X_SCALE              PS_COLOR
GRID_CROSS_SIZE_PRIMARY     TICK_LENGTH
GRID_CROSS_SIZE_SECONDARY   TICK_PEN
GRID_PEN_PRIMARY            TIME_FORMAT_PRIMARY
GRID_PEN_SECONDARY          TIME_FORMAT_SECONDARY
GRIDFILE_FORMAT             UNIX_TIME_FORMAT
GRIDFILE_SHORTHAND          UNIX_TIME_POS
HEADER_FONT_SIZE            UNIX_TIME
HEADER_FONT                 VECTOR_SHAPE
HEADER_OFFSET               VERBOSE
HISTORY                     WANT_LEAP_SECONDS
HSV_MAX_SATURATION          X_ORIGIN
HSV_MAX_VALUE               XY_TOGGLE
HSV_MIN_SATURATION          Y_AXIS_TYPE
HSV_MIN_VALUE               Y_ORIGIN
INPUT_CLOCK_FORMAT          Y2K_OFFSET_YEAR
INPUT_DATE_FORMAT