1 / 38

Generic Mapping Tools (GMT)

Generic Mapping Tools (GMT). Maps. Final script for Example 1. #!/bin/bash psbasemap -R0/100/0/10 -JX5i/3i -Ba20g10/a2f1g2WSne - X2 -Y14 -K -P > plot.ps sample1d -I1 << END | nawk '{print $1, sqrt($1)}' > input.dat 0 100 END psxy input.dat -R -JX -W5t15_15:0 -K -O >> plot.ps

avi
Download Presentation

Generic Mapping Tools (GMT)

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Generic Mapping Tools (GMT) Maps

  2. Final script for Example 1 #!/bin/bash psbasemap -R0/100/0/10 -JX5i/3i -Ba20g10/a2f1g2WSne \ -X2 -Y14 -K -P > plot.ps sample1d -I1 << END |\ nawk'{print $1, sqrt($1)}' > input.dat 0 100 END psxyinput.dat -R -JX -W5t15_15:0 -K -O >> plot.ps sample1d input.dat -I10 |\ psxy -R -JX -St0.5 -G255/0/0 -W5/0/255/0 -O >> plot.ps gsplot.ps

  3. Correction on the –W flag psxyinput.dat –R –JX -W5t15_15:0 -K –O >> plot.ps The –W command can be quite complicated and is specified in order –W[width],[color],[texture] In this example, the first 5 represents line width, the t signifies a texture follows, the 15_15 specifies the dash and space width, and the :0 specifies that a dash is used. Alternately, -W5,0,15_15:- does the same thing.

  4. Common command options on first, and possibly subsequent, calls • Need on all calls • -R Define region for plot – will need on first call and at least “–R” on subsequent • -J Define projection for plot – will need this on all calls if need to define region

  5. Common command options on first, and possibly subsequent, calls • (Generally) Need on first call only • -B Borders -- annotation, frame, grid. Only need on first (or a single) call. • -P Switch between landscape and portrait modes • -X Shift X axis • -Y Shift Y axis

  6. Common command options on first, and possibly subsequent, calls. • Need when needed. • -K Don’t close PostScript (showpage), use when more will follow • need on all but last GMT call • -O Don’t initialize PostScript, use when appending to pre-existing file • need on all but first GMT call • use both –K and –O when putting a large number of GMT call outputs together

  7. Common command options on first, and possibly subsequent, calls. • Need when needed. • -V Verbose (prints out stuff to standard error for user). • -H Header records (tells GMT to skip first H lines of ascii input file)

  8. GMT Defaults • There are about 100 parameters which can be adjusted individually to modify the appearance of plots or affect the manipulation of data. Each as a default value. • GMT defaults are kept in a file called ~/.gmtdefaults4. There are tons of them and you can find out what they are and what the mean reading the man page for gmtdefaults.

  9. When a program is run, it initializes all parameters to the GMT defaults, then tries to open the file .gmtdefaults4 in the current directory. If not found, it will look for that file in a sub-directory ~/.gmt, and finally in your home directory itself. If successful, the program will read the contents and set the default values to those provided in the file. • To view your current gmtdefault settings %gmtdefaults –L • To view the list of options for each default parameter %man gmtdefaults

  10. Plotting Defaults example of start of .gmtdefaults4 # GMT-SYSTEM 4.2.1 Defaults file #-------- Plot Media Parameters -- PAGE_COLOR = 255/255/255 PAGE_ORIENTATION = landscape PAPER_MEDIA = letter #--- Basemap Annotation Parameters -- ANNOT_MIN_ANGLE = 20 ANNOT_MIN_SPACING = 0 ANNOT_FONT_PRIMARY = Helvetica ANNOT_FONT_SIZE = 14p ANNOT_OFFSET_PRIMARY = 0.075i

  11. Changing the defaults • You can edit your local copy of .gmtdefaults4 using nedit or vim • You can explicitly reset a default within a script using the command gmtset #!/bin/sh gmtset PAPER_MEDIA letter MEASURE_UNIT cm gmtset OUTPUT_DEGREE_FORMAT +D gmtset PLOT_DEGREE_FORMAT +D

  12. How about making prettyMAPS?

  13. Map projections available in GMT

  14. List of “standard” command line options. The –J option sets the “projection” One has to look at the man page for each one as “different things vary”

  15. pscoast -R-90/-70/0/20 -JM6i -P -B5g5 -G180/120/60 > map1.ps “All” gmt programs plot “maps” through the projection command line option or switch (even the x-y plot).

  16. pscoast -R-90/-70/0/20 -JM6i -P -B5g5 -G180/120/60 > map1.ps All projections give you two selections for specifying the scale (note GMT takes the mapmakers attitude that a map has to have a predetermined/known scale – nicely filling the page does not cut it.)

  17. pscoast -R-90/-70/0/20 -JM6i -P -B5g5 -G180/120/60 > map1.ps -Jmparameters (Mercator [C]). Specify one of: -Jmscaleor -JMwidth Give scale along equator (1:xxxx or UNIT/degree).

  18. pscoast -R-90/-70/0/20 -JM6i -P -B5g5 -G180/120/60 > map1.ps -Jmlon0/lat0/scaleor -JMlon0/lat0/width Give central meridian, standard latitude and scale along parallel (1:xxxx or UNIT/degree, UNIT = number inches or cms).

  19. Mercator Projection: One way to address plotting sphere on a plane (which is whole ‘nother subject) Conformal (maintains shapes) Cylindrical projection

  20. pscoast -R-130/-70/24/52 -JB-100/35/33/45/6i -B10g5:."Conic\ Projection": -N1/2p -N2/0.25p -A500 -G200 -W0.25p -P >! map.ps Region is “rectangle” on the spherical earth. -N for boundaries (international, US/Canadian/Mexican state boundaries “built in”), rivers.

  21. pscoast -R-130/-70/24/52 -JB-100/35/33/45/6i -B10g5:."Conic\ Projection": -N1/2p -N2/0.25p -A500 -G200 -W0.25p -P >! map.ps -A to get rid of small water/island features Albers projection (b/B) – need to know something (center and standard parallels).

  22. pscoast -R-130/-70/24/52 -JB-100/35/33/45/6i -B10g5:."Conic\ Projection”: -N1/2p -N2/0.25p -A500 -G200 -W0.25p -P >! map.ps -Jblon0/lat0/lat1/lat2/scale or -JBlon0/lat0/lat1/lat2/width Albers - Give projection center, two standard parallels, and scale (1:xxxx or UNIT/degree).

  23. Albers Also conformal (maintains/conserves shape) Conical projection

  24. pscoast -R0/360/-90/90 -JG280/30/6i -Bg30/g15 -Dc -A5000 \ -G255/255/255 -S150/50/150 -P >! map.ps azimuthal Orthographic projection (g/G) mimics looking at earth from infinite distance.

  25. pscoast -R0/360/-90/90 -JG280/30/6i -Bg30/g15 -Dc -A5000 \ -G255/255/255 -S150/50/150 -P >! map.ps New option -Dc Controls resolution of coastline f full h high l low c crude Helps manage file sizes.

  26. Example 2Simple location mapwith focalmechanisms

  27. gmtset: To change individual GMT default parameters • pscoast: Plot coastlines, filled continents, rivers, and political borders • psxy: Plot symbols, polygons, and lines in 2-D • pstext: Plot text strings • psmeca: Plot focal mechanisms on maps

  28. #!/bin/sh ###Create a basemap for the Sumatra-Andaman Island region #Set global variables Plots=../BasicInfo CPT=/usr/local/GMT4.2.0/share/cpt LATMIN=-6.2 LATMAX=16 LONMIN=90 LONMAX=106 OUTFILE=$0.ps#What does the $0 signify on this line? #Explicitly set some GMT default values gmtsetMEASURE_UNIT cm ANOT_FONT_SIZE 12 LABEL_FONT_SIZE 12 gmtsetBASEMAP_TYPE fancy

  29. #Map • # Create the basic basemap using GMT coastline data • pscoast-JM15 -R$LONMIN/$LONMAX/$LATMIN/$LATMAX -Ba2f1WNes -Y5 -X3 -K -P \ • -Dh -A100 -N1 -W1 -G155 > $OUTFILE • -D Selects the resolution of the data set ((f)ull, (h)igh, (i)ntermediate, (l)ow, (c)rude) • -A Features with an area smaller than min_area in km^2 • -W and –G Control the fill color and line color of the land regions • -N Plots political boundaries. 1 is for national boundaries, 2 is for states in the Americas

  30. #Plot basic tectonic information, including subduction zones and volcanoes • psxy$Plots/trench.right.gmt -R -JM -M -W5/0 -Sf0.3i/0.08irt -G0 -O -K >> $OUTFILE • psxy$Plots/trench.left.gmt -R -JM -M -W5/0 -Sf0.3i/0.08ilt -G0 -O -K >> $OUTFILE • psxy$Plots/trench.other.gmt -R -JM -M -Wt30_30:1p5/0 -O -K >> $OUTFILE • psxy$Plots/transform.gmt -R -JM -M -Wt10_10:1p5/0 -O -K >> $OUTFILE • psxy$Plots/ridge.gmt -R -JM -M -W1/0 -O -K >> $OUTFILE • psxy$Plots/volcanoes.simkin+siebert.gmt -R -JM -St.35 -W1/0 -G0 -O -K >> $OUTFILE • -M Multiple segment file. Segments are separated by a record whose first character is flag. [Default is '>']. Example of trench.right.gmt • > • 92.278 6.948 • 91.903 7.663 • >

  31. #Information for the 2004 Sumatra event • awk '{print $7, $6, $9/30}' $Plots/122604_032705.pde |\ #prints lat, lon,mag/30 • psxy-R -JM -Sc -W1/0 -G255/255/100 -O -K >> $OUTFILE • When no value is specified on the –S flag (ie 0.20), then the size of the circle (denoted by the c) is controlled by the third column of info sent to the psxy command. • So we are scaling the circle size by what earthquake parameter in this example? • Example of 122604_032705.pde from http://neic.usgs.gov/neis/epic/epic_rect.html • PDE-W 2004 12 26 005853.45 3.30 95.98 30 9.00 MwHRV 9C M STS...M • PDE-W 2004 12 26 011710.33 4.94 94.27 30 5.50 mb GS .. . .......

  32. psxy-R -JM -Sa.89 -W1/0 -G255/255/0 -: -O -K << END >> $OUTFILE • 3.3 95.98 • END • -Sa creates a star • -: indicates that the data is in y-x format rather than x-y format. This is useful because we speak of earthquakes in terms of latitude (y) and longitude (x) and therefore tend to write geographic coordinates in this order as well.

  33. #Global CMT solutions for events above Mw ~6.6 • psmeca -R -JM -Sm0.5 -G255/255/100 -T -O -K << END >> $OUTFILE • #lon lat depth mrrmttmppmrtmrpmtpiexp name • 94.26 3.09 29 1.04 -0.43 -0.61 2.98 -2.40 0.43 29 X Y 122604A • 92.79 6.61 14 5.26 -0.84 -4.41 3.95 -2.91 2.10 26 X Y 122604B • 92.45 8.58 12 0.94 -0.12 -0.81 0.02 -0.27 0.32 26 X Y 122604C • 92.20 4.99 12 -0.02 -0.64 0.66 -0.11 -0.16 -0.80 26 X Y 010105A • 95.38 2.84 12 0.26 -0.17 -0.09 0.90 -0.91 0.16 26 X Y 022605A • END • psmeca creates lower hemisphere projections, aka beachballs. It is a rather complicated little command. • -S Selects the meaning of the columns in the data file. • -Sa Focal mechanisms in Aki and Richard convention. • -Sc Focal mechanisms in Harvard CMT convention. • -Sm Seismic moment tensor (Harvard CMT, with zero trace). • -Sd plot only the double couple part of moment tensor. • -Sz plot anisotropic part of moment tensor (zero trace). • http://www.globalcmt.org/CMTsearch.html • -Sp Focal mechanisms given with partial data on both planes. • -Sx Principal axis. • -G Controls fill of the compressive quadrant of the moment tensor • -T Plots the nodal planes

  34. pstext-R -JM -: -O << END >> $OUTFILE • -3.5 90.5 14 0 0 0 Yellow: 2004 Mw 9.2 series • -5.75 90.5 14 0 0 0 CMT solutions for Mw >= 6.6 • END • format is x, y, size, angle, font#, justify, text • %pstext –L will tell you the font # and font name • you can use [T|M|B][L|C|R] to specify justification • (top/middle/bottom/left/center/right) • gs$OUTFILE #gs indicates ghostscript if you are using vim to edit this script :!script.gmtwill run the script, open the resulting file using gs, and you can see what you have done. If you don’t like it, type quit, and continue editing the file.

  35. nawk commands as input to GMT nawk '{print $9, $8, $11}' EBH.HDF You can put this into GMT several ways If this is the only file you want to plot – this would work nawk '{print $9, $8, $11}' EBH.HDF | pxsy …

  36. If you had a number of files that needed conversion you could do it this way (only need one psxycall) psxy … << END … … `nawk '{print $9, $8, $11}' EBH.HDF` … END Converting each file on the fly.

  37. If you want to do the same thing to a list of files filelist=“$SAMDATA/eq-rupt-1995.dat $DEM/eq-rupt-1960.dat” for FILE in $filelist do psxy -R -$PROJ$SCALE -M$ -: $CONTINUE -W$LINETHICK/$PURPLE $FILE \ $VBSE >> $OUTPUTFILE done Other ways to make list (notice the different kinds of quotes: “,’ and `) filelist=`ls -1 $ROOT/dem/topocontours/andes_3000_*` contourlist=‘1 2 3 4’

  38. Some other nawk tricks – doing math and passing variables to nawk(quote heaven) SCALE=`echo $STNDTMLON | \ nawk ‘{print ($1>=0?$1:360+$1)”/”’${jTRESCALE}_1’*’$FACTOR’}’ `

More Related