africa.m

用matlab实现世界地图调用和显示

set input0 = TOPO5.grd
set input = africa.grd
set inputint = africa_intens.grd
set input2 = FPS.dat
set cpt = africa.cpt
set range = -20/60/-40/40
set proj = M12
set output = africa.ps
set frame = a10f5/a10f5NEWS
gmtset ANOT_FONT_SIZE 14 LABEL_FONT_SIZE 14;
#resample from TOPO5.GRD
#grdsample $input0 -G$input -R$range -I10m/10m
#make a color table
makecpt -Crainbow -T-6000/6000/50 -Z -V > $cpt
#compute data from output file
grdgradient $input -G$inputint -A0/180 -M -Nt0.75 -V
#plot map
grdimage $input -I$inputint -R$range -J$proj -B$frame -C$cpt -X4 -Y6 -K -V -P > $output
#contour
grdcontour $input -R -J$proj -C500 -G2 -E180/90 -K -O -V >> $output
#plot KM scale
psbasemap -B -J$proj -R -Lf50/-48/1/2000 -K -O >> $output 
#plot color scale
psscale -C$cpt -D6/-3/12/0.1h -B2000:Topography:/:m: -K -O -V >> $output
#plot coast line, river and political boundary
pscoast -R -J$proj -Dh -I1/2/0/0/255 -N1/2 -A2 -W2 -K -O -V >> $output
#plot beach ball
awk '{print $1,$2,$3,$4,$5,$6,$7,$8,$9,$10,$1,$2}' $input2 |psmeca -J$proj -R -Sm0.2 -G255/0/0 -K -O -V >> $output
#text
pstext -R -J$proj -G0/0/0 -N -K -O -V << END >> $output
20 60 30 0 10 CM Africa Topography
END

#symbol
#psxy -J$proj -R -Sc0.1 -G255/0/0 -N -O -K -V << END >> $output
#20 -42
#END

#text
pstext -R -J$proj -G0/0/255 -N -O -V << END >> $output
20 50 15 0 10 CM Earthquake 1975/01/01 ~ 2005/07/15
END

ghostview $output&
convert $output africa.jpg