(22) World-wide seismicity the last month¶
The next example uses the command-line tool curl to obtain a data
file from a specified URL [1]. In the example script this line is
commented out so the example will run even if you do not have wget
(we use the remote file usgs_quakes_22.txt
which normally would be created by curl);
remove the comments to get the actual current seismicity plot using the
live data. The main purpose of this script is not to show how to plot a
map background and a few circles, but rather demonstrate how a map
legend may be composed using the module
legend. Some scripting is used to
pull out information from the data file that is later used in the
legend. The legend will normally have the email address of the script
owner; here that command is commented out and the user is hardwired to
“GMT guru”. The USGS logo, taken from their web page and converted to a
Sun raster file, is used to spice up the legend.
The script produces the plot in Figure, giving the URL where these and similar data can be obtained.
#!/usr/bin/env bash
# GMT EXAMPLE 22
#
# Purpose: Automatic map of last month of world-wide seismicity
# GMT modules: coast, plot, legend
# Unix progs: cat, sed, awk, wget|curl
#
gmt begin ex22
# Get the data (-s silently) from USGS using the curl
# Hardwired here to the month of October, 2017
# SITE="https://earthquake.usgs.gov/fdsnws/event/1/query.csv"
# TIME="starttime=2017-09-01%2000:00:00&endtime=2017-10-01%2000:00:00"
# MAG="minmagnitude=3"
# ORDER="orderby=magnitude"
# URL="${SITE}?${TIME}&${MAG}&${ORDER}"
# curl -s $URL > usgs_quakes_22.txt
# Count the number of events (to be used in title later. one less due to header)
file=$(gmt which @usgs_quakes_22.txt -G)
n=$(gmt info $file -h1 -Fi -o2)
# Pull out the first and last timestamp to use in legend title
first=$(gmt info -h1 -f0T -i0 $file -C --TIME_UNIT=d -I1 -o0 --FORMAT_CLOCK_OUT=-)
last=$(gmt info -h1 -f0T -i0 $file -C --TIME_UNIT=d -I1 -o1 --FORMAT_CLOCK_OUT=-)
# Assign a string that contains the current user @ the current computer node.
# Note that two @@ is needed to print a single @ in gmt text:
#set me = "$user@@$(hostname)"
me="GMT guru @@ GMTbox"
# Create standard seismicity color table
gmt makecpt -Cred,green,blue -T0,100,300,10000 -N
# Start plotting. First lay down map, then plot quakes with size = magnitude * 0.015":
gmt coast -Rg -JK180/22c -B45g30 -B+t"World-wide earthquake activity" -Gburlywood -Slightblue -A1000 -Y7c
gmt plot -C -Sci -Wfaint -hi1 -i2,1,3,4+s0.015 $file
# Create legend input file for NEIS quake plot
cat > neis.legend <<- END
H 16p,Helvetica-Bold $n events during $first to $last
D 0 1p
N 3
V 0 1p
S 0.25c c 0.25c red 0.25p 0.5c Shallow depth (0-100 km)
S 0.25c c 0.25c green 0.25p 0.5c Intermediate depth (100-300 km)
S 0.25c c 0.25c blue 0.25p 0.5c Very deep (> 300 km)
D 0 1p
V 0 1p
N 7
V 0 1p
S 0.25c c 0.15c - 0.25p 0.75c M 3
S 0.25c c 0.20c - 0.25p 0.75c M 4
S 0.25c c 0.25c - 0.25p 0.75c M 5
S 0.25c c 0.30c - 0.25p 0.75c M 6
S 0.25c c 0.35c - 0.25p 0.75c M 7
S 0.25c c 0.40c - 0.25p 0.75c M 8
S 0.25c c 0.45c - 0.25p 0.75c M 9
D 0 1p
V 0 1p
N 1
END
# Put together a reasonable legend text, and add logo and user's name:
cat <<- END >> neis.legend
G 0.25l
P
T USGS/NEIS most recent earthquakes for the last month. The data were
T obtained automatically from the USGS Earthquake Hazards Program page at
T @_https://earthquake.usgs.gov@_. Interested users may also receive email alerts
T from the USGS.
T This script could be called monthly to update the latest information.
G 0.9c
# Add USGS logo
I @USGS.png 1i RT
G -0.75c
L 12p,Times-Italic LB $me
END
# OK, now we can actually run gmt legend. We center the legend below the map.
# Trial and error shows that 1.7i is a good legend height:
gmt legend -DJBC+o0/1c+w18c/4.2c -F+p+glightyellow neis.legend --FONT_ANNOT_PRIMARY=10p,Helvetica
rm neis.legend usgs_quakes_22.txt
gmt end show