#!/bin/tcsh -f
#
# plot a Cassini projected map of Alaska with a user-defined scale
#  This is a bit different from the lecture:
#     The map scale is hard-coded instead of passed as an argument
#     We use the projection code "C" instead of "c", which changes how
#       the scale is interpreted. With "C", the scale is computed by GMT
#       so that the final map size is given by the scale parameter
#       (-JC$pro_lon/$pro_lat/6i = 6 inches wide)

set workingdir = `pwd`
gmtset `cat $workingdir/gmt.plot.parameters`

set scale = 6i
#   The lecture example used $1 to allow the user to set scale at runtime

set w_lon = 180
set e_lon = 225
set s_lat = 50
set n_lat = 72
set map_box = -R$w_lon/$e_lon/$s_lat/$n_lat

set filename = alaska.ps

set pro_lon = `echo $w_lon $e_lon | awk '{printf"%.2f",($1+$2)/2.0}'`
set pro_lat = `echo $n_lat $s_lat | awk '{printf"%.2f",($1+$2)/2.0}'`
set map_pro = -JC$pro_lon/$pro_lat/$scale

gmt pscoast $map_pro $map_box -Di -N1 -G125/255/125 -K -P > $filename

gmt psxy $map_pro $map_box -Sc0.25 -G255/125/125 -K -O >> $filename << END
210 70
210 65
210 60
210 55
END

gmt psbasemap $map_pro $map_box -Bg5a10 -O >> $filename
