因为我想把最近三年随同实验室出差的地方、路线一一记录下来，所以我写了这篇博客，目的是记录我用GMT绘制中国地形图并添加不同标识的过程。

Mac OSX 10.11.2

GMT version

GMT 4.5.13

data resolution

1-minute grids

Script language

perl

数据

对于gmt而言，只要有网格数据就能绘图。地球科学上常见的数据主要有两种：海岸线数据（e.g., netcdf）和地形数据（e.g., DEM）。

海岸线数据

一般来说在安装GMT的同时都会额外下载安装海岸线数据或者世界各国家边界线数据。

Seisman在GMT命令之用pscoast绘图海岸线如何用gmt绘制含海岸线的地图，这个跟MATLAB的m_map工具箱功能差不多。

地形图

对于地形图数据而言，所区分的就在于数据精度。分辨率越高，数据量越大，绘图过程越缓慢、内存消耗越多。Seisman在全球地形起伏数据总结中分别介绍了不同精度的全球地形数据，从5 minutes到1 second。之前用过SRTM的3弧秒的数据(.hgt format)绘制$3^o x 1^o$跨度的地形图，现在绘制全中国尺幅必须选用小精度的地形数据否则电脑跑不动啊。

为了省事，我选择了UCSD提供的1 mintue ASCII XYZ 70/135/15/55 地形图数据China_Topo.cgi。

下载的China_Topo.cgi是ASCII XYZ数据需要转换为GMT的grd网格数据

system "surface China_Topo.cgi -GChina_topo.grd -I0.5m -R$R -V ";

绘制地形图

# 计算梯度文件
system  "grdgradient China_topo.grd -Gintense.grd -A0/270 -Nt1.0 -V";
# 制作cpt文件
system  "makecpt -Cglobe -T-10000/10000/100 -Z -V >topo.cpt";
# grd文件成图
system  "grdimage -R$R -B$B -J$J -X4.5 -Y5.5  China_topo.grd -Iintense.grd -Ctopo.cpt  -K  -V >$PS";
# 绘制海岸线，国界边界数据
system  "pscoast -R$R -J$J  -O -K -Dh -A300 -I1 -W0.1p -N1/1p  >>$PS";

添加城市标识

添加城市标识就是利用psxy绘制地震台站点位,然后利用pstext在相应的位置添加文本信息，注意：文本信息位置应该稍稍偏离台站点位坐标点否则会有重叠。

# add location position
open (PSXY, "| psxy -J$J -R$R -B$B -St0.3c -W1p/black -Gblack -K -O >> $PS");
    print PSXY "   109.19  30.48   \n"; #enshi
    print PSXY "   120.17  36.01   \n"; #qingdao
    print PSXY "   114.17  30.35   \n"; #wuhan
    print PSXY "   113.4   34.46   \n"; #zhengzhou
    print PSXY "   84.51   45.36   \n"; #karamary
    print PSXY "   93.28   42.5    \n"; #hami
    print PSXY "   87.36   43.45   \n"; #urumqi
    print PSXY "   116.24  39.55   \n"; #beijing
    print PSXY "   121.29  31.14   \n"; #shanghai
    print PSXY "   103.51  36.04   \n"; #lanzhou
    print PSXY "   108.57  34.17   \n"; #xian
    print PSXY "   102.79  24.9    \n"; #kuming
    print PSXY "   115.98  39.49   \n"; #BGP
    print PSXY "   93.03   43.60   \n"; #balikun
    print PSXY "   120.16  30.25   \n"; #hangzhou
    print PSXY "   121.09  32.08   \n"; #suzhou
    print PSXY "   120.60  31.30   \n"; #nantong
    # print PSXY "   -122.40  37.78   \n";
close(PSXY);
# STRESS WUHAN CITY
open (PSXY, "| psxy -J$J -R$R -Sa0.4c -W1p/red -Gred -K -O >> $PS");
    print PSXY "    114.17  30.35                       \n";
close(PSXY);
# add location Name
open(PSTEXT, "| pstext -R$R -J$J -B$B -K -O -N  >> $PS");
    print PSTEXT "  125  31.14   10  0   2   2   Shanghai  \n";
    print PSTEXT "  109.19 30.80    10  0   2   2   Enshi  \n";
    print PSTEXT "  120.17 36.40    10  0   2   2   Qingdao  \n";
    print PSTEXT "  114.17  30.90   10  0   2   2   Wuhan  \n";
    print PSTEXT "  113.4   34.90   10  0   2   2   Zhengzhou  \n";
    print PSTEXT "  84.51   45.80   10  0   2   2   Karamary    \n";
    print PSTEXT "  93.28   41    10  0   2   2   Hami  \n";
    print PSTEXT "  87.36   43.99   10  0   2   2   Urumqi    \n";
    print PSTEXT "  116.8  41   10  0   2   2   BeiJing  \n";
    print PSTEXT "  103.51  36.99   10  0   2   2   Lanzhou  \n";
    print PSTEXT "  107.57  34.99   10  0   2   2   Xi'an  \n";
    print PSTEXT "  102.79  25.80    10  0   2   2   Kuming  \n";
    print PSTEXT "  115.0  39.90   10  0   2   2   BGP  \n";
    print PSTEXT "  93.03 43.960   10  0   2   2   Balikun  \n";
    print PSTEXT "  120.16 29   10  0   2   2   Hangzhou  \n";
    print PSTEXT "  119 31.31   10  0   2   2   Suzhou  \n";
    print PSTEXT "  121.8 32.8   10  0   2   2   Nantong  \n";
    # print PSTEXT "  -122.40  37.78  10  0   2   2   SanFranciso  \n";
close(PSTEXT);

添加旅行路线

添加旅行路线就像是添加断层信息

将分段路线坐标点位信息存储在文本2015.xy中

114.17  30.35
121.29  31.14
-122.40  37.78
>
114.17  30.35
121.29  31.14
120.16  30.25
>
109.19  30.48
114.17  30.35
>
120.17  36.01
114.17  30.35
>
114.17  30.35
103.51  36.04
93.28   42.5
93.03   43.60
>
116.24  39.55
114.17  30.35
>
102.79  24.9
114.17  30.35
>
115.98  39.49
114.17  30.35

绘制分段线段

# -m 操作选项通过'>'分隔符来判断分段标识
system  "psxy 2015.xy -J$J -R$R -W2p,blue -K -O -m'>' >> $PS";

GMT绘制中国地形图

Contents

数据