A script full of GMT commands

Sorry, this one ain't standalone so you gotta be at the seminar for the explanation. It makes an image that looks like this.


#!/bin/csh

# CERI_NM 
#
# plot dem raster
# plot NMSZ epicenters
# plot CERI NMSZ stations
# plot telemetry paths

#gmtdefaults I changed (man gmtdefaults for explanation)
# PAPER_MEDIA = a4+
# MEASURE_UNIT = inch
# DEGREE_FORMAT = 5

# map corners
set LONMIN = -90.7
set LONMAX = -89.0
set LATMIN = 35.0
set LATMAX = 37.0

#inches per degree
set SCALE = 3.5
 
# used to center on page
set XOFFSET = 1.50
set YOFFSET = 1.0

# postscipt file goes here
set OUTPUTFILE = ../ps/$0.ps

#illum direction determines where the "sun" is
set ILLUM = 320

# params for color file CPTMASTER is the color table and CSAMP is min/max/increment
set CPTMASTER = rainbow
set CSAMP = -S-50/180/1.0

# other more unsightly colormaps
#set CPTMASTER = sealand
#set CSAMP = -S-100/150/0.1
#set CPTMASTER = topo
#set CSAMP = -S-10/200/0.1
#set CPTMASTER = globe
#set CSAMP = -S-200/300/0.1

# limits when make colors
#set CLIM = -L1/300
set CLIM =

# use -P for portait, else landscape
set ORIENT = -P

# handy shapes
set HEXAGON = -Sh0.08
set TRIANGLE = -St0.2
set SQUARE = -Ss0.2
set CIRCLE = -Sc0.2

# handy fonts
set STA_SIZE = 14
set STA_ANG = 0
set STA_FONT = 5
set STA_JUST = 11

#resolution of gmt dataset
set RESOL = -Df

#standard switched for more to follow
set CONTINUE = "-K -O"

# type of map projection
set PROJ = Jm

# where to put scratch files
set SCRATCH = ../scratch

# default colors
set WHITE = 255
set LTGRAY = 192
set VLTGRAY = 225
set EXTGRAY = 250
set GRAY = 128
set BLACK = 0
set RED = 250/0/0
set DKRED = 196/50/50
set BLUE = 0/0/255
set GREEN = 0/255/0
set YELLOW = 255/255/50
set ORANGE = 255/192/50
set PURPLE = 255/50/255
set CYAN = 50/255/255
set LTBLUE = 192/192/250
set VLTBLUE = 225/250/250
set LTRED = 250/225/225
set PINK = 255/225/255
set BROWN = 160/64/32
set LTBROWN = 200/150/150

# convert corners to a region gmt understands
set REGION = $LONMIN/$LONMAX/$LATMIN/$LATMAX

# datasets defined in the grdraster.info file which should be
# found in the directory pointed to by $GMT_GRIDDIR
# as of this writing 9 is US 30 second USGS DEM
set DATASET = 9

# defines the resolution of the grid (30 by 30 seconds)
set DATAGRID0 = -I30c/30c

# now make the grd file
echo performing grdraster of dataset $DATASET
grdraster $DATASET -G$SCRATCH/$0_topo0.grd $DATAGRID0 -R$REGION

# resample
echo resampling grid
#set DATAGRID = -I10c/10c
set DATAGRID = -I30c/30c
grdsample $SCRATCH/$0_topo0.grd -G$SCRATCH/$0_topo.grd $DATAGRID -F -R$REGION

#clean up the increments for GMT 3.1 and older files
echo fixing grd file
grdedit -A $SCRATCH/$0_topo.grd


# what norm to use when creating intensities
set NORM = -Nt
#set NORM = -Ne

# get the intensities
echo working on illumination for $CPTMASTER pallette.
grdgradient $SCRATCH/$0_topo.grd -A$ILLUM -G$SCRATCH/$0_topo.intns $NORM

# now make the color palette
grd2cpt $SCRATCH/$0_topo.grd -C$CPTMASTER $CLIM $CSAMP  >! $SCRATCH/$0.cpt

# grid the dem (use -M for monochrome)
echo "working on dem"
grdimage $SCRATCH/$0_topo.grd -I$SCRATCH/$0_topo.intns -C$SCRATCH/$0.cpt -R$REGION -$PROJ$SCALE -Ba1:." ":/a1neSW -K -X$XOFFSET -Y$YOFFSET $ORIENT >! $OUTPUTFILE

#can we start drawing maps now please?
pscoast -R$REGION$SCALE -$PROJ -C$BLUE -S$BLUE -I1/10/$BLUE -I2/10/$BLUE -I3/10/$BLUE -I4/10/$BLUE -N2/10/$BLACK -N1/10/$BLACK -Df -W5 -O -K $ORIENT >> $OUTPUTFILE

#plot epicenters
set MYDIR = ../data
echo less than 2
nawk '{print ( $10 <= 2.0 ? $7 " " $8 " " $10/50.0 : ">" ) }' $MYDIR/NewMad2k.cat | \
  psxy -R -: -$PROJ -Sc -K -O -L -G$LTBROWN -W1/$LTGRAY -V >> $OUTPUTFILE

echo 2 to 3
nawk '{print ( $10 <= 3.0 && $10 > 2.0 ? $7 " " $8 " " $10/50.0 : ">" ) }' \
  $MYDIR/NewMad2k.cat | \
  psxy -R -: -$PROJ -Sc -K -O -L -G$LTBROWN -W1/$LTGRAY >> $OUTPUTFILE

echo 3 to 4
nawk '{print ( $10 <= 4.0 && $10 > 3.0 ? $7 " " $8 " " $10/50.0 : ">" ) }' \
  $MYDIR/NewMad2k.cat | \
  psxy -R -: -$PROJ -Sc -K -O -G$LTBROWN -L -W2/$LTGRAY >> $OUTPUTFILE

echo 4 to 5
nawk '{print ( $10 <= 5.0 && $10 > 4.0 ? $7 " " $8 " " $10/50.0 : ">" ) }' \
  $MYDIR/NewMad2k.cat | \
  psxy -R -: -$PROJ -Sc -K -O -G$BROWN -W3/$GRAY -L >> $OUTPUTFILE

echo greater than 5
nawk '{print ( $10 > 5.0 ? $7 " " $8 " " $10/50.0 : ">" ) }' $MYDIR/NewMad2k.cat | \
  psxy -R -: -$PROJ -Sc -K -O -G$BLACK -W3/$BLACK -L >> $OUTPUTFILE

# central rx
echo adding central receives
psxy $MYDIR/cnmsn-cr-op.dat -R -Jm -: $HEXAGON -K -O -G$PURPLE -W3/$BLACK -L >> $OUTPUTFILE
 
# add stns
echo adding stations
psxy $MYDIR/bb.dat -R -Jm $SQUARE -K -O -: -G$ORANGE -W3/$BLACK -L >> $OUTPUTFILE
psxy $MYDIR/lnxt.node.dat -R -Jm $TRIANGLE -K -O -: -G$BLACK -W3/$BLACK -L >> $OUTPUTFILE
psxy $MYDIR/nmad.node.dat -R -Jm $TRIANGLE -K -O -: -G$BLACK -W3/$BLACK -L >> $OUTPUTFILE
psxy $MYDIR/mkta.node.dat -R -Jm $TRIANGLE -K -O -: -G$BLACK -W3/$BLACK -L >> $OUTPUTFILE
psxy $MYDIR/ceri.node.dat -R -Jm $TRIANGLE -K -O -: -G$RED -W3/$BLACK -L >> $OUTPUTFILE
psxy $MYDIR/slu.dat -R -Jm $TRIANGLE -K -O -: -G$BLUE -W3/$BLACK -L >> $OUTPUTFILE

#add station names
echo adding station names

foreach namefile (ceri.node nmad.node lnxt.node mkta.node bb)
  nawk '{print $2 " " $1 " 12 0 5 1 " $3}' $MYDIR/$namefile.names | \
  pstext -R -Jm -G$BLACK -O -K >> $OUTPUTFILE
end
  
#pstext -R -Jm -G$BLACK -O -K <> $OUTPUTFILE
#  `nawk '{print $2 " " $1 " 12 0 5 1 " $3}' $MYDIR/ceri.node.names`
#  `nawk '{print $2 " " $1 " 12 0 5 1 " $3}' $MYDIR/nmad.node.names`
#  `nawk '{print $2 " " $1 " 12 0 5 1 " $3}' $MYDIR/lnxt.node.names`
#  `nawk '{print $2 " " $1 " 12 0 5 1 " $3}' $MYDIR/mkta.node.names`
#  `nawk '{print $2 " " $1 " 12 0 5 1 " $3}' $MYDIR/bb.names`
#END

#add telemetry paths
echo adding telemetry paths
psxy $MYDIR/nmad.backbone.dat -: -R -$PROJ -M -W8/$BLACK -K -O >> $OUTPUTFILE
psxy $MYDIR/ceri.backbone.dat -: -R -$PROJ -M -W8/$BLACK -K -O >> $OUTPUTFILE
psxy $MYDIR/lnxt.backbone.dat -: -R -$PROJ -M -W8/$BLACK -K -O >> $OUTPUTFILE
psxy $MYDIR/mkta.backbone.dat -: -R -$PROJ -M -W8/$BLACK -K -O >> $OUTPUTFILE
psxy $MYDIR/bb.backbone.dat -: -R -$PROJ -M -W12/$BLACK -K -O >> $OUTPUTFILE
psxy $MYDIR/digi.backbone.dat -: -R -$PROJ -M -W12/$RED -K -O >> $OUTPUTFILE

# memphis
echo plotting Memphis
echo " -89.941 35.122 " | \
  psxy -R -Jm -Sc0.16 -K -O -W2/$BLUE -G$BLACK -L >> $OUTPUTFILE
echo " -89.941 35.122 " | \
  psxy -R -Jm -Sc0.08 -K -O -W1/$BLACK -G$BLUE -L >> $OUTPUTFILE
echo " -89.85 35.1 14 0 5 5 Memphis" | \
  pstext -R -Jm -G$BLUE -K -O >> $OUTPUTFILE

#clean up
\rm $SCRATCH/$0*

# make a title
echo " -89.85 37.1 24 0 29 2 CERI NMSZ Seismic Network " | \
  pstext -R -N -$PROJ -G$BLACK $ORIENT -O >> $OUTPUTFILE