GMTについて

2011/02/26 筑波大学農林工学系 奈佐原(西田)顕郎

はじめに

参考になるページ:

GMTのインストール (Ubuntu Linux 10.04)

# インストール
$ sudo apt-get install gmt gmt-coast-low gv imagemagick csh

# 海岸線データのダウンロードと配置
$ wget ftp://ftp.soest.hawaii.edu/gmt3/GMT_full.zip
$ wget ftp://ftp.soest.hawaii.edu/gmt3/GMT_high.zip
$ unzip GMT_full.zip
$ unzip GMT_high.zip
$ sudo mv GMT/share/*cdf /usr/share/gmt/coast/
注意: この海岸線データは, 古いバージョンのGMTのためのものである。Ubuntuに自動インストールされるGMYはバージョンが古いので, これらが必要である。最新バージョンの海岸線データは使えない。

# 環境変数の設定
## /etc/profileまたは~/.bashrcに, 以下を記述する:
export NETCDFHOME=/usr/lib
export GMTHOME=/usr/lib/gmt
export PATH=$PATH:$GMTHOME/bin
export MANPATH=$MANPATH:$GMTHOME/man
そしてシェルを再起動。

## 印刷時のスケールをインチからcmに変える:
$ gmtdefaults -L | sed 's/inch/cm/' > ~/.gmtdefaults

使い方

GMTでは, 投影法を設定したり海岸線を描いたりマークを入れたり色を塗ったりスケールを入れたり文字を入れたり, というようなことを, ひとつひとつ, コマンドで行ない, その実行結果を, リダイレクトでポストスクリプト形式のファイルにどんどん書きこむ。コマンドはコマンドラインからぽちぽちと打ってもよいのだけど, それではめんどくさいし再利用ができないので, コマンドを一括してシェルスクリプトにまとめることが多い(GMTではなぜかCシェル(csh)のスクリプトがよく使われる。別にBシェル(shとかbash)のスクリプトでも同じことができるのだが, 慣例なので, このウェブサイトではCシェルで記述する)

GMTから吐きだされたポストスクリプトファイルは, そのまま(ポストスクリプト対応の)プリンターに流せばきれいな絵として印刷される。論文やプレゼンでもそのままで使うこともできるが, 多くの場合は, これをPNGなどの画像形式に変換するほうが, パワーポイントに張り付けたりウェブに張り付けたりできて便利である。そのためには, ImageMagickという画像処理ソフトの中のconvertというコマンドを使うのが, 最も手っ取り早い。


使ってみよう(1)

いきなりだが, 以下のコマンドを走らせてみよう:

$ pscoast -R0/360/-90/90 -JG135/35/10c -G255/168/0 -S34/34/220 -A1000 -P > globe.eps
$ convert globe.eps -crop 290x290+70+485 globe.png
$ display globe.png
左のような地球の絵が表示されればOK。

うまくいかないときは, 以下をチェックしよう。


使ってみよう(2)

以下のコマンドを走らせてみよう:

$ pscoast -R-30/30/-40/40 -Jm0.1i -B5 -I1/1p,blue -I2/0.25p,blue -N1/0.25p,- -W0.25p,white -Ggreen -Sblue  -P > test.ps
$ gv test.ps

↑このように, アフリカ大陸の西側7割ほどが緑色で表示され, 国境線が点線, 河川と海が青線で描かれていればOK。これはGMTのマニュアルに載っている例です。

上のコマンドで, 「-R-30/30/-40/40」の部分を, 「-R100/150/0/50」と変えて, やり直してみよう。すると,

↑このように, 東南アジア・中国・日本あたりが表示されればOK。これでわかったように, -Rオプションは描画範囲を指定する。

では, 上のコマンドで, -Rオプションの部分を, 「-R130/135/30/35」と変えて, やり直してみよう。九州・四国あたりが表示されるが, なんだかとても小さい。もっと拡大したい。
そこで, 「-Jm0.1i」を, 「-Jm0.8i」と変えて, やり直してみよう。九州・四国あたりがだいぶ拡大されて表示されればOK。でも, 海岸線がラフ過ぎて汚い。

そこで, 「-Jm0.8i」のあとに(ほんとはどこでもいいのだが), 「-Df」というオプションを付け加えて, やり直してみよう:

pscoast -R130/135/30/35 -Jm0.8i -Df -B5 -I1/1p,blue -I2/0.25p,blue -N1/0.25p,- -W0.25p,white -Ggreen -Sblue  -P > test.ps
$ gv test.ps
きれいになった。でも海岸線の白い線が多少うるさいので, これを無くしてしまおう。そのためには-Wオプションをとってしまえばよい:
pscoast -R130/135/30/35 -Jm0.8i -Df -B5 -I1/1p,blue -I2/0.25p,blue -N1/0.25p,- -Ggreen -Sblue  -P > test.ps
$ gv test.ps
これで↓のようなきれいな九州・四国の地図ができたはず。


このように, GMTでは, コマンドのオプションがたくさんあって, それを適当に調節することで, 望みの地図を描くことができる。


メモ

geographic projection (plate caree)で描くには? ... -Jxオプションを使えば良い!


GMTむけの県境データの整備

日本の地図なら県境があるほうがいいので, 県境のベクターデータを, ここから頂いて解凍。

すると, map.japanというテキストファイルが現れる。

この前半部の海岸線のデータはGMTの海岸線データと微妙にずれるので重ねると格好が悪い。そこで, この部分だけ削除する。

map.japanの, 海岸線データと県境データのさかいめ部分:

   25.8762  131.2539
   25.8765  131.2240
   999.9900    0.0   
   999.9900  999.9900
   999.9900    0.0   
   999.9900    0.0   
   999.9900    0.0   
   999.9900    0.0   
   999.9900    0.0   
   999.9900    0.0   
   999.9900    0.0   
   999.9900    0.0   
   40.4265  139.9323
   40.4288  139.9510
   40.4319  139.9615
awkなどで自動削除するには, 15900行まで削除すればよい。

つぎに, 県境データどうしのさかいめを

>   5
などに置換。

最後に, map.prefという名前で保存すればよい。

これらの作業を自動的にやるのが, 以下のコマンド:

	wget ftp://ftp.eri.u-tokyo.ac.jp/pub/data/map/kotake.map.japan.tar.Z
	tar zxvf kotake.map.japan.tar.Z
	cat map.japan | awk 'NR>15900 {print $2,$1}' | sed 's/0.0 999.9900/>   5/' > map.pref
こうして作ったmap.prefを, $GMTHOME/share/などに入れておくのがよいかと。


GMTにおける標高データの整備

デフォルトではGTOPO30などの標高データの情報がないので,

wget http://home.hiroshima-u.ac.jp/kojiok/wingmt/grdraster.info
sudo cp grdraster.info /usr/local/GMT4.0/share/dbase/
次に, GTOPO30をダウンロードして解凍・変換:
wget ftp://edcftp.cr.usgs.gov/pub/data/gtopo30/global/*.gz
for i in *.gz; do tar zxvf $i; done
for i in *0*.DEM; do o=`echo $i|sed 's/DEM/dem/'|sed 's/E/e/'|sed 's/W/w/'|sed 's/N/n/'|sed 's/S/s/'`; mv $i $o; done
mv ANTARCPS.DEM antarcps.dem
for i in *0.dem; do dd if=$i of=${i%.dem}s.dem conv=swab; done
dd if=antarcps.dem of=antarcpss.dem conv=swab
sudo mv *0s.dem $GMTHOME/share/dbase/
sudo mv antarcpss.dem $GMTHOME/share/dbase/
rm *0*.SRC *0*STX *0*SCH *0*PRJ *0*DMW
grdraster 19 -R100/140/-10/40 -Ge100n40.grd
grdraster 23 -R140/180/-10/40 -Ge140n40.grd
grdraster 20 -R100/140/40/90 -Ge100n90.grd
grdraster 24 -R140/180/40/90 -Ge140n90.grd
grdpaste e100n40.grd e140n40.grd -Gd1.grd
grdpaste e100n90.grd e140n90.grd -Gd2.grd
grdpaste d1.grd d2.grd -GEastAsia.grd
rm d1.grd d2.grd e100n40.grd e140n40.grd e100n90.grd e140n90.grd
sudo mv EastAsia.grd $GMTHOME/share/dbase/
grdgradient $GMTHOME/share/dbase/EastAsia.grd -Nt0.7 -A240 -M -GEastAsia_I.grd
sudo mv EastAsia_I.grd $GMTHOME/share/dbase/
- 参考
- GMTでGTOPO30のデータを使って地形表示するには, まず, 必要な部分について, 該当する.demファイルから.grdというGMT専用のラスター形式のファイルに変換する(grdrasterコマンド)。もし対象領域が複数の.demファイルにまたがっている場合は, モザイクする(grdpasteコマンド)。 - 上の処理では, とりあえず我々がよく使う東アジアの領域について, EastAsia.grdというファイルを作り, $GMTHOME/share/dbase/に置いている。
- EastAsia_I.grdというのは, 地形に陰影をつけるためのデータ。
- EastAsia.grdの外については, 別個に同様の処理が必要(grdraster以下)。
- grdrasterコマンドの第一引数は, grdraster.infoに記述してある各地形データのID。


川上演習林位置 (シンプル)
#!/bin/csh
# map of the Kawakami Experimental Forest
psbasemap -Jm4.0 -R137.5/139.5/34.9/37.2 -Ba1f1 -X5.0 -Y5.0 -P -K -V >  Kawakami.ps
pscoast -Jm -R -Dh -P -V -K -W3 -Na -S200/230/250 -O >> Kawakami.ps
echo 138.5 35.917 0.5 | psxy -Jm -Sc -G250/0/0 -R -P -V -K -O >> Kawakami.ps
psxy $GMTHOME/share/map.pref -Jm -R -M -W3 -O >> Kawakami.ps

川上演習林位置 (ゴージャス)
スクリプト 参考

北海道の気象台の分布図
スクリプト

インドシナ半島中部の地形
スクリプト

地上到達PARの分布図
スクリプト


大学演習林地図と標高分布
スクリプト