! ! バイナリファイルから netCDF データを生成するプログラム ! ! バイナリ : /global/pr_preci/2000/monthly/gpr00100.zip ! ! 方針 : 格子点, 次元の取り方等, 構造はそのまま ! netCDF コンベンションは Gtool4 - NetCDF Convention に ! 従うものとする. ! 2002/10/19 塚原 大輔 作成 ! program test_netcdf ! モジュールの呼び出し use netcdf !netcdfモジュール呼び出し implicit none ! 元データ(TRMM/GANAL cdrom)から変数を取り出すための ! 容れ物変数の宣言 ! character(len = 80) :: fname !ファイル名(ex.gpr00100) !ファイルに格納される各変数 character(len = 80) :: header(7) !ファイルの簡単な情報が書かれている変数. integer(2) :: pr(720,153) !地上における層状性, 対流性降雨量の合計 integer(2) :: pr_2(720,153) !高度2kmにおける層状性, 対流性降雨量の合計 integer(2) :: pr_4(720,153) !高度4kmにおける層状性. 対流性降雨量の合計 integer(2) :: pr_6(720,153) !高度6kmにおける層状性, 対流性降雨量の合計 integer(2) :: cpr_2(720,153) !高度2kmにおける対流性降雨量 integer(2) :: cpr_4(720,153) !高度4kmにおける対流性降雨量 integer(2) :: cpr_6(720,153) !高度6kmにおける対流性降雨量 integer(2) :: spr_2(720,153) !高度2kmにおける層状性降雨量 integer(2) :: spr_4(720,153) !高度4kmにおける層状性降雨量 integer(2) :: spr_6(720,153) !高度6kmにおける層状性降雨量 ! netcdf 出力用変数の宣言 integer :: ncid, i integer :: dim_lonid, dim_latid, dim_altid !次元ID. 左から緯度, 経度, 高度. integer :: var_lonid, var_latid, var_altid !変数ID. 同上. integer :: var_prid, var_cprid, var_sprid !降雨量. 左から合計, 対流性, 層状性. integer :: stat !エラー判断 integer :: prmid(3) !降雨量の次元IDベクトル. real, dimension(720,153) :: & rhValues = 0 !全ての格子において値が0である配列. !対流性, 層状性降雨量が共に地上での値がないので !高度の次元をそろえるために用意. !データ読み取り fname='gpr00100' open(80,file=fname,form='unformatted',status='old') read(80) header read(80) pr read(80) pr_2 read(80) pr_4 read(80) pr_6 read(80) cpr_2 read(80) cpr_4 read(80) cpr_6 read(80) spr_2 read(80) spr_4 read(80) spr_6 close(80) ! netCDFファイルの生成 命令文 stat = nf90_create('gpr00100.nc',nf90_noclobber,ncid) if(stat /= nf90_noerr) call handle_err(stat) ! 次元の定義 stat = nf90_def_dim(ncid, 'lon', 720,dim_lonid) !経度の定義 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_def_dim(ncid, 'lat', 153, dim_latid) !緯度の定義 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_def_dim(ncid,'alt',4,dim_altid) !鉛直方向ベクトルの定義 if(stat /= nf90_noerr) call handle_err(stat) ! 変数の定義 stat = nf90_def_var(ncid,'lon',nf90_double,& &dim_lonid,var_lonid) !緯度座標変数の定義 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_def_var(ncid,'lat',nf90_double,& &dim_latid,var_latid) !緯度座標変数の定義 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_def_var(ncid,'dim',nf90_double,& &dim_altid,var_altid) !鉛直座標変数の定義 if(stat /= nf90_noerr) call handle_err(stat) prmid = (/dim_lonid,dim_latid,dim_altid/) ! 降雨量の次元配列 stat = nf90_def_var(ncid,'pr',nf90_double,& &prmid,var_prid) ! 降雨強度の定義 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_def_var(ncid,'pr_c',nf90_double,& &prmid,var_cprid) ! 降雨強度の定義 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_def_var(ncid,'pr_s',nf90_double,& &prmid,var_sprid) ! 降雨強度の定義 if(stat /= nf90_noerr) call handle_err(stat) ! 属性定義 !! 大域属性 stat = nf90_put_att(ncid,nf90_global,'Conventions',& &'http://www.gfd-dennou.org/arch/gtool4/conventions/')! お決まり if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid, nf90_global,'title',& & 'TRMM/PR Precipitation ( Monthly, Jan. 2000 ), Global ') ! タイトル if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,nf90_global,'source',& &'For study of global climate change::TRMM Rain Data Set')!出典 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,nf90_global,'history',& &'2002-11-27T04:28:20+0900 tsukahara')!履歴 if(stat /= nf90_noerr) call handle_err(stat) !! 各変数の属性@単位 stat = nf90_put_att(ncid,var_lonid,'units','degrees_east')!緯度 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_latid,'units','degrees_north')!経度 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_altid,'units','km') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_prid,'units','mm') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_cprid,'units','mm') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_sprid,'units','mm') if(stat /= nf90_noerr) call handle_err(stat) !! 各変数の属性@longname stat = nf90_put_att(ncid,var_lonid,'long_name','longitude') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_lonid,'modulo','360.0') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_latid,'long_name','latitude') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_altid,'long_name','altitude') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_prid,'long_name',& & ' precipitaion ') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_cprid,'long_name',& &'convective precipitation ') if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_att(ncid,var_sprid,'long_name',& & 'stratiform precipitation') if(stat /= nf90_noerr) call handle_err(stat) ! 定義モードから抜ける stat = nf90_enddef(ncid) if(stat /= nf90_noerr) call handle_err(stat) ! 変数にデータを代入する. stat = nf90_put_var(ncid,var_lonid,(/((i-1)*0.5,i=1,720)/) )!緯度方向の座標値 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_latid,(/(((i-1)*0.5-38),i=1,153)/) )!経度方向の座標値 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_altid,(/(i,i=0,6,2)/) ) !高度方向の座標値 if(stat /= nf90_noerr) call handle_err(stat) ! 通常降雨 stat = nf90_put_var(ncid,var_prid,pr,& & start=(/1,1,1/)) !通常降雨 地表面 if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_prid,pr_2,& & start=(/1,1,2/)) !通常降雨 地上から2km if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_prid,pr_4,& & start=(/1,1,3/)) !通常降雨 地上から4km if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_prid,pr_6,& & start=(/1,1,4/)) !通常降雨 地上から6km if(stat /= nf90_noerr) call handle_err(stat) ! 対流性降雨 ! stat = nf90_put_var(ncid,var_cprid,rhValues,& & start=(/1,1,1/)) !対流性降雨 地表面はなし if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_cprid,cpr_2,& & start=(/1,1,2/)) !対流性降雨 地上から2km if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_cprid,cpr_4,& & start=(/1,1,3/)) !対流性降雨 地上から4km if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_cprid,cpr_6,& & start=(/1,1,4/)) !対流性降雨 地上から6km if(stat /= nf90_noerr) call handle_err(stat) ! 層状性降雨 ! stat = nf90_put_var(ncid,var_sprid,rhValues,& & start=(/1,1,1/)) !対流性降雨 地表面はなし if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_sprid,spr_2,& & start=(/1,1,2/)) !層状降雨 地上から2km if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_sprid,spr_4,& & start=(/1,1,3/)) !層状降雨 地上から4km if(stat /= nf90_noerr) call handle_err(stat) stat = nf90_put_var(ncid,var_sprid,spr_6,& & start=(/1,1,4/)) !層状降雨 地上から6km if(stat /= nf90_noerr) call handle_err(stat) ! 作成したnetCDFファイルを閉じる stat = nf90_close(ncid) if(stat /= nf90_noerr) call handle_err(stat) ! エラー表示用のサブルーチン部分. 杉山氏作成. contains subroutine handle_err(status) integer, intent (in) :: status if(status /= nf90_noerr) then print *, trim(nf90_strerror(status)) stop "Stopped" end if end subroutine handle_err end program test_netcdf