The earthquake map is a Lambert Conformal Conical projection, with standard
parallels of 34.5N and 36.5N. It is centered at 35.5N and 98.62W. Either
GMT or the Lambert Conformal Conical projection (or both) require projection
from a spheroid rather than an ellipsoid.

The map is plotted with Generic Mapping Tools by Paul Wessel and Walter
H. S. Smith. The county lines are plotted with an Oklahoma County file
made by Ken Luza of the OGS. A Perl program read the catalog html files
(for example www.okgeosurvey1.gov/level2/okeqcat/okeqcat.2000.html) and
produced GMT x-y files giving the location, color (red or black), and 
diameter of the circle. The postscript output files were imported, cropped,
and animated by ImageMagick.

The following PERL script converts Oklahoma Earthquake yearly catalog pages
in html to yearly files for input to GMT.


#! /usr/local/bin/perl -w
#mod 20010415 to allow "R" after latitude and/or longitude
#mod and corrections 2000 OCT09
#ok2gmt.perl
#mod from ok2vg.perl to ok2gmt.perl 2000JUL21 Jim L.
#to write files for GMT psxy
#uses /wf2/local/share/GMT_hot.cpt
#0=black(0,0,0), red=(255,0,0)
#circle diameter is in cm.

#original format version

#ok2geotoolmapoverlay.nawk modified 1997 JUN27 from ok2llnl.awk by Jim
#ok2llnl.awk reads in okeqcat format and produces LLNL map format
#w/Jim Lawson, Oklahoma Geological Survey Observatory, Leonard, 1995 JUL19.
#typical command line: awk -f ok2llnl.awk okeqmap.1995 > ok1995.eve
#mod 1999 MAR16 to only write latitude and longitude
#mod to make circle size variable date unknown
#mod 1999 APR26 to use command line infile and outfile names and option base 0
#mod 1999 APR27 to ok2vgeo.perl. REGEX used to accept lines and capture
#parameters.
#mod 1999 APR29 to make circles of felt events red and to write trailer
#comment information. Renamed ok2vg.perl.
#mod 1999 JUL20 for new format catalog.
#mod REGEX to allow '-' longitude
#mod REGEX for no required space after longitude
#mod REGEX for one (was 2) spaces after mm
#mod to give two notices that OUTFILE is being appended to,
#and was appended to if it existed before this program started

if($#ARGV != 1) {print "USAGE ok2gmt.perl infile[.html] outfile[.gmt.dat]\n";
                 exit}
print STDOUT $#ARGV, "  ",$ARGV[0], "  ", $ARGV[1], "\n";

open INFILE, '<' . $ARGV[0];
$outfile_exists = 0;
if(-e $ARGV[1]){
               print "\007TAKE NOTICE: " , $ARGV[1],
               " exists, ok2gmt.perl is APPENDING to it.\n";
               $outfile_exists = 1;
               }
open OUTFILE, '>>' . $ARGV[1];

#print OUTFILE "sources\n"; #required first line of map file
$source = 0;
$felt = 0;
while () {
    chop;       # strip record separator
print STDOUT $_, "\n";
if($_ =~ m!(?xo)
    (?#year)     ^\>?[12][890][0-9][0-9]\040
    (?#month )   (?:JAN|FEB|MAR|APR|MAY|JUN|JUL|AUG|SEP|OCT|NOV|DEC)\040
    (?#day)      \d\d.{28}
    (?mm)        ([1-9F]|\040)\040
    (?#m3hz)     (\d\.\d|\040{3})\040
    (?#mblg)     (\d\.\d|\040{3})\040
    (?#mdur)     (\d\.\d|\040{3})\040
    (?#lat)      (3[1-9]\.[0-9]{2,4})R?\040+
    (?#lon)      (-?(?:10|9)[0-9]\.[0-9]{2,4})R?!) {
$mm   = $1;
$m3hz = $2;
$mblg = $3;
$mdur = $4;
$lat  = $5;
$lon  = $6; #earlier versions had $lon = -$6 as old format
            #used pos numbers for west longitude

#comment out next two lines except for spreading older map
#epicenters, esp where there are many with same approx location
#$lon = $lon + (rand(0.245)) - 0.1225;
#$lat = $lat + (rand(0.3)) - 0.15;

$symbol_color = 0;  #not known to be felt
if($mm =~ m/[1-9F]/o) {$symbol_color = 0.375;  #known to be felt
    $felt += 1}

$source += 1;
#$mag=9.9; #This serves as an alarm for an event without a magnitude
if ($mdur ne "   "){$mag=$mdur;} #notice priority of magnitudes
if ($m3hz ne "   "){$mag=$m3hz;} #$mblg, then $m3hz, then $mdur
if ($mblg ne "   "){$mag=$mblg;} #The magnitudes are NOT averaged

$m2 = ((1 + $mag) * (1 + $mag))*0.025; #0.025 cm (.25mm) in GMT:psxy Sc Option
$m2 = 1.25 * $m2;  #mod 20010415 as some circles almost invisible on map


printf STDOUT "%8.5f %10.5f %5.3f %6.2f\n",
       $lat, $lon, $symbol_color, $m2;
printf OUTFILE "%8.5f %10.5f %5.3f %6.2f\n",
       $lat, $lon, $symbol_color, $m2;

} #end if regex block
} #end read block

# '!' is required at the beginning of a geotool map file comment
print STDOUT "! filein= ", $ARGV[0], "   fileout=", $ARGV[1], "\n";
print STDOUT "! number of events= ", $source, " of which ", $felt,
" were felt\n\n";

if($outfile_exists == 1){
               print "\007TAKE NOTICE: " , $ARGV[1],
               " existed, ok2gmt.perl APPENDED to it.\n"
                        }

close INFILE;
close OUTFILE;
exit



The following GMT script plots each year's map.

#!/usr/bin/csh
set YEAR = $1
#mod to ok11.cmd Jim Lawson 2000 OCT12
#ok9.cmd GMT map command. Jim Lawson 2000 JUL23
# 
#map 1 Oklahoma#################################################################
#######
#Lambert conical conformal parallels 34.5 36.5
#center 35.5 -98.62
#edges 33.5 37.25 -94.0 -103.25
#NB -M\# is a cshell escape. Without "\" csh sees # as end of command
#
#plot Oklahoma counties file derived from Ken Luza's file
#
#*******USE inverseletter**********

psxy -JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 \
-Ba1f.5/a1f.5WSen  \
-M\# -: -X+1.7 -Y+7 ok.county.dat -K >! ok11.ps
#
#plot 50 km scale
#
psbasemap -JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 \
 -Lf-100/33.9/slat35.5/50k -O -K >> ok11.ps
#
#plot all rivers full resolution file
#
pscoast -JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 \
-Df -Ia/W2/0/0/255ta -O -K >> ok11.ps
#
#plot earthquakes
#
psxy -JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 \
-Cmy_hot.cpt -Sc -: $YEAR.gmt.dat -K -O >> ok11.ps
#
#plot OKLAHOMA GEOLOGICAL SURVEY
#
echo "-103 34.5 24 0 9 LM OKLAHOMA"|pstext \
-JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 -K -O >> ok11.ps
echo "-103 34.25 24 0 9 LM GEOLOGICAL"|pstext \
-JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 -K -O >> ok11.ps
echo "-103 34 24 0 9 LM SURVEY"|pstext \
-JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 -K -O >> ok11.ps
#
#plot JL
#
echo "-94.4 33.7 20 0 9 LM JL"|pstext \
-JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 -K -O >> ok11.ps
#
#plot $YEAR, which is a parameter from calling ok.a.csh
#
echo "-102 36 34 0 9 LM $YEAR"|pstext \
-JL-98.62/35.5/34.5/36.5/22.8c -R-103.25/-94/33.5/37.25 -O >> ok11.ps

gs ok11.ps

#The output map ok11.ps is displayed in postscript. It is then imported
#into gif by ImageMagick. When all the years are made, ImageMagick
#produces an animated gif.

OGS HOME