GMTによる作図

Last update: 2020.7.6


1. 世界の地震の分布図

世界の地震の震源分布図を作成する。
震源を表す丸はMにあわせて半径を変え、塗る色は半透明にする。
psxyの「-t」オプションで透明度をしてして半透明にすることができる。
ただし、PSファイルだと半透明とならずベタ塗りになってしまうので、PDFやJPEGに変換する必要がある。

準備
・地震カタログ
USGSのSearch Earthquake Catalogを利用して、希望の期間、マグニチュードの地震を検索する。
CSV形式でダウンウロードできる。
それを元にExcelや自作のシェルスクリプト、プログラムなどで、「経度・緯度・表示する震源の丸の大きさ」の３列のデータに変換する。
たとえば、次のようにする。
表示する震源の丸の大きさは、マグチュードの数値を元に適当な大きさに置き換えておく。
今回の図では「M」に対して「0.01*3^(M-5)」となるようにしている。
＜例：eq_sample.txt＞
longitude latitude mag_GMT
140.000  35.000  0.810
130.000  30.000  0.270

・標高データ ETOPO1
NOAA ETOPO1 Global Relief Model から
ETOPO1_Ice_g_gmt4.grd
をダウンロードする。


・標高データの配色用のカラーパレット
今回はGMTのソースの中にある「GMT_globe.cpt」を作業ディレクトリにコピーして利用する。

作図のスクリプト
gmt_world_earthquakes_map.sh
スクリプトの例
#!/bin/sh
#
# GMT script for world earthquakes map
# 2020.7.03 Takumi Murakoshi
#
gmtset PS_MEDIA A4           # A4 size
gmtset MAP_FRAME_TYPE PLAIN  # 地図の外枠の設定

F_EQ="eq_sample.txt"            # Earthquake Catalog (longitude latitude size_for_psxy)
F_PS="world_earthquakes_map.ps" # PS file
F_GRD="ETOPO1_Ice_g_gmt4.grd" # ETOPO1 grd file (https://www.ngdc.noaa.gov/mgg/global/global.html)

F_CPT="GMT_globe.cpt" # Color Palette file

COL1="255/0/0"      # Red
COL2="255/255/255"  # White
COLT="70"           # 透明度(%)

W1="0.1"    # psxyで描く円の線の太さ

jmap="Q"    # 正距円筒図法
scale="25c" # 地図の大きさ

Ropt="-30/330/-90/90" # 地図の範囲の指定
Aopt="5000"           # 5000 km^2 以下の面積の地形の非表示

PL="  "   # Landscape (-P = Portrait)

######################################################################
# grdimage ETOPO
######################################################################
grdimage ${F_GRD} -J${jmap}${scale} -R${Ropt} -C${F_CPT} -K ${PL}  > ${F_PS}

######################################################################
# psxy earthquakes
######################################################################
psxy ${F_EQ} -J${jmap} -R -Sc -G${COL1} -W${W1}/${COL2} -O    $PL -t${COLT} >> ${F_PS}

######################################################################
# Convert PS to PDF & JPEG
######################################################################
psconvert ${F_PS} -Tf
psconvert ${F_PS} -Tj


作図
カレントディレクトリに
gmt_world_earthquakes_map.sh
eq_sample.txt
ETOPO1_Ice_g_gmt4.grd
GMT_globe.cpt
のファイルがある状態で、
$ sh gmt_world_earthquakes_map.sh
とすると次の３種類の図のファイルができる。
world_earthquakes_map.ps
world_earthquakes_map.pdf
world_earthquakes_map.jpg


作図した例
・サンプルの震源データの場合
震源分布図のサンプル

・1990-2019年のM5以上の震源データの場合
震源分布図(1990-2019,>M5)


「計算機環境」へ戻る