2018/01/09

地形モデル

遠地津波計算をするためにGEBCOのデータを変換する手順は以前書きましたので,今回は近地津波の計算をするための直交座標系での地形データの作成手順を書いておきます.

元データは内閣府の南海トラフの巨大地震モデル検討会で作成された地形データがG空間情報センターからダウンロード出来ます.無償にはなっていますが,ダウンロードするためにはアカウントを作成する必要があります.

例として,JAGURS で使える形式に変換してみます.
手順は大きく分けで3つ.(1)まずはデータをダウンロードして,(2)まともなGISフォーマットに変換して,(3)最後に必要なフォーマットへの変換をします.
地形データは,ここからダウンロードできます.利用規約はもちろん確認が必要ですが,「計算入力データ(地形,粗度,堤防データ)の説明」を確認して下さい.ここに何系がどの領域に対応しているか書かれています.例えば四国を計算する場合は第04系になるのでその地形データをダウンロードします.ダウンロードしたデータにはメタ情報が何も書かれていないので,別途「計算範囲設定」を確認する必要があります.下記の ncols, nrows, xllcorner, yllcorner の値は,エクセルのメッシュ個数(X方向,Y方向),南西端の位置(X座標,Y座標)に対応します.

データは 10f8.2 で記述してあるということなので,Fortranで読み込むのは簡単ですが, f8.2 の場合,標高(陸域)が1000mを超える場合に,符号も含めると -1000.00 となり8桁になってしまうため,前の数字とスペースなしで連結されてしまいます.Fortranなら読めるでしょって言ってしまえばそれまでなのですが,一般的なアスキーフォーマットと考えた場合,スペースで区切られていないと使えません.というわけで, 10f8.2read して, 10f9.2write します.さらに,アスキーフォーマットとしてGISデータになるようにヘッダを付与します.プログラム(reformat.f)も掲載しておきます.
     ncols 720
     nrows 540
     xllcorner -730200
     yllcorner -857700
     cellsize 2430
     NODATA_value -9999

ヘッダの中身は,さほど難しくなく,ncolsはX方向のメッシュ数,nrowsはY方向のメッシュ数,xllcorner は南西端のX座標,yllcorner は南西端のY座標,cellsize は解像度(上記例は 2430m),NODATA_value は不要ですがお決まりの -9999 を設定しています.

reformat.f:
      PROGRAM reformat
      real*8 data(10)
      character(len=128) :: fname
      read(*,'(a)') fname
      open(10,file=fname)
      do
         read(10,'(10f8.2)',end=100) data
         write(*,'(10f9.2)') data
      end do
 100  END

作成したアスキーファイル(bathy.SD01.asc)をJAGURS用のcf形式に変換します.
% grdreformat bathy.SD01.asc=gd bathy.SD01.grd=cf
これで,ネスティングしない場合は計算できるようになりました.ネスティングする場合は,メッシュの細かい調整が必要なので後日掲載します(一応,手元では動くデータは作れましたがこれで良いか確認します).

2017/12/20

JAGURS on Mac (High Sierra)

JAGURSをMacで使う手順は,JAGURS on Mac でまとめましたが,手元がHigh Sierraになったので仕切り直し.

ライブラリ関係は,MacPorts で入れてます.
  • NetCDF
    • netcdf @4.4.1.1_1+dap+netcdf4 (active)
    • sudo port install netcdf +dap +netcdf4
  • netcdf-fortran @4.4.4_1+gcc7 (active)
    • sudo port install netcdf-fortran +gcc7
  • PROJ4
    • proj @4.9.3_0 (active)
    • sudo port install proj
  • FFTW3
    • fftw-3 @3.3.5_1+gfortran (active)

netcdf-fortran の variants を +gcc7 にしました.

HPC for Mac OS X は,High Sierra 版がないので,Sierra 対応な 7.1 を入れました.大事なことを書き忘れていた気がします.JAGURSは,V0407 です.Makefile前回の記事から変更していません.


2017/12/14

Mac OS X の mail コマンドでメールを送る

外向けの 25 番ポートが塞がれているというのはよくある話かと思います.その場合,Postfix の設定ファイルである /private/etc/postfix/main.cfrelayhost を設定してあげれば動くようになります.

relayhost = smtp.example.org

確認は, /var/log/mail.log で出来ます.

あとは,,,認証が入っていることが殆どだと思いますので,その場合は,
relayhost = smtp.gmail.com:587

# GMail SMTP confguration
smtp_sasl_auth_enable = yes
smtp_sasl_password_maps = hash:/etc/postfix/sasl_passwd
smtp_sasl_security_options = noanonymous
smtp_sasl_tls_security_options = noanonymous
smtp_sasl_mechanism_filter = plain
smtp_use_tls = yes
のように設定することが出来ます(GMailの例).

パスワードを設定している /etc/postfix/sasl_passwd は,
smtp.gmail.com:587 username:password
のように記述します.

ほぼ自分用のメモでした...

2017/12/12

Emacs利用のためのiTerm2設定

iTerm2内のEmacsでMewを快適に使うためのキーバインディングを書いてましたが,自分でも新しいマシンに同じ設定をするのが面倒なのでまとめ.

⌘-d ⇒ M-d
⌘-x ⇒ M-x
⌘-< ⇒ M-<
⌘-> ⇒ M->
⌘-u ⇒ M-u
⌘-w ⇒ M-w
⌘-q ⇒ M-q
⌘-% ⇒ M-%

(副作用)
⌘-w ・・・ Windowのcloseが出来なくなります.が,ターミナルの場合は,exit (C-d)でclose出来るので弊害はない.
⌘-q ・・・ iTerm2の終了が出来なくなりますが,exitで各windowをcloseしてあれば,⌘-qがiTerm2の終了として動作します.

(禁忌)
⌘-v ・・・ こいつはペーストなのでやらない
⌘-$ ・・・ スクリーンショットの呼び出しと衝突する

具体的な手順は,iTerm2内のEmacsでMewを使うに書いてあります.

(修正履歴)
2018/02/01 ⌘-q を忘れていたので追記

2017/12/06

JAGURSのデータをQGISで表示する

JAGURSで扱うデータは,NetCDFのちょっと古いバージョンですが標準フォーマットなので,GMTだけではなくQGISでも表示できるはずです.ということで,やってみました.

まずは,JAGURS用に変換した標高データ(gebco2min.grd)です.レイヤーの追加で,座標参照システムとして,EPSG:4326(WGS84な緯度経度)を選択します.


それだけで,モノクロで表示されます.ただし,水深が正の値なので海が白くなってしまいます.


色を付けるには,プロパティ→スタイルで,レンダータイプを「単バンド疑似カラー」にします.色は,RdYIBuにしてみました.負の値が赤に割り当てられていますが,今は水深が正の値なので反転はさせないでおきます.


可視化結果が以下の通り.色は,,,これが良いと言うわけではないので(お薦めを見つけたら書き直します)自由に選んで下さい.


ここに,計算結果のSD01.zmax.grdをオーバーレイします.レイヤーの追加で同じくEPSG:4326を指定すれば読み込み完了です.このデータは正のデータだけとりあえず見たいので,0から1の間で白から赤にします.



これだと,既に読み込んである地形データが白で上塗りされてしまってもったいない状態です.プロパティ→透過性で,「画面から値を追加」で陸域を透過設定します.また,追加のデータなし値に0を追加して波が届かなかった場所を背景の地形にします.


結局,最大波高の色をYIOrRdにしています.

見やすいとは言い難いですが,以下のようになりました.