GMTで北海道海岸線+等雨量線図を描いてみる

環境研究室 更新 20051114

北海道の白図ができたら、等雨量線図にもチャレンジしてみましょう。まず、以下1〜4の準備をしてください。

参考画像
  1. 作業フォルダを設定しよう。
    ここでは\GMTの下に、\gmtsvript、\topodata、\raindatフォルダを作ってください。別に名前はなんでもいいのですが、あとで分かりやすくするためこうしています。
  2. \gmtscript フォルダに下記スクリプトをcontour.bat の形でファイルしてください。(contour.batはここからダウンロード可能)また、観測地点毎の総雨量を面的に補完するため、surfaceコマンドで作ったグリッドデータtempr19920808-09.grd(約20MBあります)もここのフォルダに置いておきます。
  3. \topodata フォルダにhokkaido.riv を置きます。(hokkaido.rivはここからダウンロード可能)※hokkaido.rivの作り方はここをクリックしてください。
  4. \raindata フォルダに全道アメダスデータo920808.csvを置きます。(ここからダウンロード可能)。このデータには1992年8月8日〜9日の全道アメダス観測所の観測所コード、緯度、経度と総雨量の値が入っています。ここで用いている全道アメダスのデータには、欠測のある観測所のデータは除いています。アメダス観測所の緯度経度などの諸量は以下サイトからも入手できます。
    http://www.data.kishou.go.jp/meta/amedas/ame_master.pdf
  5. 2のフォルダ\gmtscriptに戻り、 contour.bat をダブルクリックして実行すると、con_mask_hok19920808-09.ps というポストスクリプト形式の北海道海岸線+等雨量線図が出来上がっています(いるはず!)。
    これを見るには、イラストレータかgsviewというビューワが必要です。上の様な絵ができているのを見ることができるはずです。(結果のポストスクリプトファイルはここからダウンロードできます)
  6. タイトルや河川名などを日本語で修正する場合、GMTで作ったポストスクリプトファイルをイラストレータで開き、編集したほうが便利かもしれません。(20051114)
  7. スクリプトの説明 基本的な文法はWINDOWSのコマンドスクリプトです。 REMは注釈文です。 set文は、ある変数になにかの値を代入する命令。数字でも文字列でもいいようです。 set文で定義された変数は、別のコマンド中で%記号に挟まれて使用できます。これは便利! gawk文 そもそもUNIX系のコマンドawk(オークと読むらしい)をwindowsで使えるようにしたコマンド。ここではファイル中のデータ列の何番目と何番目を抽出する、といった機能しか使ってませんが、他にも機能があるようです。

現段階ではこんないい加減な説明しかできず、心苦しい限りです。
あとは変数をいろいろいじってみてどうなるか、各種サイトで公開されている用例を見ながら、実践あるのみです!
上図で使用したデータは環境研究室で平成16年度に整理したデータをベースにしています。アメダスのデータであればWEB上で公開されていますので、自分でもデータを入手して等雨量線図を作ることができますよね?

北海道海岸線+等雨量線図用スクリプト

本スクリプトは\gmtscript\フォルダに入っているものとする。

フォルダの位置関係(例)
c:\gmt\gmtscript\  ここに 920808-09.grd,con_920808-09.ps,GMT_rainbow.cpt 等がファイル
c:\gmt \topodata\   ここに hokkaido.riv がファイルされている
c:\gmt\raindat\    ここに o920808.csv がファイルされている

REM ----------------------------------------------
REM 北海道の海岸線に等雨量線をプロットするスクリプト 2008.4.23
REM hokriv:道内河川データ、raindat:H04.8豪雨の観測所別総雨量
REM ..はフォルダの相対位置関係を示し、1個で自己フォルダの下位、2個で同じ階層のフォルダを参照している
set hokriv=..\topodata\hokkaido.riv
REM ↓読み込みファイル名は8文字以内じゃないとうまく読み込んでくれないようです。
set raindat=..\raindata\o920808.csv
set grddata=920808-09.grd
REM of:最終的な出力イメージ(ポストスクリプト)ファイル名を格納する変数。このpsファイルにイメージが重ねられていく。
set of=con_920808-09.ps
REM cp:カラーパレット
set cp=GMT_saru_19920808-09.cpt
REM 色パレットの変換。0mmから250?まで50?ピッチで色を割り当てている。-zで連続的な色の変化を指定。
makecpt -CGMT_rainbow.cpt -T0/250/50 -Z > %cp%
REM
REM set waku=141.66666/143.33334/42.33333/43.33334 :この枠は沙流川流域用
set waku=139/146/41.0/46.0
REM ↑wakuは描画の範囲(緯度経度)を示す。
set bp=a2.0f1.0g0.0
REM ↑a:軸目盛数字の間隔(度)、b:白黒帯の間隔(度)、g:線の間隔(度)。分のときは*.*m、秒のときは*.*c
set r=2
REM r:1度あたり何cmで書くかという意味。ここでは1度あたり2cm
REM ↓基本の地図の枠を作る。この枠のイメージをベースに図が追加されていく。
psbasemap -Jm%r% -R%waku% -B%bp%/%bp% -P -K > %of%
REM −Oは追加のとき必須。>>は追加書き。
REM ----------------------------------------------
REM ↓雨量地点の緯度、経度、総雨量(?)をgawk文で並び替えし、rdat920808.txtに変換
REM ↓surfaceコマンドで面的に補間し、グリッドデータ化。このグリッドデータが等雨量線図に必要。
REM 一旦グリッドデータを作ってしまえば、以下のコマンド2行をREM文にして実行しない。(コメントアウト)
REM ----------------------------------------------
rem gawk "{print $3 , $2 , $4 }" %raindat% > rdat920808.txt
rem surface rdat920808.txt -G%grddata% -I11.25c/7.5c -R%waku% -T0.25 -Ll0.0
REM ----------------------------------------------
REM ※surfaceコマンドはランダム配列のxyzデータをグリッドに変換するGMTコマンド
REM ----------------------------------------------
REM 等雨量線図を線でなく面的に描く場合はこのコマンドを実行。
REM grdimage %grddata% -Jm%r% -R%waku% -B%bp% -C%cp% -O -K >> %of%
REM
REM ↓雨量観測点のプロット( gawkが重要な仕事をする)
REM ----------------------------------------------
gawk "{print $3,$2}" %raindat% | psxy -R%waku% -Jm%r% -Sc0.1 -G0 -O -K >> %of%
REM ↑ここでは%raindat%の3列目$3と2列目のデータ$2をパイプ"|"コマンドでpsxyに渡している。
REM なぜ3列目と2列目かというと、このデータは3列目が経度(X軸)、2列目が緯度(Y軸)であるから。
REM ----------------------------------------------
REM ↓全道河川(%hokriv%)を描く。ここでもgawkとパイプとpsxyが活躍!
REM ----------------------------------------------
gawk "{print $1,$2}" %hokriv% | psxy -Jm%r%ハ -W2/0/0/255 -M -R -P -V -O -K >> %of%
REM ----------------------------------------------
REM 等雨量線を描く。Surfaceコマンドでグリッド化した雨データを使っている。
grdcontour %grddata% -Jm%r% -R%waku%ハ -Wa4/255/0/ -Wc2/255/0/0 -C10 -A100 -L10/500 -O -K>> %of%
REM --------------------------------------4は線の太さ、255/0/0は色(ここでは赤)
REM Waは数値付等高線、Wcは数値無し等高線、−A100は数値の間隔
REM 海岸線の描画 Sコマンドで指定しているのは水域のベタ塗りの色で、ここでは白(255)を指定。
REM これを指定しなければ海の上にも等雨量線図がうねうね出てしまい、シマラない絵になってしまう。
pscoast -Jm%r% -R -B -W2 -Dh -P -V -S255 -O -K >> %of%
REM

以上、ランダムに配置されている雨量観測所のデータをグリッドデータにするあたりが面倒ですが、それさえクリアすれば、それほど面倒ではないと思います。各コマンドで定義している変数の意味はほかにもっと丁寧な説明がされているサイトがありますので、調べて見てください。もしかするともっとエレガントな使い方があるかも。

文責:水環境保全チーム 村上