時系列データの作図
GMTは習うより慣れろ。例題と課題でGMTの使い方を理解しよう。
今回は時系列データを作図します。
例題
以下のようなグラフをGMT(5.4.3)で作成する。
GMTのスクリプト例を次に示す。
plot_exam04.sh
#!/bin/bash
#--------------------------------------------------
# plot for GMT5
# 1次元時系列グラフ
#--------------------------------------------------
#----- 前処理・設定 ------------------------------
# GMTパラメーター
gmt gmtset FONT_ANNOT_PRIMARY=9p,Helvetica # 主目盛りのフォント
gmt gmtset FONT_LABEL=10p,Helvetica # 軸タイトルのフォント
# 日付関係(=前後にスペースを入れない)
gmt gmtset FORMAT_DATE_MAP=mm/dd # 日付目盛り(yy[yy]=年,mm=月,dd=日,o=月名,jjj=日数)
gmt gmtset FORMAT_TIME_PRIMARY_MAP=a # 主目盛りの月・週表示：[f|a|c]
# input
in01=chiba_wea.txt # 入力ファイル
in01_h=1 # ヘッダー行数
in01_c=0,1 # 入力する列(0始まり）
# output
ps=plot_exam04.ps # 出力psファイル
# グラフのサイズ
mp=X # 地図の投影法(X:xy, M:メルカトル図法)
gx=10 # グラフの幅[cm]
gy=4 # グラフの高さ[cm]
#log_x=l # x軸を対数軸にするならl
#log_y=l # y軸を対数軸にするならl
jopt=$mp${gx}c${log_x}/${gy}c${log_y} # -Jオプション
# グラフの範囲(range)（-Rオプション）
rw=2017-01-01T00:00 # x軸の最小値(west)
re=2017-01-05T00:00 # x軸の最大値(east)
rs=980 # y軸の最小値(south)
rn=1030 # y軸の最大値(north)
# 軸・目盛りの設定
ti_x="Date (2017)" # x軸タイトル
ti_y="Pressure [hPa]" # y軸タイトル
tics_x=a1Df6h # x軸の目盛り(a=主,f=副)
grid_x= # x軸のグリッド(g)(*)
tics_y=a10f5 # y軸の目盛り(a=主,f=副)
grid_y=g10 # y軸のグリッド(g)(*)
axis=WeSn # 軸表示（大文字：軸＋目盛り、小文字：軸）
# 線の設定
line01=0.8p,black # 線：太さ,色,線種
#----- plot ------------------------------
echo "plot $ps"
# basemap(グリッド線のみ)
gmt psbasemap -J$jopt -R$rw/$re/$rs/$rn -Bx$grid_x -By$grid_y -B+n -P -K > $ps
# psxy(line)
gmt psxy $in01 -h$in01_h -i$in01_c -J -R -W$line01 -K -O >> $ps
# basemap(枠線)
gmt psbasemap -J -R -Bx$tics_x+l"$ti_x" -By$tics_y+l"$ti_y" -B$axis -O >> $ps
#----- 後処理 ------------------------------
# 画像変換(psconvert)の設定
pc_fmt=g # g=png, b=bmp, j=jpg, t=tiff, e=eps, f=PDF, s=SVG
pc_res=300 # 画像の解像度(DPI)
pc_mgn=0.1c # グラフ外側のマージン
pc_aa=2 # アンチエイリアシングの設定(1,2,4)
# png変換
gmt psconvert -A$pc_mgn -E$pc_res -T$pc_fmt -Qt$pc_aa -Qg$pc_aa $ps
# 中間ファイルの削除
rm -f gmt.history gmt.conf
使用したデータファイル。
なお、データは気象庁HPのデータを加工して作成した（実際の観測値を保証するものではありません）。
chiba_wea.txt
Date pres[hPa]
2017-01-01T00:00 1021.8
2017-01-01T01:00 1021.7
2017-01-01T02:00 1021.8
2017-01-01T03:00 1021.7
2017-01-01T04:00 1021.6
2017-01-01T05:00 1021.7
2017-01-01T06:00 1021.8
2017-01-01T07:00 1022.4
2017-01-01T08:00 1022.9
2017-01-01T09:00 1023.4
2017-01-01T10:00 1023.4
2017-01-01T11:00 1022.7
2017-01-01T12:00 1021.7
2017-01-01T13:00 1021.0
2017-01-01T14:00 1020.5
2017-01-01T15:00 1021.1
2017-01-01T16:00 1021.7
2017-01-01T17:00 1022.2
2017-01-01T18:00 1022.5
2017-01-01T19:00 1022.9
2017-01-01T20:00 1022.9
2017-01-01T21:00 1022.6
2017-01-01T22:00 1022.3
2017-01-01T23:00 1022.3
2017-01-02T00:00 1021.8
2017-01-02T01:00 1021.5
2017-01-02T02:00 1020.8
2017-01-02T03:00 1020.2
2017-01-02T04:00 1020.1
2017-01-02T05:00 1020.1
2017-01-02T06:00 1020.1
2017-01-02T07:00 1019.4
2017-01-02T08:00 1020.1
2017-01-02T09:00 1019.8
2017-01-02T10:00 1019.0
2017-01-02T11:00 1017.5
2017-01-02T12:00 1016.1
2017-01-02T13:00 1015.2
2017-01-02T14:00 1014.3
2017-01-02T15:00 1013.8
2017-01-02T16:00 1013.6
2017-01-02T17:00 1013.9
2017-01-02T18:00 1013.7
2017-01-02T19:00 1013.7
2017-01-02T20:00 1013.5
2017-01-02T21:00 1012.7
2017-01-02T22:00 1012.3
2017-01-02T23:00 1012.1
2017-01-03T00:00 1011.7
2017-01-03T01:00 1011.5
2017-01-03T02:00 1011.5
2017-01-03T03:00 1011.3
2017-01-03T04:00 1011.2
2017-01-03T05:00 1011.2
2017-01-03T06:00 1011.8
2017-01-03T07:00 1013.0
2017-01-03T08:00 1013.6
2017-01-03T09:00 1014.4
2017-01-03T10:00 1014.6
2017-01-03T11:00 1013.7
2017-01-03T12:00 1012.7
2017-01-03T13:00 1012.3
2017-01-03T14:00 1012.8
2017-01-03T15:00 1013.2
2017-01-03T16:00 1014.0
2017-01-03T17:00 1014.5
2017-01-03T18:00 1015.4
2017-01-03T19:00 1016.0
2017-01-03T20:00 1016.5
2017-01-03T21:00 1016.8
2017-01-03T22:00 1016.8
2017-01-03T23:00 1016.8
2017-01-04T00:00 1016.6
2017-01-04T01:00 1016.6
2017-01-04T02:00 1016.3
2017-01-04T03:00 1016.0
2017-01-04T04:00 1015.5
2017-01-04T05:00 1015.3
2017-01-04T06:00 1015.1
2017-01-04T07:00 1014.7
2017-01-04T08:00 1014.7
2017-01-04T09:00 1014.6
2017-01-04T10:00 1013.8
2017-01-04T11:00 1011.9
2017-01-04T12:00 1010.6
2017-01-04T13:00 1009.7
2017-01-04T14:00 1009.8
2017-01-04T15:00 1010.2
2017-01-04T16:00 1010.6
2017-01-04T17:00 1011.6
2017-01-04T18:00 1012.5
2017-01-04T19:00 1013.1
2017-01-04T20:00 1013.5
2017-01-04T21:00 1013.4
2017-01-04T22:00 1013.8
2017-01-04T23:00 1014.1
2017-01-05T00:00 1013.9
解説
時系列データのフォーマット
- 時系列データは、上記のデータのように、「YYYY-MM-DDThh: mm: ss」とする（日付と時刻の間に「T」を挟む）（絵文字に変換されるので、「:」のあとにスペースを挿入しているが、実際には不要）。
- この書式は、ISO8601形式の拡張形式と呼ばれる。
- 使用しない時間部分は省略できる。例えば、上記のデータでは秒の部分（:ss）を省略している（0秒とみなされる）。
目盛りの設定
- 時系列データの目盛りは、(1)-Bオプション（psbasemap）と、(2)GMTパラメータ（FORMAT_DATE_MAP、FORMAT_TIME_PRIMARY_MAP）で設定する。
-Bオプションの日付設定
-
軸の目盛り(annotation, tics, grid)の設定方法は、-Ba[目盛り数値の間隔]t[目盛り線の間隔]g[グリッド線の間隔]、という風に設定する（普通のデータと同様）。
-
時系列データの場合、[間隔]の設定には時間の単位も付ける。例えば、「-Bf1d」では、「目盛り線(tics)が1日(day)間隔」という意味になる。
-
間隔の単位は、y=年、o=月、d=日、h=時、m=分、s=秒、u=週（月曜はじまり）となる。月はmではなくoなので注意。
-
間隔の単位(y,m,d等)は大文字と小文字で、annotation（目盛り数値）の表記が異なることに注意！
-
間隔の単位を小文字にすると、annotationは数値で表記される。
-
間隔の単位を大文字にすると、annotationの表記は、パラメータのFORMAT_DATE_MAPの設定値になる。
-
例えば、-Bx1D6h（annotation部分を大文字）とした場合、x軸は次のようになる。
- これを、-Bx1d6h（annotation部分を小文字）とした場合は次のようになる。
- annotation部分を小文字にした場合は、FORMAT_DATE_MAPの設定に関わらず、目盛り数値が「01, 02」と日数表記になっている。
日付表示の設定（FORMAT_DATE_MAP）
- annotationなどの間隔の設定は、-Bオプションで設定するが、annotationの表記は、FORMAT_DATE_MAPで設定する。
- 表記のフォーマットは、yyyy=年（4桁）、yy=年（下2桁）、mm=月（数値）、o=月（月名）、dd=日、jjj=年初からの日数、を使って指定する。（※）週の表記はww=週数なのだが、自分の環境では表示できなかった（GMTのバグと思われる）。
- 例えば、「FORMAT_DATE_MAP=mm/dd」とすれば、1月10日は"01/10"と表記される（-Bオプションの間隔単位が大文字であること）。
- フォーマットの最初に、-（ハイフン）を付けると、2桁目の0がなくなる。例えば、1月10日を、「-mm/dd」とすれば"1/10"となる。
- 月を数値ではなく、月名（January, Februaryなど）で表示する場合は、パラメータのFORMAT_TIME_PRIMARY_MAPを設定する必要がある。
月名表示の設定（FORMAT_TIME_PRIMARY_MAP）
-
月名表示を設定するパラメータFORMAT_TIME_PRIMARY_MAPの設定値は、3種類ある。（※）週名も同様に設定できるらしいが、うまくいかなかったので、省略。
-
f(=full)：「January」のように、省略なしで表示する。
-
a(=abbreviated)：「Jan」のように、3文字の短縮形で表示する。
-
c(=character)：「J」のように、1文字の頭文字で表示する。
-
f, aを大文字で指定した場合は、単語が全部大文字になる。Fの場合は「JANUARY」のようになる。
-
例えば、FORMAT_DATE_MAP=-o.dd、FORMAT_TIME_PRIMARY_MAP=a、とした場合は、次のようになる。
時刻表示の設定（FORMAT_CLOCK_MAP）
- 日付の表示はFORMAT_DATE_MAPだが、時刻の表示はFORMAT_CLOCK_MAPで設定する。
- 表記のフォーマットは、hh=時、mm=分、ss=秒、を使って指定する。
- 例えば、「FORMAT_CLOCK_MAP=hh:mm」とすれば、2時15分は"02:15"となる。
- なお、日付と時刻を同じ軸で同時に表示することはできない（"MM/DD hh:mm"のようなフォーマットは無理）。どちらかを副軸を利用するなど、工夫が必要となる。
課題
課題(1)
- 以下に示すtemp.txtのデータを使って、下記のような気温のグラフを作成せよ。
- データはtab区切りで、日付(YYYY/MM/DD hh:mm)と気温[℃]のデータとなっている。
- 日付はGMTの形式（ISO8601形式）と異なるので、変換が必要。
- 温度記号「℃」は、"@~\260@~C"で出力できる。
- なお、データは気象庁HPのデータを加工して作成した（実際の観測値を保証するものではありません）。
temp.txt
Date-time Temp[℃]
2017/08/01 00:00 27.3
2017/08/01 01:00 27.2
2017/08/01 02:00 27.2
2017/08/01 03:00 26.8
2017/08/01 04:00 26.9
2017/08/01 05:00 26.8
2017/08/01 06:00 27.0
2017/08/01 07:00 27.6
2017/08/01 08:00 27.6
2017/08/01 09:00 28.4
2017/08/01 10:00 29.1
2017/08/01 11:00 28.9
2017/08/01 12:00 29.8
2017/08/01 13:00 30.6
2017/08/01 14:00 29.9
2017/08/01 15:00 27.7
2017/08/01 16:00 25.6
2017/08/01 17:00 25.1
2017/08/01 18:00 24.6
2017/08/01 19:00 24.2
2017/08/01 20:00 24.2
2017/08/01 21:00 24.2
2017/08/01 22:00 24.1
2017/08/01 23:00 23.6
2017/08/02 00:00 23.0
2017/08/02 01:00 22.7
2017/08/02 02:00 22.5
2017/08/02 03:00 21.9
2017/08/02 04:00 21.7
2017/08/02 05:00 21.6
2017/08/02 06:00 22.3
2017/08/02 07:00 22.7
2017/08/02 08:00 23.3
2017/08/02 09:00 23.9
2017/08/02 10:00 24.7
2017/08/02 11:00 24.5
2017/08/02 12:00 24.3
2017/08/02 13:00 24.7
2017/08/02 14:00 24.3
2017/08/02 15:00 24.8
2017/08/02 16:00 24.8
2017/08/02 17:00 24.4
2017/08/02 18:00 23.9
2017/08/02 19:00 23.4
2017/08/02 20:00 23.1
2017/08/02 21:00 22.7
2017/08/02 22:00 22.3
2017/08/02 23:00 21.7
2017/08/03 00:00 21.3
課題(2)
- 以下のようなグラフを作成せよ（テスト用として枠線のみとし、データのプロットは省略）。
- x軸は、annotationが月名（頭文字のみ）で、ticsが10日間隔としている。
解答
課題(1)
- 日付データの変換をExcelで行う場合、セルの表示形式を「yyyy-mm-ddThh:mm」とすればGMT形式となる。
- 以下の例では、awkで変換を行う場合のスクリプトを記述しておきます。
- awkの-Fオプションで区切り文字をtabに変更し、gsub関数で文字列を置換している。
plot_assign04-01.sh
#!/bin/bash
#--------------------------------------------------
# plot for GMT5
# 1次元時系列グラフ
#--------------------------------------------------
#----- 前処理・設定 ------------------------------
# GMTパラメーター
gmt gmtset FONT_ANNOT_PRIMARY=9p,Helvetica # 主目盛りのフォント
gmt gmtset FONT_LABEL=10p,Helvetica # 軸タイトルのフォント
# 日付関係
gmt gmtset FORMAT_DATE_MAP=mm/dd # 日付目盛り(yy[yy]=年,mm=月,dd=日,o=月名,jjj=日数)
gmt gmtset FORMAT_TIME_PRIMARY_MAP=a # 主目盛りの月・週表示：[f|a|c]
# input
in01=temp.txt # 入力ファイル
#in01_h=1 # ヘッダー行数
#in01_c=0,1 # 入力する列(0始まり）
# output
ps=plot_assign04-01.ps # 出力psファイル
# グラフのサイズ
mp=X # 地図の投影法(X:xy, M:メルカトル図法)
gx=10 # グラフの幅[cm]
gy=4 # グラフの高さ[cm]
#log_x=l # x軸を対数軸にするならl
#log_y=l # y軸を対数軸にするならl
jopt=$mp${gx}c${log_x}/${gy}c${log_y} # -Jオプション
# グラフの範囲(range)（-Rオプション）
rw=2017-08-01T00:00 # x軸の最小値(west)
re=2017-08-03T00:00 # x軸の最大値(east)
rs=0 # y軸の最小値(south)
rn=40 # y軸の最大値(north)
# 軸・目盛りの設定
ti_x="Date (2017)" # x軸タイトル
ti_y="Temperature [@~\260@~C]" # y軸タイトル
tics_x=a1Df6h # x軸の目盛り(a=主,f=副)
grid_x=g0 # x軸のグリッド(g)(*)
tics_y=a10f2 # y軸の目盛り(a=主,f=副)
grid_y=g10 # y軸のグリッド(g)(*)
axis=WeSn # 軸表示（大文字：軸＋目盛り、小文字：軸）
# 線の設定
line01=0.8p,black # 線：太さ,色,線種
#----- plot ------------------------------
echo "plot $ps"
# basemap(グリッド線のみ)
gmt psbasemap -J$jopt -R$rw/$re/$rs/$rn -Bx$grid_x -By$grid_y -B+n -P -K > $ps
# psxy(line)
awk -F"\t" 'NR>1 {gsub("/", "-", $1); gsub(" ","T",$1); print $1, $2}' $in01 | \
gmt psxy -J -R -W$line01 -K -O >> $ps
# basemap(枠線)
gmt psbasemap -J -R -Bx$tics_x+l"$ti_x" -By$tics_y+l"$ti_y" -B$axis -O >> $ps
#----- 後処理 ------------------------------
# 画像変換(psconvert)の設定
pc_fmt=g # g=png, b=bmp, j=jpg, t=tiff, e=eps, f=PDF, s=SVG
pc_res=300 # 画像の解像度(DPI)
pc_mgn=0.1c # グラフ外側のマージン
pc_aa=2 # アンチエイリアシングの設定(1,2,4)
# png変換
gmt psconvert -A$pc_mgn -E$pc_res -T$pc_fmt -Qt$pc_aa -Qg$pc_aa $ps
# 中間ファイルの削除
rm -f gmt.history gmt.conf
課題(2)
- psabasemapの-Bオプションに加え、FORMAT_DATE_MAP, FORMAT_TIME_PRIMARY_MAPの設定値を変更する。
plot_assign04-02.sh
#!/bin/bash
#--------------------------------------------------
# plot for GMT5
# 1次元時系列グラフ
#--------------------------------------------------
#----- 前処理・設定 ------------------------------
# GMTパラメーター
gmt gmtset FONT_ANNOT_PRIMARY=9p,Helvetica # 主目盛りのフォント
gmt gmtset FONT_ANNOT_SECONDARY=9p,Helvetica # 副目盛りのフォント
gmt gmtset FONT_LABEL=10p,Helvetica # 軸タイトルのフォント
# 日付関係
gmt gmtset FORMAT_DATE_MAP=o # 日付目盛り(yy[yy]=年,mm=月,dd=日,o=月名,jjj=日数)
gmt gmtset FORMAT_TIME_PRIMARY_MAP=c # 主目盛りの月・週表示：[f|a|c]
# input
#in01=temp.txt # 入力ファイル
#in01_h=1 # ヘッダー行数
#in01_c=0,1 # 入力する列(0始まり）
# output
ps=plot_assign04-02.ps # 出力psファイル
# グラフのサイズ
mp=X # 地図の投影法(X:xy, M:メルカトル図法)
gx=10 # グラフの幅[cm]
gy=4 # グラフの高さ[cm]
#log_x=l # x軸を対数軸にするならl
#log_y=l # y軸を対数軸にするならl
jopt=$mp${gx}c${log_x}/${gy}c${log_y} # -Jオプション
# グラフの範囲(range)（-Rオプション）
rw=2017-01-01T00:00 # x軸の最小値(west)
re=2018-01-01T00:00 # x軸の最大値(east)
rs=0 # y軸の最小値(south)
rn=40 # y軸の最大値(north)
# 軸・目盛りの設定
ti_x="Date (2017)" # x軸タイトル
ti_y="Temperature [@~\260@~C]" # y軸タイトル
tics_x=a1Of10d # x軸の目盛り(a=主,f=副)
grid_x=g0 # x軸のグリッド(g)(*)
tics_y=a10f2 # y軸の目盛り(a=主,f=副)
grid_y=g10 # y軸のグリッド(g)(*)
axis=WeSn # 軸表示（大文字：軸＋目盛り、小文字：軸）
# 線の設定
line01=0.8p,black # 線：太さ,色,線種
#----- plot ------------------------------
echo "plot $ps"
# basemap(グリッド線のみ)
gmt psbasemap -J$jopt -R$rw/$re/$rs/$rn -Bx$grid_x -By$grid_y -B+n -P -K > $ps
# basemap(枠線)
gmt psbasemap -J -R -Bx$tics_x+l"$ti_x" -By$tics_y+l"$ti_y" -B$axis -O >> $ps
#----- 後処理 ------------------------------
# 画像変換(psconvert)の設定
pc_fmt=g # g=png, b=bmp, j=jpg, t=tiff, e=eps, f=PDF, s=SVG
pc_res=300 # 画像の解像度(DPI)
pc_mgn=0.1c # グラフ外側のマージン
pc_aa=2 # アンチエイリアシングの設定(1,2,4)
# png変換
gmt psconvert -A$pc_mgn -E$pc_res -T$pc_fmt -Qt$pc_aa -Qg$pc_aa $ps
# 中間ファイルの削除
rm -f gmt.history gmt.conf