# 熱, 水収支を計算します # # 履歴 # 2010/03/26 納多哲史物理量の名前を文字変数で与えられるようにした # 2009/10/02 納多哲史北半球と南半球で領域が重複していたバグを修正 # 2009/10/01 納多哲史とりあえず完成 # 2009/09/24 納多哲史作成 require "numru/ggraph" include NumRu ana_start_time = 1000 ana_end_time = 2000 name_SSR = "SSRA" name_SLR = "SLRA" name_OSR = "OSRA" name_OLR = "OLRA" name_Sens = "SensA" name_Evap = "EvapA" name_Rain = "Rain" # 各変数について # * 調べる範囲 # * ある点の量 # * 半球平均 (東西南北) # * 計算する量 # * 大気の熱収支 # * 大気の水収支 # 解析対象となるファイルの場所 prefix = '.' # debug #prefix = '/home/noda/work/dcpam/dcpam5/practice/090629_test_reprod-morikawa_SR090_1Omega_nonSR' #prefix = '/home/noda/work/dcpam/dcpam5/practice/090605_test_reprod-morikawa_SR090_0.5Omega' def help() print < 0 case ARGV[0] when "--range_time" ARGV.shift ana_start_time, ana_end_time = get_range(ARGV[0]) when "--prefix" ARGV.shift prefix = ARGV[0] when "-h", "--help" help exit 0 else puts "ERROR: unknown argument!\n" help exit 20 end ARGV.shift end # HTML のヘッダを出力 print < Budget

Budget

Time range: #{ana_start_time} - #{ana_end_time} (day)
Unit is W m-2.
Upward is positive.
(Heat Budget) = (OSR + OLR) + (SSR + SLR) + Sens + Rain
(Water Budget) = EVap - Rain
(Surface Budget) = SSR + SLR + Sens + Evap
(Top Budget) = OSR + OLR

EOM puts "" puts "" for varname in var_list puts "" end puts "" puts "" puts "" puts "" puts "" path = File.join(prefix, var_list[0] + ".nc") lon_na = GPhys::IO.open(path, 'lon').data.val lat_na = GPhys::IO.open(path, 'lat').data.val #lon_weight_na = GPhys::IO.open(path, 'lon_weight').data.val lat_weight_na = GPhys::IO.open(path, 'lat_weight').data.val lat = GPhys::IO.open(path, 'lat') lat_weight = GPhys::IO.open(path, 'lat_weight') nlon = lon_na.length nlat = lat_na.length ntime = GPhys::IO.open(path, 'time').cut("time"=>ana_start_time..ana_end_time).val.length # 範囲の平均を取る時の場所を [経度開始, 経度終了, 緯度開始, 緯度終了] で指定. # 注意: Gphys は指定された値がないときは近い値を取ってくるので, # 半球平均などの際のダブルカウントに注意. ranges = [ [ 0, lon_na[nlon-1] , -90, 90 ], # 全球 [ 0, lon_na[nlon-1] , lat_na[nlat/2], 90 ], # 北半球 [ 0, lon_na[nlon-1] , -90, lat_na[nlat/2-1] ], # 南半球 [ 0, lon_na[nlon/2-1], -90, 90 ], # 東半球 [ lon_na[nlon/2], lon_na[nlon-1] , -90, 90 ] # 西半球 ] var = Hash.new # ファイルの読み込み for varname in var_list path = File.join(prefix, varname + ".nc") # GPhys クラスの変数をハッシュに格納 var[varname] = GPhys::IO.open( File.join(prefix, varname + ".nc"), varname).cut("time"=>ana_start_time..ana_end_time) end ave = Hash.new stddev = Hash.new # 領域平均 for range in ranges puts "" x_start = range[0] x_end = range[1] y_start = range[2] y_end = range[3] puts "" for varname in var_list tmp = var[varname].cut("lon"=>x_start..x_end, "lat"=>y_start..y_end) tmp = tmp.mean("lon") tmp_ave = 0.0 for y in lat.cut("lat"=>y_start..y_end).val weight = lat_weight.cut("lat"=>y).val.to_f # weight = lat_weight.cut("lat"=>y).val tmp_ave += tmp.cut("lat"=>y) * weight end # 全球平均でも {南, 北} 半球平均でも重みが同じになるように補正 if lat.cut("lat"=>y_start..y_end).val.size == nlat tmp_ave /= 2.0 end ave[varname] = tmp_ave.mean("time") # stddev[varname] = tmp_ave.stddev("time") puts "" end total_heat = (ave[name_OSR] + ave[name_OLR]) \ + (ave[name_SSR] + ave[name_SLR]) + ave[name_Sens] + ave[name_Rain] # total_heat = (ave[name_SLR] - ave[name_OLR]) + ave[name_Sens] + ave[name_Rain] total_wat = ave[name_Evap] - ave[name_Rain] total_surf = ave[name_SSR] + ave[name_SLR] + ave[name_Sens] + ave[name_Evap] total_top = ave[name_OSR] + ave[name_OLR] puts "" puts "" puts "" puts "" puts "" puts "" end # 指定された点での解析 for point in points puts "" point_lon = point[0] point_lat = point[1] puts "" for varname in var_list tmp = var[varname].cut("lon"=>point_lon, "lat"=>point_lat) ave[varname] = tmp.mean("time") # stddev[varname] = tmp.stddev("time") puts "" end total_heat = (ave[name_OSR] + ave[name_SSR]) \ + (ave[name_SLR] + ave[name_OLR]) + ave[name_Sens] + ave[name_Rain] # total_heat = (ave[name_SLR] - ave[name_OLR]) + ave[name_Sens] + ave[name_Rain] total_wat = ave[name_Evap] - ave[name_Rain] total_surf = ave[name_SSR] + ave[name_SLR] + ave[name_Sens] + ave[name_Evap] total_top = ave[name_OSR] + ave[name_OLR] puts "" puts "" puts "" puts "" puts "" puts "" end # html を閉じる puts "

" + "Point or Range" + "	" + varname + "	" + "Heat Budget" + "	" + "Water Budget" + "	" + "Surface Budget" + "	" + "Top Budget" + "
(" + sprintf("%#.1f", x_start) + ":" + sprintf("%#.1f", x_end) + ", " + sprintf("%#.1f", y_start) + ":" + sprintf("%#.1f", y_end) + ")	" + sprintf("%#.3f", ave[varname].val) + "	" + sprintf("%#.3f", total_heat.val) + "	" + sprintf("%#.3f", total_wat.val) + "	" + sprintf("%#.3f", total_surf.val) + "	" + sprintf("%#.3f", total_top.val) + "
(" + point_lon.to_s + ", " + point_lat.to_s + ")	" + sprintf("%#.3f", ave[varname].val) + "	" + sprintf("%#.3f", total_heat.val) + "	" + sprintf("%#.3f", total_wat.val) + "	" + sprintf("%#.3f", total_surf.val) + "	" + sprintf("%#.3f", total_top.val) + "

" puts "" __END__