英国でデータサイエンスを学ぶ

30代になってから海外で統計学・機械学習・プログラミングを勉強

MENU

GMTで地図を描く その1 基本的な使い方(メルカトル図法)

地図を描画するにはいろいろなソフトがありますけれど、個人的にはGMTで作成する地図が好きですね。
なんとなく綺麗な印象がありますので。



GMTとは何か、を詳しく知りたい方は、ハワイ大学のこちらのページが参考になります。
gmt.soest.hawaii.edu



以下は基本的な操作を備忘的に。



[目次]

1. GMTを用いたメルカトル図法での地図描画

1-1. イギリスの白地図を作成してみる(最小限の設定)

せっかくイギリスに留学していますので、イギリスの白地図を作成するところから。
世界地図だと広すぎますしね。
必要なコマンドは以下の通り。
ターミナルを開いて以下のコマンドを実行します。
1行目でポストスクリプトファイルを作成して、
2行目でポストスクリプトファイルをjpegファイルに変換しています。

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di > map1.ps
gmt psconvert map1.ps -E124 -A

たったこれだけのコマンドで次のような白地図が作成されます。


f:id:detailed-balance:20180408235626j:plain

こんなの簡単だろうと思う方にはこれから先の記事は役に立たないと思います。
ちなみに私の場合は、なんだこれは。。。。理解ができんと途方にくれました。

ということで詳しい中身を、自分が忘れないためにも備忘メモです。
計算するためのプログラム言語と違ってGMTって使用頻度少ないですからね。


2. GMTを用いたメルカトル図法での地図描画(海岸線の描画)

GMTでの海岸線の描画 : pscoast

GMTで海岸線を描画するには、

gmt pscoast

というコマンドを使用します。
ただ、これだけでは不十分ですので、付け加えるオプションを以下に説明していきます。

GMTでの描画範囲の指定 -R

まず、地図を描画する場合は、どの範囲での地図にするのかを決める必要があります。
具体的には、日本近郊の地図なのか、アメリカ大陸の地図なのか、または中央アジアの地図なのか、といったように範囲を指定してあげる必要があります。
GMTでは、

-R西経/東経/南緯/北緯

で指定することになります。
例えば、

-R-10/10/-30/45

と範囲指定すれば、これは西経10度東経10度南緯30度北緯45度で囲まれた範囲の地図を書くことになります。

GMTでの図法およびサイズの指定 -J

今回はメルカトル図法を用いています。
図法とサイズの指定は

-J図法サイズ

で行います。
例えば、

-JM16c

とすると、メルカトル図法で16cmの地図を作ることになります。

GMTで海岸線の太さを指定する -W

GMTで海岸線を描くときに太さを指定するには

-W太さ

で指定することになります。
太さ"2"であれば、

-W2

となります。

同一の画像で太さを1から5まで変更させた時の図は以下のようになります。
太すぎると少し見辛いですが、これは後述する海岸線の荒さ(解像度)と関連するので、それらを併せて調整するといいと思います。
f:id:detailed-balance:20180907015516j:plain

GMTの図を90度回転させる -P

GMTで描画すると、なぜか横方向に南北が走り縦方向に東西が走るような図になってしまいます。
なので、それを90度回転させて縦方向に南北が走り、横方向に東西が走るような図を作るため

-P

という指示をしてあげる必要があります。

GMTで海岸線の解像度を指定する -D

GMTで海岸線を描画する際には、どの程度海岸線の詳細を描画するか指定が可能です。
具体的には

-D解像度の指定

で行います。
今回は解像度をintermidiateで記載しているので

-Di

とします。

海岸線の解像度は

1. full
2. high
3. intermediate
4. low
5. crude

の5種類になります。

それぞれ図示してみると以下のようになります。
流石に-Dlや-Dcはちょっと荒すぎる気がしますが、-Df、-Dh、-Diならそれほど大きな差はない気がします。
f:id:detailed-balance:20180907020521j:plain

ここまでをまとめると

1. GMTで海岸線を描画し
2. 西経12度東経2度北緯49度北緯60度で囲まれた範囲を指定し
3. メルカトル図法で大きさ16cmの図とし、
4. 海岸線の太さを1とし
5. 90度回転させ
6. 海岸線の荒さをintermidiateにする
場合は、以下のコマンドで実行できます。

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di > map1.ps

GMTで作成したポストスクリプトファイルをjpegに変換

ポストスクリプトファイルはそれはそれでいいのですが、jpegの方が何かと使い勝手が良いので、jpegファイに変換してみます。
変換にはpsconvertというコマンドを使用します。
解像度の指定は

-E解像度

で行います。
例えば、124dpiにするのであれば

-E124

というように指定します。
さらに、このままでは余白が残ってしまうので、余白をカットするために

-A

というオプションを使用します。
これらをまとめると

gmt psconvert map1.ps -E124 -A

となります。


3. GMTで陸地と海に色を塗る -G, -S


f:id:detailed-balance:20180409001151j:plain
上記のような図を作成するには、以下のコマンドで実行します。

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue > map8.ps
gmt psconvert map8.ps -E36 -A

ここで、

-S色

で海の色を指定し、

-G色

で陸地の色を指定しています。



4. GMTで国境線を描画


f:id:detailed-balance:20180409002113j:plain
上記のような図を作成するには、

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 > map9.ps
gmt psconvert map9.ps -E36 -A

とコマンドを実行します。
国境線のオプションは

-N国境線のオプション

で指定し、国境線オプションは以下の4通りがあります。

  • 1 = National boundaries
  • 2 = State boundaries within the Americas
  • 3 = Marine boundaries
  • a = All boundaries (1-3)

5. GMT で細かな構造を省略

小さい構造物を省略したい時(例えば、小さな湖を表示しない時)には、-Aコマンドを使用します。

-A省略する規模(km)

例えば、1,000km未満の構造物を省略する時には、以下のようにコマンドを実行します。

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 > map10.ps
gmt psconvert map10.ps -E36 -A


f:id:detailed-balance:20180409003434j:plain

省略する規模でどういった違いがあるかの目安としては以下の図を参考にしてください。
f:id:detailed-balance:20180919095121j:plain

6. GMTで地図に目盛り・フレーム・グリッド線を入れる

GMTで作成した地図に緯度や経度の目盛り、グリッド線を入れるには

-Ba目盛りの刻みfフレームの刻みgグリッド線の刻み

とします。

例えば、5度ずつの刻みにする場合には、

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5 > map1.ps
gmt psconvert map1.ps -E36 -A

としてやればOKです。


f:id:detailed-balance:20180409004306j:plain



7. GMTで解像度を下げたことによる文字の粗さの修正

GMTで作成したポストスクリプトファイルをJpegに変換した際に、かなりファイルサイズを圧縮しています。
そこで、解像度を下げたことによる文字の粗さを修正するコマンドとして、

-Qt

をしようすると、文字の粗さを補完してくれます。

具体的には、

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5 > map12.ps
gmt psconvert map12.ps -E36 -A -Qt

とすると、次のような画像になります。


f:id:detailed-balance:20180409010121j:plain


8. 解像度を下げたことによる画像の粗さを修正

上記の文字の粗さを修正した場合と同様に、画像の粗さを修正する方法として

-Qg

を用います。

具体的には、

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5 > map13.ps
gmt psconvert map13.ps -E36 -A -Qt -Qg

とすると、次のような画像になります。


f:id:detailed-balance:20180409010442j:plain



9. GMTで画像に表題を付ける

GMT作成した画像に表題を付けるためには、

:."タイトル":

というコマンドを用います。

例えば、GMTで作成した地図にFig1という表題をつけたいとします。
その場合、次のようなコマンドを実行すれば、表題を付けることが可能です。

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5:."Fig1":> map14.ps
gmt psconvert map14.ps -E36 -A -Qt -Qg



f:id:detailed-balance:20180409012106j:plain


10. GMTで河川を描画する

GMTで河川を描画するには、

-I河川オプション

を用います。

オプションとしては、

1 = Permanent major rivers
2 = Additional major rivers
3 = Additional rivers
4 = Minor rivers
5 = Intermittent rivers - major
6 = Intermittent rivers - additional
7 = Intermittent rivers - minor
8 = Major canals
9 = Minor canals
10 = Irrigation canals
a = All rivers and canals (1-10)
r = All permanent rivers (1-4)
i = All intermittent rivers (5-7)
c = All canals (8-10)
があります。

ここでは例として、rオプション(All permanent rivers)を用いてGMTで河川を描画してみると、次のコマンドを用いれば

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5:."Fig1": -Ir > map15.ps
gmt psconvert map15.ps -E36 -A -Qt -Qg

以下のような画像が作成されます。



f:id:detailed-balance:20180409015551j:plain

GMTで河川の色を指定する

河川の色を指定したい場合には、

-I河川のオプション/色

で指定可能です。
例えば、青色にしたい場合には、

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5:."Fig1": -Ir/0/0/250 > map16.ps
gmt psconvert map16.ps -E36 -A -Qt -Qg

としてやると、次のように青色で河川が描画されます。


f:id:detailed-balance:20180919151126j:plain

色を変えた場合には、例えば、


f:id:detailed-balance:20180919152502j:plain


11. GMTでスケールを設置する

GMTの地図上でどのくらいの距離なのかを示すスケールを設置するには、

-Lスケール位置の経度/スケール位置の緯度/スケールの距離を計算する緯度/スケールの幅

のコマンドを用います。

例えば、今はイギリスの地図を描いているので、経度0度、緯度48度の位置にスケールを設置し、緯度48での距離の表示を行い、その距離は250km分を表示するとします。
その場合は、

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5:."Fig1": -Ir -Lf0/48/48/250 > map16.ps
gmt psconvert map16.ps -E36 -A -Qt -Qg

というコマンドになります。
この結果作成された図は以下のようになります。


f:id:detailed-balance:20180409015611j:plain

(注意)
ここで注意が必要なのは、2点間の距離というのは緯度が異なれば大きく変わってくるということです。
地球が楕円体であるために、上図の場合は北緯48度での250kmを表示しましたが、これが北緯0度では大きく変わることになります。


12. GMTで作図した地図にシンボルマークを描画する

GMTで作図した地図中に、例えば都市を示すマークを描画したい場合には、都市の緯度・経度を指定したファイルを用意し、それを用いて

-Sシンボルマークの種類・シンボルマークのサイズ

オプションで描画します。

例えば、ロンドンに星マークを、グラスゴーに丸マークを描画したい場合には以下のステップを踏む必要があります。

プロットするシンボルマークの経度・緯度をtxtファイルに保存

ロンドンの経度・緯度をlondon.txtファイルに保存し、

-0.07 51.3

グラスゴーの経度・緯度をgrasgow.txtファイルに保存します。

-4.15 55.5

ファイル名は自分の管理しやすいもので構いません。

GMTで図に重ねて描画する -O -K オプション

GMTでシンボルマークを描画するには、作図した図にシンボルマークを重ねて描画すると手続きが理解しやすくなります。
そこで、
最初に全体図を作成し-Kオプションを用いて次にシンボルマークを重ねることを宣言します。
次に、-Oオプションを用いてロンドンのシンボルマークを重ね、さらにグラスゴーのシンボルマークを重ねるために-Kオプションを再度追加します。
最後に、-Oオプションを用いてグラスゴーのシンボルマークを重ねます。

これらのステップを全て踏むと、以下のコマンドを用いることになります。

gmt pscoast -R-12/2/49/60 -JM16c -W1 -P -Di -GForestGreen -SLightBlue -N1 -A1000 -Ba5f5g5:."Fig1": -Ir -Lf0/48/48/250 -K > map181.ps
gmt psxy london.txt -JM16c -R-12/2/49/60 -Sa1 -Gred -O -K >> map181.ps
gmt psxy glasgow.txt -JM16c -R-12/2/49/60 -Sc1 -Gblue -O >> map181.ps
gmt psconvert map181.ps -E36 -A -Qt -Qg
-K:作製した図にさらに追記する場合に記載
-O:すでにある図に追記する場合に記載
london.txtにロンドンの経度・緯度を記載
glasgow.txtにグラスゴーの経度・緯度を記載
-Sa1:星マークでサイズ1
-Sc1:丸マークでサイズ1


f:id:detailed-balance:20180410224605j:plain