Cross of Sumatra

Keterangan sudah dicantumkan dalam command program

#!/bin/csh
# & rubah menjadi <
gmtset ANNOT_FONT_SIZE_PRIMARY 10p HEADER_FONT_SIZE 18p PLOT_DEGREE_FORMAT ddd:mm:ssF
set region = 90/107/-10/7
set size = M11c
set psfile = 4cross.ps
set grdfile=./sumatra.nc
set cptfile=color.cpt
makecpt -Cglobe -Z > $cptfile
set AZ = 45

#Cross box near subduction
set boxlon = 99.0
set boxlat = -4.0
set boxdep = 0
echo $boxlon $boxlat $boxdep >! box.d

#Draw bathymetry and basemap
grdgradient $grdfile -A30/270 -Gintens.grd -Nt0.30 -V
grdimage $grdfile -R$region -J$size -C$cptfile -Iintens.grd -K -P -Ba4f1NWse -X2 -Y17 >! $psfile

gmtset ANOT_FONT_SIZE 9
pscoast -R$region -J$size -C$cptfile -Dl -W1 -P -Lf92.0/-9/40/60 -Tf92.0/-8.65/10i/3 -O -K  >> $psfile

#Draw cross centre
awk -F, '{ print $1, $2}' box.d | psxy -J$size -R$region -Sx.5/.5 -W2::0 -P -O -K  >> $psfile

# Draw tectonical structures
awk -F, '{ print $1, $2}' ../focalproposal/trench.gmt | psxy -R$region -J$size -Gdarkblue -Sc0.05 -O -K >> $psfile
awk -F, '{ print $1, $2}' ../focalproposal/transform.gmt | psxy -R$region -J$size -Gdarkblue -Sc0.05 -O -K >> $psfile
awk -F, '{ print $1, $2}' ../focalproposal/ridge.gmt | psxy -R$region -J$size -Gdarkblue -Sc0.05 -O -K >> $psfile

#Draw usgs data
awk -F, '{ print $1, $2}' ../focalproposal/gempa_kecil.lst | psxy -J$size -R$region -Sc0.1c -W0.1 -H -P -O -K >> $psfile
awk -F, '{ print $1, $2}' ../focalproposal/gempa_sedang.lst | psxy -J$size -R$region -Sc0.15c -W0.15 -H -P -O -K >> $psfile
awk -F, '{ print $1, $2}' ../focalproposal/gempa_besar.lst | psxy -J$size -R$region -Sc0.2c -W0.2 -H -P -O -K >> $psfile


#Projection with the center
awk -F, '{ print $1, $2, $3}' ../focalproposal/trench.gmt  | project -Q -C$boxlon/$boxlat -A$AZ -Fxyzpqrs -L-10/550 -W-50/50 >! xyzpqrs_trench
awk -F, '{ print $1, $2, $3}' ../focalproposal/transform.gmt | project -Q -C$boxlon/$boxlat -A$AZ -Fxyzpqrs -L-10/550 -W-50/50 >! xyzpqrs_transform
awk -F, '{ print $1, $2, $3}' ../focalproposal/ridge.gmt  | project -Q -C$boxlon/$boxlat -A$AZ -Fxyzpqrs -L-10/550 -W-50/50 >! xyzpqrs_ridge

awk -F, '{ print $1, $2, $3}' ../focalproposal/gempa_kecil.lst  | project -Q -C$boxlon/$boxlat -A$AZ -Fxyzpqrs -L-10/550 -W-50/50 >! xyzpqrs_kecil
awk -F, '{ print $1, $2, $3}' ../focalproposal/gempa_sedang.lst  | project -Q -C$boxlon/$boxlat -A$AZ -Fxyzpqrs -L-10/550 -W-50/50 >! xyzpqrs_sedang
awk -F, '{ print $1, $2, $3}' ../focalproposal/gempa_besar.lst  | project -Q -C$boxlon/$boxlat -A$AZ -Fxyzpqrs -L-10/550 -W-50/50 >! xyzpqrs_besar


#Plot projected earthquakes  (red dots) and mainshock (yellow star) on map
awk '{print $1, $2}' xyzpqrs_kecil | psxy -J$size -R$region -Sc0.1c -W0.1/red -P -O -K >> $psfile
awk '{print $1, $2}' xyzpqrs_sedang | psxy -J$size -R$region -Sc0.15c -W0.15/red  -P -O -K >> $psfile
awk '{print $1, $2}' xyzpqrs_besar | psxy -J$size -R$region -Sc0.2c -W0.2/red  -P -O -K >> $psfile

#Mainshock will appear in front of all earthquakes
#BE CAREFULL !!  MAINSHOCK FILE MUST CONSIST OF LON LAT DEPTH
awk '{print $1, $2}' ../focalproposal/cobi.d | psxy -J$size -R$region -Sa0.4c -Gyellow  -Wthin -P -O -K  >> $psfile
project -Q -C$boxlon/$boxlat -A$AZ ../focalproposal/cobi.d -Fxyzpqrs -L-10/550 -W-50/50 >! xyzpqrs_cobi

#Focalmec lines taken from true epicenter
sed -n '2,3p' ../focalproposal/test.txt | psxy -R$region -J$size -Ggray -W1 -O -K >> $psfile
sed -n '5,6p' ../focalproposal/test.txt | psxy -R$region -J$size -Ggray -W1 -O -K >> $psfile
sed -n '8,9p' ../focalproposal/test.txt | psxy -R$region -J$size -Ggray -W1 -O -K >> $psfile
....................................

#Modified location of Global CMT's focalmec
psmeca ../focalproposal/focaltop_cmt.gmt -R$region -J$size -Sm0.25c/-1 -Gred -W1 -O -K >> $psfile
pstext ../focalproposal/cibi.txt -R$region -J$size -Gblack -O -K >> $psfile

#Drawing box and line of cross section on map
psxy box.txt -J$size -R$region -W0.5 -P -O -K >> $psfile
awk '{print $1, $2}' line.txt | psxy -J$size -R$region -W0.5 -P -O -K >> $psfile
awk '{print $1, $2, 14, 0, 1, "LT", $3 }' line.txt | pstext -J$size -R$region -Gdarkblue -P -O -K >> $psfile

#DRAW CROSSEC DIAGRAM
psbasemap -JX9.5/-5 -R0/550/0/300 -B50/30NWes -P -O -K  -Y-7 >> $psfile
pstext -JX -R  -Gdarkblue -P -O -K <> $psfile
0    0   14 0 1 LT C
520  0   14 0 1 LT C'
EOF

#Plot stars of structures, earthquakes, and mainshocks on crossec diagram
awk '{print $4, $3}' xyzpqrs_trench | psxy -JX -R  -Gdarkblue -Sc0.2 -P -O -K >> $psfile
awk '{print $4, $3}' xyzpqrs_transform | psxy -JX -R  -Gdarkblue -Sc0.2 -P -O -K >> $psfile
awk '{print $4, $3}' xyzpqrs_ridge | psxy -JX -R  -Gdarkblue -Sc0.2 -P -O -K >> $psfile

awk '{print $4, $3}' xyzpqrs_kecil | psxy -JX -R  -Sc0.1c -W0.1/red -P -O -K >> $psfile
awk '{print $4, $3}' xyzpqrs_sedang | psxy -JX -R  -Sc0.15c -W0.15/red -P -O -K >> $psfile
awk '{print $4, $3}' xyzpqrs_besar | psxy -JX -R  -Sc0.2c -W0.2/red -P -O -K >> $psfile

awk '{print $4, $3}' xyzpqrs_cobi | psxy -JX -R  -Sa0.4c -Gyellow -Wthin -P -O -K >> $psfile

psscale -C$cptfile -Iintens.grd -D0.4c/1/10/0.15 -B1000 -P -O -K -X11.5c -Y12c >> $psfile

# Create legend of magnitude size
cat << EOF >! magnitude.legend
N 3
S 0.025i c 0.04i red 0.010p 0.1i   mag 3.5-4.5
S 0.020i c 0.06i red 0.015p 0.1i   mag 4.5-5.5
S 0.015i c 0.08i red 0.020p 0.1i   mag 5.5-9.5
&
EOF

pslegend magnitude.legend -Dx0.1i/-0.2i/6.8i/1.5i/TC -Jx0.4i -R0/8/0/8 -O -X-4c -Y-5 >> $psfile

\rm -f *.d *.xy xyzpqrs_*
gs -sDEVICE=x11 $psfile