Murakami's Memorandum

Home Page > My Memo               災害科学ホワイトボード   災害科学コース   理学部

(GMT) ETOPO1からデータを抽出

1分毎の高度データ ETOPO1から必要な領域のデータを切り出すスクリプト

# ====
# GMT command
# ubuntu 18.04LTS
# gmt ver.5.3.0
# ====
src=ETOPO1_Ice_g_gmt4.grd
range=132/135/33/35   # W/E/S/N (range of drawing)
size=15         # width of drawing (cm)
file=topo_shikoku     # part of drawing name
fig=fig_$file.eps      # output image file name
gmt grdcut $src -R$range -G$file.grd
gmt grd2xyz $file.grd > $file.dat   #ASCII形式で出力 (経度 緯度 高度)
#
#---------------------
gmt makecpt -Chaxby -T-500/2000/10 > $file.cpt
#gmt makecpt -Chaxby -T-2000/2000/500 > $file.cpt
gmt grdimage $file.grd -R$range -JM$size -C$file.cpt -P -K > $fig
#gmt grdcontour $file.grd -R -J -C500 -W0 -L0/5000 -A500 -O -K >> $fig
gmt grdcontour $file.grd -R -J -C200 -W0 -L0/5000 -A- -O -K >> $fig
gmt psscale -Ba2500g500f500 -C$file.cpt -D6c/-1c/6c/0.3ch -O -K >> $fig
gmt pscoast -R -J -N1/2p,#000000,2_2:0 -I1/1p,#0000ff -I2/1p,#0000ff -O -K >> $fig
gmt psbasemap -R -J -Bxa5f1 -Bya5f1 -O >> $fig
#---------------------
# sort data
sort -n -k 1 $file.dat > temp.dat
sort -n -k 2 temp.dat |gawk '{print $1,$2,$3}' >$file.dat2
#
#メルカトル図法で経度の範囲を15cmとしたときの (x,y)
gmt mapproject $file.dat2 -JM15 -R >$file.xyz.dat
#
# (経度、緯度) ===> (方位、距離) 基準(133.5,34)
#gmt mapproject $file.dat2 -G133.5/34+ue -Af > test.dat

最終更新時間:2020年06月28日 19時21分15秒


HomePage > FrontPage