gmt mapprojectの使い方

mapproject は地図上の緯度経度と紙面上の座標を行き来できるコマンドです。
任意の点が紙面の原点から何cmの場所に位置するか、などを計算することができます。
また測量計算もすることができます。

もっと詳しく

たとえば、地図を描く時を思い出してください。
メルカトル図法で日本地図を描く際は、以下のようなコマンドを打つと思います。

#!/bin/bash
gmt begin japan jpg
    gmt basemap -JM12 -R120/150/20/50 -Bafg -BWSNE
    gmt coast -Df -W0.25 -LJBR+jTR+o0/2+c20+w500+f
gmt end

このとき、この地図のタテの長さは何cmでしょうか？

地図のヨコの長さは12㎝ですが、タテの長さは分かりません(=指定できない)。
そこで、タテの長さを求めるには mapproject を使います。

これはイメージしやすいと思います。

任意の2点の直線距離を求めることができます。

紙面上/地図上の座標を取得

mapproject <経度・緯度データ> -J[投影法][ヨコの長さ] -R[描画範囲] [-I]
※ – I をつけると紙面座標から緯度経度を求めることができる
2点間の距離を求める

mapproject <経度・緯度データ> -G[基準経度]/[基準緯度]+u[単位]
※ 指定できる単位は以下の通り。
d : 度
m : 分
s : 秒
e : メートル
f : フィート
k : キロメートル
方位角を求める

mapproject -A[f/b] <固定点> <方位角を求めたいデータ>
※ -Afは固定点からの方位角、-Abは固定点への方位角

#!/bin/bash
gmt begin japan jpg
    gmt basemap -JM12 -R120/150/20/50 -Bafg -BWSNE
    gmt coast -Df -W0.25 -LJBR+jTR+o0/2+c20+w500+f
gmt end

たとえば先ほどの日本地図の右上の座標は東経150度・北緯50度でした。
なので、ここの(x,y)座標は

#!/bin/bash
gmt mapproject -JM12 -R120/150/20/50 << END
150 50
END

とコマンドを打つと、12 14.9303627441 と出力されるはずです。
すなわち、この地図はヨコ12㎝・タテ14.9303627441cmとなります。

逆に地図の左下から右に5cm・上に6cmの緯度・経度は、

#!/bin/bash
gmt mapproject -JM12 -R120/150/20/50 -I << END
5 6
END

とコマンドを打つと、132.5 33.4242416154と出力されます。
すなわち、右に5cm・上に6cmの緯度・経度は東経132.5度・北緯33.4242416154度となります。

たとえば、東経130度北緯40度～東経135度北緯42度の直線距離は、

#!/bin/bash
gmt mapproject -G135/42+uk << END
130 40
END

とコマンドを打つと、130 40 475.426828239と出力されます。
すなわち、東経130度北緯40度～東経135度北緯42度の直線距離は約475.42kmとなります。

さらに単位を変えてみてメートルにしてみると、

#!/bin/bash
gmt mapproject -G135/42+ue << END
130 40
END

とコマンドを打つと、130 40 475426.828239と出力されます。
確かに +uk の時に比べて1000倍されていることが分かります。