スポンサーリンク
スポンサーリンク

gmt coupeの使い方

スポンサーリンク
スポンサーリンク
スポンサーリンク
スポンサーリンク
coupe — GMT 6.6.0 documentation

gmt coupe は任意の断面の発震機構解を描くモジュールです。
gmt mecaを理解していることが前提ですので、gmt mecaについてよく分からない人は下記を見てください。

meca — GMT 6.6.0 documentation
gmt meca の使い方
GMTのmecaモジュールについて説明しています。

gmt coupe <データファイル名> -A[断面の切り方の設定] -S[発震機構解の種類][スケールの大きさ] [-E【伸張部分の色】] [-G【圧縮部分の色】] [-Fa] [-C] -JX[ヨコの長さ]/[タテの長さ] -B[x|y|z][a][f][g][+l”label”] -B[WSEN]
※ スケールの大きさはM5の時の大きさ。
※ “-E”を指定しない場合は白色で塗りつぶす。
※ “-G”を指定しない場合は黒色で塗りつぶす。
※ “-Fa”を指定するとT軸・P軸もプロットされる。
※ カラーパレットを用いて色を付ける場合は“-C”を付ける。色は3列目が参照される。

スポンサーリンク
スポンサーリンク

gmt meca

復習です。気象庁のサイトに飛んで発震機構解のデータを取得しましょう。
今回は2021年10月の発震機構解データを取得してみます。

気象庁|地震・津波|初動発震機構解
初動発震機構解のページです
経度 緯度 深さ 走向 傾斜 すべり角 マグニチュード
137.2416667 37.515 13 47 34 103 3.3
141.8383333 38.18 52 90 48 161 4
139.5416667 35.98666667 118 270 48 51 3.7
142.2533333 40.00666667 56 83 46 125 5.9
142.24 40.00333333 53 145 48 32 4.1
141.7266667 38.18 51 192 23 86 3.6
139.9216667 36.075 46 244 14 111 3.7
140.1033333 35.59 75 171 27 76 5.9
140.095 35.56166667 71 193 11 74 3.5
141.9383333 42.38 97 69 74 166 4.5
142.245 39.98666667 54 68 43 113 4.7
142.2383333 40.02 52 129 79 -7 3.8
143.11 42.355 51 250 23 137 4.7
138.405 34.99333333 25 68 59 36 2.7
132.29 33.59833333 50 309 30 -114 3.3
141.5466667 37.685 52 324 24 20 3.6
134.675 33.73666667 37 188 74 -173 3.5
139.0416667 35.475 20 352 57 23 3.6
135.2366667 34.02333333 5 209 32 115 3.3
136.0983333 34.77166667 56 330 69 4 3.6
135.6183333 35.07 12 342 43 63 3.9
142.2166667 40.94166667 61 262 36 -100 4.6
141.6916667 38.17166667 52 190 23 83 4.1
135.075 33.73166667 40 121 81 -179 3.5
137.2166667 37.51166667 13 20 31 67 4
142.9166667 41.53 35 267 40 -171 3.6
135.7166667 35.63166667 10 153 49 17 3.4
141.97 38.88833333 49 187 19 74 3.8
138.45 32.18 356 112 21 163 5.6
140.5866667 36.81833333 7 197 44 -69 3.3
138.9283333 35.50666667 22 187 34 18 3.8
142.915 41.55166667 36 238 47 173 4
140.15 34.28166667 36 73 15 85 3.5
140.975 36.52166667 45 214 17 132 3.5
138.9266667 35.505 22 185 19 10 3.1
141.4016667 37.18 46 213 16 106 4
141.835 39.75833333 59 204 21 90 3.8
141.8633333 38.26 63 64 57 151 3.5
139.9083333 36.06333333 45 242 18 107 4.5
141.71 37.755 54 204 50 133 4
140.0616667 35.65666667 60 97 57 -23 4.2

まずはこれをプロットしてみます。

#!/bin/bash
gmt begin hogehoge png
    gmt basemap -JM12 -R120/150/20/50 -Bafg -BWSNE
    gmt coast -Df -W0.25 -Slightcyan -Glightgray
    gmt makecpt -Cseis -T0/100/1 -Z
    gmt meca "./mecha.csv" -Sa0.5 -C
    gmt colorbar -DJBR+jBL+o1.5/0+w-5/0.2 -Baf
gmt end

この結果が、

となります。ここまでは簡単。

発震機構解断面図を作成する

今回は発震機構解の断面図を描いてみます。

下半球投影

気象庁|地震・津波|震源球と断層面の関係
震源球と断層面の関係の解説ページです

上記のサイトを読んでみてください。

平面図では、発震機構解は下半球投影で描くのが慣例ですが、断面図を切る場合は横断図と縦断図にて投影方法が変えたい時もあると思います。下図を見てください

たとえば、図は左横ずれ断層ですが、断層の下面、横、そしてディスプレイの手前から奥側から各々投影する場合は発震機構解が変わってくるわけです。

これに対応するのがgmt coupeです。

gmt coupe

gmt coupe <データファイル名> -A[断面の切り方の設定] -S[発震機構解の種類][スケールの大きさ] [-E【伸張部分の色】] [-G【圧縮部分の色】] [-Fa] [-C] -JX[ヨコの長さ]/[タテの長さ] -B[x|y|z][a][f][g][+l”label”] -B[WSEN]
※ スケールの大きさはM5の時の大きさ。
※ “-E”を指定しない場合は白色で塗りつぶす。
※ “-G”を指定しない場合は黒色で塗りつぶす。
※ “-Fa”を指定するとT軸・P軸もプロットされる。
※ カラーパレットを用いて色を付ける場合は“-C”を付ける。色は3列目が参照される。

gmt coupeは断面図を描くモジュールですが、-Rを必須としません。
GMTが勝手に範囲を設定し断面図を作成します。
(もちろんこだわりたければ-Rを使うことも可能です)

断面の切り方には色々あります。たとえば

  • 任意の2地点を結ぶ線(-Aa)
  • ある1点を基準として方向角と長さを指定する方法(-Ab)

などがありますが、ここでは-Aaについて紹介します。

-Aa

-Aa[始点経度/始点緯度/終点経度/終点緯度/断面傾斜角/許容距離/最小深度/最大深度+f]

の順番で記載します。

ついでに、公式ドキュメントの書き方ではうまく行きませんでした。よく分かりません。

ということで、先ほどの平面図の断面を描いてみます。
今回は全体の断面を図示してみます。

#!/bin/bash
gmt begin hogehoge png
    gmt basemap -JM12 -R125/148/30/45 -Bafg -BWSNE+t"Kanto"
    gmt coast -Df -W0.25 -Slightcyan -Glightgray
    gmt makecpt -Cseis -T0/100/1 -Z
    gmt meca "./mecha.csv" -Sa0.5 -C
    gmt colorbar -DJBR+jBL+o1.5/0+w-5/0.2 -Baf
    gmt coupe "./mecha.csv"  -JX12/-6 -Bxafg+l"distance" -Byafg+l"depth" -BWSne -Aa125/40/148/40/90/1000/0/100+f -Sa0.5 -W0.25 -C -Y-8 
gmt end

この結果が、

となります。横軸は強制的に距離になります。

-Aa125/40/148/40/90/1000/0/100+f

  • 北緯40度東経125度から北緯40度東経148度までの、
  • 垂直(90°)で切った断面を描く。
  • 許容距離は上下1000kmまでで、
  • 0~100kmまでの深さで発生した地震を対象にする。

という意味です。なお、+fを付けないと図示されません。

-Ab

-Ab[基準経度/基準緯度/方位面線の距離/断面傾斜角/許容距離/最小深度/最大深度+f]

これはmapprojectと合わせ技なのかな、、と思っています。

gmt mapprojectの使い方
GMTのmapprojectモジュールについて説明しています。
gmt mapproject -G125/30+uk << END
> 148 30
> END

この結果が、148 30 2213.58090951なので、

gmt begin hogehoge png
    gmt basemap -JM12 -R125/148/30/45 -Bafg -BWSNE+t"Kanto"
    gmt coast -Df -W0.25 -Slightcyan -Glightgray
    gmt makecpt -Cseis -T0/100/1 -Z
    gmt meca "./mecha.csv" -Sa0.5 -C
    gmt colorbar -DJBR+jBL+o1.5/0+w-5/0.2 -Baf
    # gmt basemap -JX12/-6 -R139/141/0/100 -Bxafg+l"longitude" -Byafg+l"depth" -BWSne 
    gmt coupe "./mecha.csv"  -JX12/-6 -Bxafg+l"distance" -Byafg+l"depth" -BWSne -Ab125/38/90/2213.6/90/1000/0/100+f -Sa0.5 -W0.25 -C -Y-8 

となり、

上のような図が作成されます。
ここの横軸の距離は自身で設定した値にセットされます。

-Aaと断面図が若干異なるのは、mapprojectで超適当な位置の距離を求めているからですね。
地球は楕円体なので緯度によって距離が若干変わっています。

今回はここまで。

コメント

スポンサーリンク