2020/03/15
2018/11/09
直交座標系での計算
JAGURSで直交座標系(CARTESIAN)で計算するには,Makefileで
#CARTESIAN=ON
を有効にするために#のコメントアウトを外します.
従って,球座標(JAGURSのデフォルト)で作成したバイナリとは別になりますので時々で使い分ける場合はファイル名を変えるなり,ディレクトリを別にするなりしておく必要があります.
OpenTSUNAMI (杞憂プロジェクト)はぢめます.詳しくは後日.
#CARTESIAN=ON
を有効にするために#のコメントアウトを外します.
従って,球座標(JAGURSのデフォルト)で作成したバイナリとは別になりますので時々で使い分ける場合はファイル名を変えるなり,ディレクトリを別にするなりしておく必要があります.
OpenTSUNAMI (杞憂プロジェクト)はぢめます.詳しくは後日.
Open TSUNAMI(オープン津波)とか称してオープンソースソフトであるJAGURSとオープンデータ(広義)を使って色々やってみることにしました. https://t.co/hsjagQXTkg— 山本直孝 (@naotakayamamoto) 2017年12月21日
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.2 で read して, 10f9.2 で write します.さらに,アスキーフォーマットとして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形式に変換します.
日本海講側は,中央防災会議のものが公開されています(古いです).
元データは内閣府の南海トラフの巨大地震モデル検討会で作成された地形データが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.2 で read して, 10f9.2 で write します.さらに,アスキーフォーマットとして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
これで,ネスティングしない場合は計算できるようになりました.ネスティングする場合は,メッシュの細かい調整が必要なので後日掲載します(一応,手元では動くデータは作れましたがこれで良いか確認します).日本海講側は,中央防災会議のものが公開されています(古いです).
登録:
投稿 (Atom)