Continuing from the previous post about GMT, now we will cover some others code to produce a map, and also to get some data from the free available sources. To start it, we may try a very simple command to produce a map. Keep in mind that, in this post we work on a Linux Operating System (OS) with GMT version 4.5 (some features, i.e., to make a file executable, post script viewing (.ps) and others, you may need to change it or install some program in order to work and viewing the output).
Here is the idea: We want to produce a map that consist of considerably regional map and a study area within a single plot (one figure with two maps ~ or we can say it, as such as the famous “figure of study area”).
Well, there are couple of ways to do it so, and we will start again as we mention earlier as a very simple one. Here we go:
A. We will only use a pscoast command to make it. You may try it as follows:
#!/bin/csh
#environment setting
gmtset PAPER_MEDIA A3
gmtset LABEL_FONT_SIZE = 12p
gmtset BASEMAP_TYPE plain
gmtset ANNOT_FONT_SIZE_PRIMARY = 12p
gmtset ANNOT_FONT_SIZE_SECONDARY = 12p
gmtset HEADER_FONT_SIZE = 14p
###regional map
pscoast -R95/102/-2/6 -JM4i -Ba1f0.5NSEw -Dh -W0.5p,black -X11.7i -Y6.8i -K -V > nearshore.ps
###study area map
pscoast -R95.18/95.55/5.5/5.75 -JM10i -Ba0.05f0.01NSEW -Dh -W0.5p,black -X-11i -Y-6.5i -O -V >> nearshore.ps
gv nearshore.ps
After you write this code you may make it executable by:
After the program run, you may get a figure as shown below.
Now you get a figure as what we have in mind. Okay, now, I will try to explain the code that I wrote above so it will be meaningful for everyone. The explanation are as follows:
- The first line (LN) is a command for Linux to access the csh.
- LN 3-8 are some setting that I want to change for my gmt default (this is kind of personal stuff, so you may do it as what ever you want). You may check other options, just type in your terminal gmtdefaults -L to make a change or keep it, as that is based on your personal taste.
- LN 11 is the pscoast command that will plot you the regional map. the -R option gives you the range coverage of your map, -JM is your map projection (Mercator) with the size of 4i (i stand for inch), -Ba is your map annotation which shown the degree value (we set it to show the value in each 1 degree and in between we set a 0.5 degree as a small thick), and the NSEw will show you where you need to show the degree value (capital later) and no need to show it (small later). -D indicate the coast line resolution that you gonna use and h is high. Caution : you may check it first in your machine, is the high resolution coast line is available or not. You may check it by going to the root as is shown in figures below.
And then you may find
If you see something just like the one indicated by the red rectangle, then you already have the high resolution coast file installed in your machine. If it is not showing like that, then you only have the low resolution (please use -Dl). To install it, you may just yum install GMT-coastlines-high.noarch. Okay, we continue, -W will gives you the thickness of the coast line and we use 0.5p (p is point/pixel) and a comma sign (,) to gives the color (here we used a black color).-X and -Y were used to set the location of the map on the paper with (+) X and Y value will put your map to the right and up respectively, and opposite with (-) value. Here we want to put the regional map at the upper righthand side of the study area map. -K indicate that you will have another command after this line and -V is a verbose mode (report everything in terminal during the execution/ the program run). The greater then sign (>) will gives you the output and the output file name is nearshore.ps.
- LN 14 is the same thing as what we already explain in number 3, but here we make a smaller range (see the -R), make a larger size (see -JM), set up the location of the map on the paper (see -X and -Y), an additional of -O is to overlay this map with the previous one and two greater then symbols (>>) to print this map in the same output file as previous command.
- LN 16 is the ghost view command that will pop up the output file as soon as it is finish (to make your life easy).
Okay, we already got the figure, are we finish?not quiet. We can add the location of the study area in our regional map. This can be done by adding a psxy command after the first pscoast command. Now the GMT script is shown like this:
#!/bin/csh
#environment setting
gmtset PAPER_MEDIA A3
gmtset LABEL_FONT_SIZE = 12p
gmtset BASEMAP_TYPE plain
gmtset ANNOT_FONT_SIZE_PRIMARY = 12p
gmtset ANNOT_FONT_SIZE_SECONDARY = 12p
gmtset HEADER_FONT_SIZE = 14p
###regional map
pscoast -R95/102/-2/6 -JM4i -Ba1f0.5NSEw -Dh -W0.5p,black -X11.7i -Y6.8i -K -V > nearshore.ps
###make a rectangle indicate for the study area
psxy -R -J -W0.01i/red -K -O -V << EOF >> nearshore.ps
95.18 5.5
95.55 5.5
95.55 5.75
#> next (this gonna disconnect the points)
95.18 5.75
95.18 5.5
EOF
###study area map
pscoast -R95.18/95.55/5.5/5.75 -JM10i -Ba0.05f0.01NSEW -Dh -W0.5p,black -X-11i -Y-6.5i -O -V >> nearshore.ps
gv nearshore.ps
You may execute this program again and it will give you this figure:
If you see carefully there is a red rectangle at the upper left corner of the regional map (northern tip of Sumatra Island).
To make that rectangle see the ###make rectangle … The psxy command will plot a point on the map (first line below psxy, longitude (lon) and latitude (lat) information) and connect it to another point (second line, lon and lat) and so on. For this command we will follow the -R, -J of the previous command.
Now the information seems quiet complete. But, you may add some important information (tectonically), such as topography and trench location. To do that, we revise this script, and now we go to B much more advance than A.
B. We will use a grid file (.nc or .grd) from the free available resources. Here we use a GEBCO 30-arc second in our map. You can just follow the instruction in GEBCO website to get the .nc or .grd file. The other information that we need is a trench location (for a subduction area), the file is available at the USGS website and please see or follow the instruction which file that you need in your map at this website.
After you downloaded these two files, make sure that these files are located within the same folder as your .gmt script file located and please make sure that the grid file that you downloaded from GEBCO was using the same file name as well as the trench clip file, i.e. in the script below: gebco_1.nc and sumatra_trench.clip, respectively.
Now you may try to execute this script below.
#!/bin/csh
#environment setting
gmtset PAPER_MEDIA A3
gmtset LABEL_FONT_SIZE = 12p
gmtset BASEMAP_TYPE plain
gmtset ANNOT_FONT_SIZE_PRIMARY = 12p
gmtset ANNOT_FONT_SIZE_SECONDARY = 12p
gmtset HEADER_FONT_SIZE = 14p
###resample the grd file to fit in range of regional map
grdsample gebco_1.nc -Gsumatra.grd -R95/102/-2/6 -I0.005 -V
###set color for the grid file
makecpt -Cglobe -T-5775/3590/10 -Z -V > sumatra.cpt
###grd image
grdimage sumatra.grd -R95/102/-2/6 -JM4i -Ba1f0.5NSEw -Csumatra.cpt -X11.7i -Y6.8i -K -V > nearshore.ps
###regional map
pscoast -R -JM -B -Dh -W0.5p,black -K -O -V >> nearshore.ps
###input sumatra_trench
psxy sumatra_trench.clip -R -B -J -Sf0.5i/0.1ilt -Gred -W0.03i,red -K -O -V >> nearshore.ps
###make a rectangle indicate for the study are
psxy -R -J -W0.03i/red -K -O -V << EOF >> nearshore.ps
95.18 5.5
95.55 5.5
95.55 5.75
#> next (this gonna disconnect the points)
95.18 5.75
95.18 5.5
EOF
###study area map
pscoast -R95.18/95.55/5.5/5.75 -JM10i -Ba0.05f0.01NSEW -Dh -W0.5p,black -X-11i -Y-6.5i -O -V >> nearshore.ps
gv nearshore.ps
After the program run, you may get a figure as shown below.
Now your regional map already have the topographic view and a trench location.
Institute of Earth Sciences Coders 