Main contents

☆主なコンテンツ
1、新着論文 2、論文概説 3、コラム 4、本のレビュー 5、雑記(PC・研究関連)
6、気になった一文集(日本語English) 7、日記(日本語English) おまけTwilog

2013年11月9日土曜日

GMTで等深線を描く

(2016/10/1 追記➡︎続編:「ETOPO1の標高データを使い、色付きの等深線をGMTで描く」)

MacでGMT(Generic Mapping Tool)を使っているのですが、ここ数日で等深線の描き方について色々と学んだことについて簡単にまとめておきます。
僕自身、GMTプロの沖野研の横にいながらずぶの素人なので、初心者で同じようなエラーに悩まされる人の参考になればと、記事を書いている次第です。

GMTのバージョンは4.5.8です。

今回やろうとしていたのは沿岸部の海岸線と等深線を描くという作業なのですが、たったそれだけの作業にかなりの時間をかけてしまいました。

使おうとしたデータはまず「ETOPO5.dat
アメリカ大気海洋庁(NOAA)HPに行くとダウンロードできるようになっています。
ちょっと古いデータということですが、緯度・経度5分間隔と、多少軽いデータとなっています。

ただし、これだとどうしても深度情報(z)にエラーが出てしまい、お手上げ状態となりました。(後で詳しい人に聞いたところによると、救済できる手段もあるのだとか)

「$ head ETOPO5.txt」で最初の10行を表示しても明らかに文字化けしている。

PCにあまり詳しくない僕はこの時点でお手上げ。

 〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜

次に手をつけたのが最新バージョンでとんでもなく重い(圧縮されて約1GB)、
Bedrock grid-registeredの
ETOPO1._Bed_g_int.xyz」ファイル。こちらは緯度・経度1分間隔の細かいデータです。

一応NOAA的には推奨されるファイルということになっています。(HP

「$ head ETOPO1._Bed_g_int.xyz」や「tail ETOPO1._Bed_g_int.xyz」で中を確認すると、今度はデータそのものはちゃんと表示されていそう。
そのまま「$ grdcut ETOPO1._Bed_g_int.xyz -R/120/155/20/50 -Gjapan.grd」でグリッドデータを切り出そうとするも、「Not a supported grid format(サポートされたグリッド形式ではない)」というエラーが返されてしまう。

〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜〜

それで色々と模索したあげく、辿り着いたのが「ETOPO1_Bed_c_gmt4.grd
要するにETOPO1がすでにGMT ver.4に対応した、grdファイルになったものだということです。
こちらのHPからダウンロードできます。ソースは同じくNOAAです。

かなりシンプルに、

pscoast -R${range} -JM${size} -D${reso} -W1 -G245 -S -B${xanot}/${yanot} -K -X3 -Y3 -A30 -Ia  > ${out}.eps 
grdcut ETOPO1_Bed_c_gmt4.grd -R -G${gridfile}.grd
grdcontour ${gridfile}.grd -JM -C10 -W0.2 -A20 -L-60/0 -O -K >> ${out}.eps

と、ETOPO1_Bed_c_gmt4.grdからある指定した範囲で切り出す(-R${range})という風にすると成功。

それまで、あれこれと
grdraste.infoの中身をいじったり、
「$ grdraster」でグリッドファイルの書き出し(参考HP
を試みるもうまくいきませんでした。
GMT用の質問サイトで「ETOPO1_Bed_c_gmt4.grdにはgrdrasterは必要ない」という記述を見つけて、上の方法でやったところすんなりいきました。

以下は「Kamichu 日本地図を描いてみる」を参考に描いた日本周辺の海底地形図。

range=120/155/20/50
#size is in cm scale
size=10
reso=i
xanot=a10f5
yanot=a10f5
out=japan
gridfile=japan
xy=japan_location

pscoast -R${range} -JM${size} -D${reso} -W1 -G245 -S -B${xanot}/${yanot} -K -X3 -Y3 -A30 -Lf74/18/18/200k > ${out}.eps 
grdcut ETOPO1_Bed_c_gmt4.grd -R -G${gridfile}.grd
grdcontour ${gridfile}.grd -JM -C2000 -W0.1 -A20 -L-12000/0 -O -K >> ${out}.eps
#psxy ${xy}.txt -JM -R -W -Sc6p -Gblack -O >> ${out}.eps


あとは細かいところを変えていけばいろんな地域の海底地形図を描けます。



ちなみに、

◎grdcontourの-Lオプションで「-12000/0」としてますが、これは等深線を描く範囲を-12000mから0m(海岸線)と指定しており、失くせば等高線も描けます。

◎japan_location.txtで緯度経度情報を含めたテキストファイルを用意しておくと、地図上に点を打つこともできます。
使うときは上のコマンドの#psxyの#を外してください。

◎同期に教えてもらった方法なのですが、コマンドをテキストエディタ(TextWrangler、GNU Emacs、X codeなど)で.shで保存しておくと、何度もコマンドを打つ必要がなく楽チンです。
読み出すときは「$ sh ファイル名.sh」でOK。

◎何度かSYNTAX errorというのを返されましたが(まだよく理解していませんが)、文字コードなど.shファイルの形式にも気をつけた方が良さそうです。



海洋底科学部門なので、もっと勉強しなきゃですね…沖野先生のGMT講座もしっかり見なければ…。。