1 ! 2 ! バイナリファイルから netCDF データを生成するプログラム 3 ! 4 ! バイナリ : /global/pr_preci/2000/monthly/gpr00100.zip 5 ! 6 ! 方針 : 格子点, 次元の取り方等, 構造はそのまま 7 ! netCDF コンベンションは Gtool4 - NetCDF Convention に 8 ! 従うものとする. 9 ! 2002/10/19 塚原 大輔 作成 10 ! 11 12 program test_netcdf 13 14 ! モジュールの呼び出し 15 use netcdf !netcdfモジュール呼び出し 16 17 implicit none 18 19 ! 元データ(TRMM/GANAL cdrom)から変数を取り出すための 20 ! 容れ物変数の宣言 21 ! 22 23 character(len = 80) :: fname !ファイル名(ex.gpr00100) 24 25 !ファイルに格納される各変数 26 27 character(len = 80) :: header(7) !ファイルの簡単な情報が書かれている変数. 28 29 integer(2) :: pr(720,153) !地上における層状性, 対流性降雨量の合計 30 integer(2) :: pr_2(720,153) !高度2kmにおける層状性, 対流性降雨量の合計 31 integer(2) :: pr_4(720,153) !高度4kmにおける層状性. 対流性降雨量の合計 32 integer(2) :: pr_6(720,153) !高度6kmにおける層状性, 対流性降雨量の合計 33 34 integer(2) :: cpr_2(720,153) !高度2kmにおける対流性降雨量 35 integer(2) :: cpr_4(720,153) !高度4kmにおける対流性降雨量 36 integer(2) :: cpr_6(720,153) !高度6kmにおける対流性降雨量 37 38 integer(2) :: spr_2(720,153) !高度2kmにおける層状性降雨量 39 integer(2) :: spr_4(720,153) !高度4kmにおける層状性降雨量 40 integer(2) :: spr_6(720,153) !高度6kmにおける層状性降雨量 41 42 ! netcdf 出力用変数の宣言 43 44 integer :: ncid, i 45 integer :: dim_lonid, dim_latid, dim_altid !次元ID. 左から緯度, 経度, 高度. 46 integer :: var_lonid, var_latid, var_altid !変数ID. 同上. 47 integer :: var_prid, var_cprid, var_sprid !降雨量. 左から合計, 対流性, 層状性. 48 integer :: stat !エラー判断 49 integer :: prmid(3) !降雨量の次元変数群. 50 real, dimension(720,153) :: & 51 rhValues = 0 !全ての格子において値が0である配列. 52 !対流性, 層状性降雨量が共に地上での値がないので 53 !高度の次元をそろえるために用意. 54 !データ読み取り 55 fname='gpr00100' 56 open(80,file=fname,form='unformatted',status='old') 57 read(80) header 58 read(80) pr 59 read(80) pr_2 60 read(80) pr_4 61 read(80) pr_6 62 63 read(80) cpr_2 64 read(80) cpr_4 65 read(80) cpr_6 66 67 read(80) spr_2 68 read(80) spr_4 69 read(80) spr_6 70 close(80) 71 ! netCDFファイルの生成 命令文 72 stat = nf90_create('gpr00100.nc',nf90_noclobber,ncid) 73 if(stat /= nf90_noerr) call handle_err(stat) 74 75 ! 次元の定義 76 77 stat = nf90_def_dim(ncid, 'lon', 720,dim_lonid) !経度の定義 78 if(stat /= nf90_noerr) call handle_err(stat) 79 80 stat = nf90_def_dim(ncid, 'lat', 153, dim_latid) !緯度の定義 81 if(stat /= nf90_noerr) call handle_err(stat) 82 83 stat = nf90_def_dim(ncid,'alt',4,dim_altid) !鉛直方向ベクトルの定義 84 if(stat /= nf90_noerr) call handle_err(stat) 85 86 ! 変数の定義 87 88 stat = nf90_def_var(ncid,'lon',nf90_double,& 89 &dim_lonid,var_lonid) !緯度座標変数の定義 90 if(stat /= nf90_noerr) call handle_err(stat) 91 92 93 stat = nf90_def_var(ncid,'lat',nf90_double,& 94 &dim_latid,var_latid) !緯度座標変数の定義 95 if(stat /= nf90_noerr) call handle_err(stat) 96 97 98 stat = nf90_def_var(ncid,'dim',nf90_double,& 99 &dim_altid,var_altid) !鉛直座標変数の定義 100 if(stat /= nf90_noerr) call handle_err(stat) 101 102 prmid = (/dim_lonid,dim_latid,dim_altid/) ! 降雨量の次元配列 103 104 stat = nf90_def_var(ncid,'pr',nf90_double,& 105 &prmid,var_prid) ! 降雨強度の定義 106 if(stat /= nf90_noerr) call handle_err(stat) 107 108 stat = nf90_def_var(ncid,'pr_c',nf90_double,& 109 &prmid,var_cprid) ! 降雨強度の定義 110 if(stat /= nf90_noerr) call handle_err(stat) 111 112 stat = nf90_def_var(ncid,'pr_s',nf90_double,& 113 &prmid,var_sprid) ! 降雨強度の定義 114 if(stat /= nf90_noerr) call handle_err(stat) 115 116 ! 属性定義 117 !! 大域属性 118 119 stat = nf90_put_att(ncid,nf90_global,'Conventions',& 120 &'http://www.gfd-dennou.org/arch/gtool4/conventions/')! お決まり 121 if(stat /= nf90_noerr) call handle_err(stat) 122 123 124 stat = nf90_put_att(ncid, nf90_global,'title',& 125 & 'TRMM/PR Precipitation ( Monthly, Jan. 2000 ), Global ') ! タイトル 126 if(stat /= nf90_noerr) call handle_err(stat) 127 128 stat = nf90_put_att(ncid,nf90_global,'source',& 129 &'For study of global climate change::TRMM Rain Data Set')!出典 130 if(stat /= nf90_noerr) call handle_err(stat) 131 132 stat = nf90_put_att(ncid,nf90_global,'history',& 133 &'2002-11-27T04:28:20+0900 tsukahara')!履歴 134 if(stat /= nf90_noerr) call handle_err(stat) 135 136 !! 各変数の属性@単位 137 stat = nf90_put_att(ncid,var_lonid,'units','degrees_east')!緯度 138 if(stat /= nf90_noerr) call handle_err(stat) 139 stat = nf90_put_att(ncid,var_latid,'units','degrees_north')!経度 140 if(stat /= nf90_noerr) call handle_err(stat) 141 stat = nf90_put_att(ncid,var_altid,'units','km') 142 if(stat /= nf90_noerr) call handle_err(stat) 143 stat = nf90_put_att(ncid,var_prid,'units','mm') 144 if(stat /= nf90_noerr) call handle_err(stat) 145 stat = nf90_put_att(ncid,var_cprid,'units','mm') 146 if(stat /= nf90_noerr) call handle_err(stat) 147 stat = nf90_put_att(ncid,var_sprid,'units','mm') 148 if(stat /= nf90_noerr) call handle_err(stat) 149 150 151 !! 各変数の属性@longname 152 153 stat = nf90_put_att(ncid,var_lonid,'long_name','longitude') 154 if(stat /= nf90_noerr) call handle_err(stat) 155 stat = nf90_put_att(ncid,var_lonid,'modulo','360.0') 156 if(stat /= nf90_noerr) call handle_err(stat) 157 stat = nf90_put_att(ncid,var_latid,'long_name','latitude') 158 if(stat /= nf90_noerr) call handle_err(stat) 159 stat = nf90_put_att(ncid,var_altid,'long_name','altitude') 160 if(stat /= nf90_noerr) call handle_err(stat) 161 162 stat = nf90_put_att(ncid,var_prid,'long_name',& 163 & ' precipitaion ') 164 if(stat /= nf90_noerr) call handle_err(stat) 165 stat = nf90_put_att(ncid,var_cprid,'long_name',& 166 &'convective precipitation ') 167 if(stat /= nf90_noerr) call handle_err(stat) 168 stat = nf90_put_att(ncid,var_sprid,'long_name',& 169 & 'stratiform precipitation') 170 if(stat /= nf90_noerr) call handle_err(stat) 171 172 173 ! 定義モードから抜ける 174 175 stat = nf90_enddef(ncid) 176 177 if(stat /= nf90_noerr) call handle_err(stat) 178 179 180 ! 変数にデータを代入する. 181 182 stat = nf90_put_var(ncid,var_lonid,(/((i-1)*0.5,i=1,720)/) )!緯度方向の座標値 183 if(stat /= nf90_noerr) call handle_err(stat) 184 stat = nf90_put_var(ncid,var_latid,(/(((i-1)*0.5-38),i=1,153)/) )!経度方向の座標値 185 if(stat /= nf90_noerr) call handle_err(stat) 186 stat = nf90_put_var(ncid,var_altid,(/(i,i=0,6,2)/) ) !高度方向の座標値 187 if(stat /= nf90_noerr) call handle_err(stat) 188 189 ! 通常降雨 190 191 stat = nf90_put_var(ncid,var_prid,pr,& 192 & start=(/1,1,1/)) !通常降雨 地表面 193 if(stat /= nf90_noerr) call handle_err(stat) 194 stat = nf90_put_var(ncid,var_prid,pr_2,& 195 & start=(/1,1,2/)) !通常降雨 地上から2km 196 if(stat /= nf90_noerr) call handle_err(stat) 197 stat = nf90_put_var(ncid,var_prid,pr_4,& 198 & start=(/1,1,3/)) !通常降雨 地上から4km 199 if(stat /= nf90_noerr) call handle_err(stat) 200 stat = nf90_put_var(ncid,var_prid,pr_6,& 201 & start=(/1,1,4/)) !通常降雨 地上から6km 202 if(stat /= nf90_noerr) call handle_err(stat) 203 204 ! 対流性降雨 205 ! 206 stat = nf90_put_var(ncid,var_cprid,rhValues,& 207 & start=(/1,1,1/)) !対流性降雨 地表面はなし 208 if(stat /= nf90_noerr) call handle_err(stat) 209 stat = nf90_put_var(ncid,var_cprid,cpr_2,& 210 & start=(/1,1,2/)) !対流性降雨 地上から2km 211 if(stat /= nf90_noerr) call handle_err(stat) 212 stat = nf90_put_var(ncid,var_cprid,cpr_4,& 213 & start=(/1,1,3/)) !対流性降雨 地上から4km 214 if(stat /= nf90_noerr) call handle_err(stat) 215 stat = nf90_put_var(ncid,var_cprid,cpr_6,& 216 & start=(/1,1,4/)) !対流性降雨 地上から6km 217 if(stat /= nf90_noerr) call handle_err(stat) 218 219 ! 層状性降雨 220 ! 221 stat = nf90_put_var(ncid,var_sprid,rhValues,& 222 & start=(/1,1,1/)) !対流性降雨 地表面はなし 223 if(stat /= nf90_noerr) call handle_err(stat) 224 stat = nf90_put_var(ncid,var_sprid,spr_2,& 225 & start=(/1,1,2/)) !層状降雨 地上から2km 226 if(stat /= nf90_noerr) call handle_err(stat) 227 stat = nf90_put_var(ncid,var_sprid,spr_4,& 228 & start=(/1,1,3/)) !層状降雨 地上から4km 229 if(stat /= nf90_noerr) call handle_err(stat) 230 stat = nf90_put_var(ncid,var_sprid,spr_6,& 231 & start=(/1,1,4/)) !層状降雨 地上から6km 232 if(stat /= nf90_noerr) call handle_err(stat) 233 234 235 236 ! 作成したnetCDFファイルを閉じる 237 238 stat = nf90_close(ncid) 239 if(stat /= nf90_noerr) call handle_err(stat) 240 241 242 ! エラー表示用のサブルーチン部分. 杉山氏作成. 243 244 contains 245 246 subroutine handle_err(status) 247 integer, intent (in) :: status 248 if(status /= nf90_noerr) then 249 print *, trim(nf90_strerror(status)) 250 stop "Stopped" 251 end if 252 end subroutine handle_err 253 254 end program test_netcdf