[更新情報]

▽Go to footer

GMT でヒストグラムを描く事例です. プロットするデータは,気象庁のホームページに掲載されている,前橋市の過去 114 年間の年最大日雨量です. 並び替えを自動で行ってくれ,プロットデータも出力できるので便利です.

実行用バッチファイルと解説

rem **************************
rem ヒストグラム描画
rem **************************
set inp_org=inp_PRB_DM.txt
set out_text=out_DM_HISTO.txt
set wd=10
set xlabel="Maximum daily rainfall (mm/day)"
set fig_out=fig_DM_HISTO.eps
rem
minmax %inp_org% -H1 > %out_text%
pshistogram %inp_org% -W%wd% -IO -H1 >> %out_text%
rem ***************************************************************************
rem * minmax %inp_org% -H1 > %out_text%                                       *
rem *     ファイル%inp_org%の内容より,ファイル名,データ数,最小/最大を      *
rem *     ファイル%out_text%に新規出力                                        *
rem *     -H1 により先頭1行読み飛ばし                                         *
rem *     minmax の出力形式は以下のとおり.                                   *
rem *     (ファイル名): N = (データ数)  <最小値/最大値>                       *
rem *     事例== inp_PRB_DM.txt: N = 114	<43.2/357.4>                  *
rem *                                                                         *
rem * pshistogram %inp_org% -W%wd% -IO -H1 >> %out_text%                      *
rem *     ファイル%inp_org%の内容を値域幅%wd%でx値,y値(頻度)を               *
rem *     ファイル%out_text%に追加出力                                        *
rem *     -W%wd% === 値域幅wd(頻度を集計するx軸幅)                            *
rem *     -IO  ===x値,y値(頻度)を出力                                        *
rem *     -H1  ===先頭1行読み飛ばし                                           *
rem *     WARNING: 269046071 multiple modes found という警告がでるが          *
rem *     ファイル出力は行われる.                                            *
rem *     -IOオプションでxy値の出力は行われるが-Iオプションのみでは同様の     *
rem *     警告が出て画面出力はするがファイル出力はしない                      *
rem *     このため,最大・最小をファイルに落とすにはminmaxコマンドを          *
rem *     使用する必要有.                                                    *
rem ***************************************************************************
rem
gawk "{if(NR==1){nd=$4;str0=$5}}END{sub(//,\"\",str0);split(str0,a,\"/\");print nd,a[1],a[2]}" %out_text% > _tempdat.txt
rem ***************************************************************************
rem * ファイル%out_text%の1行目(minmaxの出力)からデータ数,最小値,最大値を   *
rem * 抜き出しファイル_tempdat.txtに書き出す                                  *
rem ***************************************************************************
rem
gawk "{printf \"0.95 0.95 12 0 0 MR n=%%g\n\",$1}" _tempdat.txt > _txtdat.txt
gawk "{printf \"0.95 0.90 12 0 0 MR min=%%g\n\",$2}" _tempdat.txt >> _txtdat.txt
gawk "{printf \"0.95 0.85 12 0 0 MR max=%%g\n\",$3}" _tempdat.txt >> _txtdat.txt
rem ***************************************************************************
rem * ファイル_tempdat.txt(1行のファイル)からデータ数,最小値,最大値の描画用 *
rem * ファイル_txtdat.txtに書き出す                                           *
rem ***************************************************************************
rem
gmtset ANOT_FONT_SIZE 14
gmtset LABEL_FONT_SIZE 14
gmtset TICK_LENGTH -0.2c
pshistogram %inp_org% -JX15/10 -B:%xlabel%:/:"Frequency":WSen -H1 -W%wd% -L0.5p -G230 -Z0 -P -K > %fig_out%
pstext _txtdat.txt -R0/1/0/1 -J -N -P -O -K >> %fig_out%
gawk "BEGIN{FS=\",\"}{if(NR==1)printf \"0.2 10.2 12 0 0 BL %%s\n\",$0}" %inp_org% | pstext -R -J -N -P -O >> %fig_out%
rem ***************************************************************************
rem * pshistogram                                                             *
rem *     -JX15/10 === 描画領域幅15cm,高さ10cm                               *
rem *     -B:"Maximum daily rainfall (mm/day)":/:"Frequency":WSen             *
rem *            === 軸ラベルと数値を描画する軸を指定                         *
rem *     -H1    === 先頭1行読み飛ばし                                        *
rem *     -W%wd% === 値域幅wd(頻度を集計するx軸幅)                            *
rem *     -L0.5p === 指定した太さのペンで棒の輪郭を描く                       *
rem *                 [デフォルトでは輪郭なし] ここでは0.5ptを指定            *
rem *     -G200  === 領域をグレースケール230で塗りつぶし                      *
rem *     -Z0    === 縦軸の選択.デフォルト(指定無しor -Z0)は個数.           *
rem *                -Z1を指定すれば頻度%                                    *
rem *     -P     === 縦長紙を指定                                             *
rem *     -Rを指定していないので軸範囲は自動選定                              *
rem *                                                                         *
rem * pstext                                                                  *
rem *     -R0/1/0/1 === pshistogramでは軸範囲はGMTにまかせたので,            *
rem *                   データ数・最小値・最大値描画時は軸範囲を再設定する    *
rem *     最後に左上にコメントを追加                                          *
rem ***************************************************************************

set inp_org=
set out_text=
set wd=
set xlabel=
set fig_out=

del .gmt*
del _*

バッチファイルのダウンロードと画像確認

ファイル名概要
GMT実行用バッチファイル
入力データファイル(年最大日雨量)
minmaxとpshistogramで作成したプロットデータ(1行目はminmaxの出力)
描画事例(ヒストグラム)



inserted by FC2 system