gmtset FRAME_PEN 1.0p
gmtset FRAME_WIDTH 0.1c
gmtset HEADER_FONT 5
gmtset LABEL_FONT 5
gmtset HEADER_FONT_SIZE 10
gmtset LABEL_FONT_SIZE 12
gmtset ANOT_FONT_SIZE 10
gmtset BASEMAP_TYPE  plain
gmtset TICK_LENGTH   -0.1c
gmtset TICK_PEN  1.0p
gmtset LABEL_OFFSET 0.1c
gmtset MAP_TITLE_OFFSET 5p
set ps=yuanzhen.ps
set scale=2.4e8
set eq=eq.txt
set lon=103.6
set lat=36.0            
rem grdmath -R-180/180/2/4 -I6/0.1 X 4 MUL PI MUL 180 DIV COS Y 2 POW MUL = dist.nc
grdmath -Rg -I1 %lon% %lat% SDIST  = dist.nc
grdmath -Rg -I1 %lon% %lat% SBAZ  = baz.nc
rem grdmath -Rg -I1 %lon% %lat% SDIST = dist.nc
rem pscoast -R-180/180/-90/90 -Je103.6/36/1:%scale% -B30g30 -Ggray -K -A2500 -Dh -X3.5i -Y1.25i > %ps%
pscoast -R0/360/-90/90 -Je103.6/36/1:%scale% -I2/0.25p,blue -N1/0.25p -W0.25p,white -Ggreen -Sblue  -K -A2500 -Dh -X3.0i -Y1.1i > %ps%
rem pscoast -R-180/180/-90/90 -JA103.6/36/6i -I2/0.25p,blue -N1/0.25p -W0.25p,white -Ggreen -Sblue  -K -A2500 -Dh -X3.0i -Y1.1i > %ps%
rem pscoast -R-180/180/2/4 -JP6i -I2/0.25p,blue -N1/0.25p -W0.25p,white -Ggreen -Sblue  -K -A2500 -Dh -X3.5i -Y1.25i > %ps%
grdcontour dist.nc -A60 -C30 -S8  -K -O -R -J  -Wathick,black -Wcthin,black>> %ps%
rem grdcontour dist.nc -B30Ns -P -C2 -S4 --PLOT_DEGREE_FORMAT=+ddd  >> %ps%
grdcontour baz.nc -C45 -S8  -K -O -R -J  -Wcthicker,white,-  >> %ps%

rem grdcontour dist.nc -A15+u+fblack -S8 -O -K -J -Wathin,black >> %ps%
gawk "{print $1,$2,$3*$3*$3/1000}" %eq% |psxy -: -R -J -O -K -Sc -Gred  -Wthinnest  >> %ps%
rem psxy -: %eq% -R -J -O -K -Sc0.08 -Gblack -Wthinnest  >> %ps%
echo %lon% %lat% | psxy -R -J -O -K -St0.7 -Gblack -Wthin >> %ps%
pause
del dist.nc .gmt*