gmt coupeの使い方

coupe — GMT 6.4.0 documentation

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

meca — GMT 6.4.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.241666737.5151347341033.3
141.838333338.185290481614
139.541666735.9866666711827048513.7
142.253333340.006666675683461255.9
142.2440.003333335314548324.1
141.726666738.185119223863.6
139.921666736.07546244141113.7
140.103333335.597517127765.9
140.09535.561666677119311743.5
141.938333342.389769741664.5
142.24539.986666675468431134.7
142.238333340.025212979-73.8
143.1142.35551250231374.7
138.40534.99333333256859362.7
132.2933.598333335030930-1143.3
141.546666737.6855232424203.6
134.67533.736666673718874-1733.5
139.041666735.4752035257233.6
135.236666734.023333335209321153.3
136.098333334.77166667563306943.6
135.618333335.071234243633.9
142.216666740.941666676126236-1004.6
141.691666738.171666675219023834.1
135.07533.731666674012181-1793.5
137.216666737.51166667132031674
142.916666741.533526740-1713.6
135.716666735.631666671015349173.4
141.9738.888333334918719743.8
138.4532.18356112211635.6
140.586666736.81833333719744-693.3
138.928333335.506666672218734183.8
142.91541.5516666736238471734
140.1534.28166667367315853.5
140.97536.5216666745214171323.5
138.926666735.5052218519103.1
141.401666737.1846213161064
141.83539.758333335920421903.8
141.863333338.266364571513.5
139.908333336.0633333345242181074.5
141.7137.75554204501334
140.061666735.65666667609757-234.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で超適当な位置の距離を求めているからですね。
地球は楕円体なので緯度によって距離が若干変わっています。

今回はここまで。

コメント

タイトルとURLをコピーしました