gmt5绘制网格数据(surface_grdimage)

任务

从文本文件（经度、纬度、高程）读取数据绘制高程图。

思路

先将文本文件转为gmt格网数据，再绘制格网

解决方法

1.数据转换

从文本（ascii）到格网数据，如果文本里面的数据是等间隔的，可以使用xyz2grd命令，如果数据是不等间隔的，可以使用surface命令，本处的数据是不等间隔的，所以使用surface，gmt命令如下：

# ubuntu gmt5
#R=320/325/45/50
R=300/345/40/55
data=hv_re.txt
grdfile=hv_re.grd

# ascii to grd
# http://modules.gmt-china.org/xyz2grd.html?highlight=xyz2grd
# http://gmt.soest.hawaii.edu/doc/5.4.1/xyz2grd.html
#gmt xyz2grd $data -G$grdfile -R$R -I15m -V

# http://gmt.soest.hawaii.edu/doc/5.4.1/surface.html
gmt surface $data -G$grdfile -R$R -I1m

2.绘图

使用grdimage命令可以绘制格网图，命令如下：

#R=320/325/45/50
R=300/345/40/55
J=M5i
ps=topo.ps
grdfile=hv_re.grd
cptfile=hv_re.cpt

# head
gmt psxy -J$J -R$R -T -K > $ps

# 制作cpt
gmt grd2cpt $grdfile -Ctopo -S-50/80/20 -Z > $cptfile

# http://modules.gmt-china.org/grdimage.html?highlight=grdimage
gmt grdimage $grdfile -R$R -J$J -B2/2 -K -P -O >> $ps

# tail
gmt psxy -J$J -R$R -T -O >> $ps

# save as pic
gmt psconvert $ps -Tg -P -E300 -A