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

【適宜更新】最近作成した図たち

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

ここでは最近GMTで作成した図を紹介します。

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

2024/4/28 オクラホマ

スクリプト

  1. 震源

    #!/bin/bash
    gmt begin oklahoma_hypo png
    gmt basemap -JM12 -R-103/-93/33/41 -Bafg -BWSNE+t"Oklahoma_2001-2023_USGS"
    gmt coast -Df -W0.25 -Glightgray -Slightcyan -LJBR+jTR+o0/1.5+w100+f -N2
    gmt makecpt -Cseis -T0/40/1 -Z
    awk -F"," '{print $3, $2, $4, $5*0.1}' query.csv | gmt plot -Sc -C -W0.25
    gmt colorbar -DJBR+jBL+o1.5/0+w-5/0.2 -Ba10f10+l"Depth (km)"
    gmt inset begin -DjBR+w2.5/1.4+o0.1c -F+gwhite+p0.25 -M0.01
        gmt set MAP_FRAME_TYPE = "plain"
        gmt basemap -JM2.5 -R-125/-65/24.6/50 -Bafg -BWSNE
        gmt coast -Df -W0.25 -Glightgray -Slightcyan -LJBR+jTR+o0/1.5+w100+f -N2
        gmt plot -W0.25,blue -L <<END
    -103 33
    -93 33
    -93 41
    -103 41
    END
    gmt inset end
    gmt basemap -JX12/-12 -R-103/-93/0/40 -Bxafg+l"longitude" -Byafg+l"depth" -BWSne -Y-15
    awk -F"," &#039;{print $3, $4, $4, $5*0.1}&#039;  query.csv | gmt plot -Sc -W0.01 -C
    gmt colorbar -DJBR+jBL+o1.5/0+w-5/0.2 -Ba10f10+l"Depth (km)"
    gmt end
  2. N-T図

    #!/bin/bash
    input="query2.csv"
    # ${input}のうち2001から2023が含まれる行の数を数える
    gmt begin count png
    awk -F"," &#039;{print $1}&#039; ${input} | gmt histogram -JX12 -R2001/2023/0/4000 -Bxa1fg+a45+l"Year" -By1000a+l"Count" -BWSne+t"Oklahoma_2001-2023_USGS_M2.5+"  -Glightgray -W -D -T1
    gmt inset begin -DjTL+w5/2.8+o0.1c -F+gwhite+p0.25 -M0.01
        gmt set MAP_FRAME_TYPE = "plain"
        gmt basemap -JM5 -R-125/-65/24.6/50 -Bafg -BWSNE
        gmt coast -Df -W0.25 -Glightgray -Slightcyan -LJBR+jTR+o0/1.5+w100+f -N2
        gmt plot -W3.5,blue -L <<END
    -103 33
    -93 33
    -93 41
    -103 41
    END
    gmt inset end
    gmt end

    謝辞

    データ出典:USGS
    記して感謝申し上げます。USGSはデータ取得しやすい!便利!!

2024/4/17 豊後水道

スクリプト

  1. 震源
#!/bin/bash
input="./h2023.csv"
gmt begin bungo_suido png
    gmt basemap -JM12 -R132.2/132.6/33.0/33.4 -Bafg -BWSNE+t"2023/01/01-2023/12/31_JMA"
    gmt coast -Df -W0.25 -Glightgray -Slightcyan
    gmt makecpt -Cseis -T0/50/1 -Z
    awk -F"," &#039;{print $11,$9,$13,$15*0.1}&#039; ${input} | gmt plot -Sc -W0.25 -C
    gmt plot -Sa -W0.25 -C << END
132.4 33.2 41 0.62
END
    gmt colorbar -DJBR+jBL+o1.5/0+w-5/0.2 -Baf+l"depth"
    gmt basemap -JX12/-12 -R132.3/132.5/0/50 -Bxafg+l"longitude" -Byafg+l"depth" -BWSne -Y-15
    awk -F"," &#039;{print $11,$13,$13,$15*0.1}&#039; ${input} | gmt plot -Sc -W0.01 -C
    gmt plot -Sa -W0.25 -C -l"202404172314" << END
132.4 41 41 0.65
END

    gmt colorbar -DJBR+jBL+o2.5/0+w-5/0.2 -Baf+l"Depth"
gmt end
  1. 発震機構解
#!/bin/bash
# ローカルのファイル先を指定したりURLを入力してください。
input="./mecha/mecha2011-2023.txt"
rm -f mecha.txt mecha_aki.txt
# データ収集(URLから取得する場合のみ)
# wget $url -O $input
## 安芸・リチャーズ
## 1,2行目は削除してカンマ区切りにする
sed -e 1,2d $input | sed -e &#039;s/  */,/g&#039; | awk -F"," &#039;{print $4,$3,$5,$8,$9,$10,$6}&#039; > mecha.txt
## 負の値の対応
sed -e &#039;s/;[0-9]*//g&#039; mecha.txt | sed -e &#039;s/-([0-9]*)-([0-9]*)/-1/g&#039; > mecha_aki.txt
# GMT
gmt begin bungo_suido_mecha png
    gmt basemap -JM12 -R132.2/132.6/33.0/33.4 -Bafg -BWSNE+t"2011-2023_F-net"
    gmt coast -Df -W0.25 -Glightgray -Slightcyan -LJBR+jTR+o0/1+c40+w10+f
    gmt makecpt -Cseis -T0/50/1 -Z
    gmt meca "mecha_aki.txt" -Sa0.5 -C
    gmt meca -Sa0.5 -C << END
132.4 33.2 41 166.1 69.6 -114 6.2 0 0 202404172314
END
    gmt colorbar -DJBR+jBL+o2/0+w-5/0.2 -Baf+l"Depth(km)"
gmt end

注意事項

  • 震源の位置はAQUAシステム メカニズム解カタログに基づくため今後変わる可能性があります。
  • 2023年の気象庁一元化震源データは今後見直される可能性があります。

謝辞

以下のデータを使用させていただきました。ここに記して感謝申し上げます。

  • 震源データとして、気象庁の「地震カタログ」を一部改変したもの。
気象庁|震源データ
気象庁が提供するページです
  • 発震機構解のデータとして、国立研究開発法人防災科学技術研究所のF-netデータ
防災科研機関リポジトリ
  • HinetPy

Dongdong Tian, 2020 , HinetPy: A Python package to request and process seismic waveform data from Hi-net. , doi: 10.5281/zenodo.3885779

2024/8/8 日向灘

スクリプト

#!/bin/bash
# 131.8593,31.704,32,0.35
rm danmen.txt danmen.csv
magnitude="2"
input="./h1998-h2024_m${magnitude}over.csv"
lon1="130"
lon2="134"
lat1="31"
lat2="34"
output="1998-2024_m${magnitude}over"
title="1998-2024 Dep<100km M${magnitude}>"
depth_max="100"
gmt set GMT_VERBOSE="long"
gmt begin $output png
    gmt basemap -R$lon1/$lon2/$lat1/$lat2 -JM12 -Bafg -BWSNE+t"$title"
    gmt coast -W0.25 -Df -LJBR+jTR+w100+f+o0/1 -Slightcyan -Glightgray
    gmt makecpt -Cseis -T0/$depth_max/1 -Z
    # 9カラム目が緯度、11カラム目が経度、13カラム目が深度、15カラム目がマグニチュード
    # awkに変数を渡して13カラム目が$depth_max以下の行を抽出
    # 15カラム目が小さい順に並びかえる
    awk -F, -v depth_max="${depth_max}" &#039;$13 <= depth_max {print $11, $9, $13, $15*0.05}&#039; "$input"| sort -t &#039; &#039; -n -k4 | gmt plot -Scc -W0.25 -C
    # 20240808
    echo "131.8593,31.704,32,0.35" | gmt plot -Ss -W0.25 -C
    #
    gmt colorbar -DJBR+jBL+o1.5/0+w-5/0.2 -Baf+l"Depth (km)" 
    # lon1 lon2 lat1 lat2の中の震源の個数を図の左下にN=で記載する
    awk -F, -v lon1="${lon1}" -v lon2="${lon2}" -v lat1="${lat1}" -v lat2="${lat2}" &#039;$11 >= lon1 && $11 <= lon2 && $9 >= lat1 && $9 <= lat2 {n++} END {print lon1, lat1, "N="n}&#039; "$input" | gmt text -F+f12p,Helvetica-Bold+jBL -D0/0.1
    gmt legend -DJTL+jTL+o0.2c -F+p2p,black+gwhite << END
N 1
S 0.3 c 0.10 red 0.5p,black 0.8 M2
G0.25
S 0.3 c 0.15 red 0.5p,black 0.8 M3
G 0.25
S 0.3 c 0.20 red 0.5p,black 0.8 M4
G 0.25
S 0.3 c 0.25 red 0.5p,black 0.8 M5
G 0.25
S 0.3 c 0.30 red 0.5p,black 0.8 M6
G 0.25
S 0.3 c 0.35 red 0.5p,black 0.8 M7
END
    # 断面図の作成
    start_lon="131"
    start_lat="32.5"
    end_lon="133"
    end_lat="31.5"
    width="15"
    echo -e "$start_lon $start_lat n$end_lon $end_lat" | gmt plot -W1,blue
    # startからendまでの距離をmapprojectで計算
    distance=(`echo "$start_lon $start_lat" | gmt mapproject -G$end_lon/$end_lat+uk`)
    awk -F, -v depth_max="${depth_max}" &#039;$13 <= depth_max {print $11, $9, $13, $15*0.05}&#039; "$input"| sort -t &#039; &#039; -n -k4 | gmt project -C$start_lon/$start_lat -E$end_lon/$end_lat -W-${width}/${width} -Lw -Q > danmen.txt
    # danmen.txtをタブ区切りからカンマ区切りに変換
    sed &#039;s/t/,/g&#039; danmen.txt > danmen.csv
    gmt basemap -JX12/-12 -R${lon1}/${lon2}/0/${depth_max} -Bxafg+l"longitude" -Byafg+l"depth[km]" -BWSne -Y-14
    awk -F, &#039;{print $7,$3,$3,$4}&#039; "danmen.csv" | gmt plot -Sc -W0.25 -C -N
    # 20240808
    echo "131.8593,32,32,0.35" | gmt plot -Ss -W0.25 -C
    #
    # danmen.csvの行数を数えて図の左下に"N="で記載する。また、widthを記載する
    awk -F, &#039;END {print "N="NR, "Width="width}&#039; width="${width}km" "danmen.csv" | gmt text -F+f12p,Helvetica-Bold+cBL -D0.1/0.1
gmt end

注意事項

  • 震源の位置はAQUAシステム メカニズム解カタログに基づくため今後変わる可能性があります。
  • 2024年の気象庁一元化震源データは今後見直される可能性があります。

謝辞

以下のデータを使用させていただきました。ここに記して感謝申し上げます。

  • 震源データとして、気象庁の「地震カタログ」を一部改変したもの。
気象庁|震源データ
気象庁が提供するページです
  • 発震機構解のデータとして、国立研究開発法人防災科学技術研究所のF-netデータ
防災科研機関リポジトリ
  • HinetPy

Dongdong Tian, 2020 , HinetPy: A Python package to request and process seismic waveform data from Hi-net. , doi: 10.5281/zenodo.3885779

2024/9/1 富士山

高いぞ高いぞ…

スクリプト

  • 震源分布図
#!/bin/bash
input="./h2014-h2024zan.csv"
gmt begin mt_fuji pdf,png
    gmt set GMT_VERBOSE=long
    gmt basemap -R138.5/139.5/35/36 -Baf -BWSNE+t"Mt. Fuji(2014-2024)"
    gmt coast -Df -W0.25
    gmt makecpt -Cgeo -T-8000/8000/200 -Z
    gmt grdcut @earth_relief_15s -R138.5/139.5/34/36 -Gdem.nc
    gmt grdgradient dem.nc -Ggrad.grd -A45 -Ne0.8
    gmt grdimage dem.nc -Igrad.grd -C
    gmt makecpt -Cseis -T0/100/1 -Z
    # 11カラム目は138-140, 9カラム目は34-36, 15カラム目の大きい順に並び変える
    awk -F"," '{if($11>=138.5 && $11<=139.5 && $9>=35 && $9<=36 && $13<=100) print $11,$9,$13,$15*0.1}' $input | sort -k4 -n | gmt plot -Scc -W0.25 -C
    # 138.9647, 35.4896 , 24.0, 4.3
    gmt plot -Sa -W0.25 -C << END
138.9647 35.4896 24 0.43
END
    gmt coast -Df -W0.25 -LJRB+jTR+o0/1.5+w50+f
    gmt colorbar -DJBR+jBL+w-5/0.2+o1.5/0 -Baf+l"Depth(km)"
gmt end
  • N-T図
#!/bin/bash
input="./h2014-h2024zan.csv"
north=35.6
south=35.3
west=139.3
east=138.9
# 9カラム目が緯度、11カラム目が経度なので、絞り込む
awk -F"," '{if($11>=138.5 && $11<=139.5 && $9>=35 && $9<=36 && $13<=100) print $0}' $input > range_hypo.csv
# 2カラム目が年、3カラム目が月、4カラム目が日なのでyyyy-mm-ddに変換
awk -F"," '{print $2"-"$3"-"$4","$0}' range_hypo.csv > range_hypo_date.csv
# 横軸に日付、縦軸にその日に発生した地震の回数をプロットする
awk -F"," '{print $1}' range_hypo_date.csv | sort | uniq -c | awk '{print $2","$1}' > date_count.csv
# 2カラム目を大きい順に並び替えて上位10行のみを取得して、range_hypo_date.csvからその日付のデータを抽出する
sort -t"," -k2 -n -r date_count.csv | head -n 10 | awk -F"," 'NR==FNR{a[$1];next}($1 in a){print $0}' - range_hypo_date.csv > top10.csv

gmt begin n-t pdf,png
    gmt set GMT_VERBOSE=long
    gmt basemap -JX24/10 -R2014-01-01/2024-08-31/0/500 -Bxa6Ofg+a45+l"Date" -Byafg+l"Count" -BWSne+t"2014-2024"
    awk  -F"," '{print $1}' range_hypo_date.csv | gmt histogram -T0.5 -W1
    gmt histogram top10.csv -T0.5 -W1 -D -N
    gmt inset begin -DjTR+w4/4+o0.1/0.1 -F+gwhite+p0.25 -M0.01
        gmt set MAP_FRAME_TYPE = "plain"
        gmt basemap -JM? -R138.5/139.5/35/36 -Bxy
        gmt coast -Df -W0.25 -Glightbrown
        gmt makecpt -Cseis -T0/100/1 -Z
        awk -F"," '{if($11>=138.5 && $11<=139.5 && $9>=35 && $9<=36 && $13<=100) print $11,$9,$13,$15*0.01}' $input | sort -k4 -n | gmt plot -Scc -W0.01 -C
        gmt plot -W0.5,blue -L << END
${west} ${north}
${east} ${north}
${east} ${south}
${west} ${south}
END
    gmt inset end
gmt end

注意事項

  • 震源の位置はAQUAシステム メカニズム解カタログに基づくため今後変わる可能性があります。
  • 2024年の気象庁一元化震源データは今後見直される可能性があります。

謝辞

以下のデータを使用させていただきました。ここに記して感謝申し上げます。

  • 震源データとして、気象庁の「地震カタログ」を一部改変したもの。
気象庁|震源データ
気象庁が提供するページです
  • 発震機構解のデータとして、国立研究開発法人防災科学技術研究所のF-netデータ
防災科研機関リポジトリ
  • HinetPy

Dongdong Tian, 2020 , HinetPy: A Python package to request and process seismic waveform data from Hi-net. , doi: 10.5281/zenodo.3885779

コメント

スポンサーリンク