#! /usr/bin/ruby # -*- coding: euc-jp -*- # # 2013-5-12 石渡正樹 # バンドごとの長波放射加熱率 require "getoptlong" # for option_parse require "numru/ggraph" include NumRu ## オプション解析 parser = GetoptLong.new parser.set_options( ### global option ### ['--lon', GetoptLong::REQUIRED_ARGUMENT], ['--lat', GetoptLong::REQUIRED_ARGUMENT], ['--time_start', GetoptLong::REQUIRED_ARGUMENT], ['--time_end', GetoptLong::REQUIRED_ARGUMENT], ['--itr', GetoptLong::REQUIRED_ARGUMENT], ['--cint', GetoptLong::REQUIRED_ARGUMENT], ['--wsn', GetoptLong::REQUIRED_ARGUMENT] ) parser.each_option do |name, arg| eval "$OPT_#{name.sub(/^--/, '').gsub(/-/, '_')} = '#{arg}'" # strage option value to $OPT_val end time_start = ($OPT_time_start||1).to_f time_end = ($OPT_time_end||9999).to_f lon = ($OPT_lon||0).to_s lat_opt = ($OPT_lat||0).to_s cint = ($OPT_cint||1e-5).to_f itr = ($OPT_itr||2).to_f # band number 軸の生成 bandnum_a = VArray.new( NArray.sfloat(10).indgen(1,1), {"long_name"=>"Band Number", "units"=>""}, "BandNum" ) bandnum_axis = Axis.new.set_pos(bandnum_a) p 'in DTempDtRadL: cint=', cint # sig 軸の生成 sig = GPhys::IO.open('DTempDtRadLBand1.nc', 'sig').val #sig_weight = GPhys::IO.open('DTempDtRadLBand1.nc', 'sig_weight').val sig_a = VArray.new( sig, {"long_name"=>"atmosphere_sigma_coordinate", "units"=>"1"}, "sig" ) sig_axis = Axis.new.set_pos(sig_a) # lat 軸情報の取得 lat = GPhys::IO.open('DTempDtRadLBand1.nc', 'lat').val lat_weight = GPhys::IO.open('DTempDtRadLBand1.nc', 'lat_weight').val nlat = lat.size lat_weight_dummy = VArray.new( NArray.sfloat(lat.size,sig.size), {"long_name"=>"lat_weigh", "units"=>""}, "lat_weight" ) for j in 0..lat.size-1 lat_weight_dummy[j,true]= lat_weight[j] end #p lat_weight_dummy # データを格納する配列を作る (まだ箱だけ) DTempDtRadL_a = VArray.new( NArray.sfloat(sig.size,10), {"long_name"=>"long wave radiative heating rate for bands", "units"=>"K s-1 m"}) ## データ切出し, Ts 軸データの作成 for n in bandnum_a.val.to_i lon = lon.to_s p lon.type varname = 'DTempDtRadLBand' + n.to_s # たとえば, DTempDtRadLBand1 ncfile = varname + '.nc' # たとえば, DTempDtRadLBand1.nc if lat_opt == 'mean' then # 緯度方向には極から極まで平均する場合 if lon == 'mean' then # NArray の和をとるので sum の引数として次元の番号を与えた. DTempDtRadL_a[true,n-1] = ( GPhys::IO.open(ncfile, varname).cut('time'=>time_start..time_end).mean('time').mean('lon').val * 0.5 * lat_weight_dummy.val ).sum(0) elsif lon.include?(':') then # 経度方向の範囲を指定した場合 lon_start = lon.split(':')[0].to_f lon_end = lon.split(':')[1].to_f # DTempDtRadL_a[true,n-1] = GPhys::IO.open(ncfile, varname).cut('time'=>time_start..time_end,'lon'=>lon_start..lon_end,'lat'=>0.0).mean('lon').mean('time').val DTempDtRadL_a[true,n-1] = (GPhys::IO.open(ncfile, varname).cut('time'=>time_start..time_end,'lon'=>lon_start..lon_end).mean('lon').mean('time').val * 0.5 * lat_weight_dummy.val ).sum(0) else # 経度には関しては 1 点が指定されているはず. lon = lon.to_f DTempDtRadL_a[true,n-1] = (GPhys::IO.open(ncfile, varname).cut('time'=>time_start..time_end,'lon'=>lon).mean('time').val * 0.5 * lat_weight_dummy.val ).sum(0) end # elsif lon == 'mean' and lat.include?(':') then # lat_start = lat.split(':')[0].to_f # lat_end = lat.split(':')[1].to_f else # lat として 1 点が指定された場合 lat_opt = lat_opt.to_f if lon == 'mean' then # 経度方向には全部平均 DTempDtRadL_a[true,n-1] = GPhys::IO.open(ncfile, varname).cut('time'=>time_start..time_end,'lat'=>lat_opt).mean('time').mean('lon').val elsif lon.include?(':') then # 経度方向の範囲を指定した場合 lon_start = lon.split(':')[0].to_f lon_end = lon.split(':')[1].to_f DTempDtRadL_a[true,n-1] = GPhys::IO.open(ncfile, varname).cut('time'=>time_start..time_end,'lon'=>lon_start..lon_end,'lat'=>lat_opt).mean('lon').mean('time').val else # 経度には関しては 1 点が指定されているはず. lon = lon.to_f DTempDtRadL_a[true,n-1] = GPhys::IO.open(ncfile, varname).cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>lat_opt).mean('time').val end end end DTempDtRadL = GPhys.new(Grid.new(sig_axis, bandnum_axis), DTempDtRadL_a) ## 作図 DCL.gropn($OPT_wsn||4) DCL.sgpset('lcntl', false) ; DCL.uzfact(0.7) GGraph.set_fig('itr'=> itr, 'viewport'=>[0.2,0.8,0.3,0.6] , "yreverse"=>true) GGraph.contour( DTempDtRadL, true, "annotate"=>false, "exchange"=>true,"interval"=>cint) DCL.grcls exit # 前バージョンのデータ格納 DTempDtRadL_a[true,0] = GPhys::IO.open('DTempDtRadLBand1.nc', 'DTempDtRadLBand1').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,1] = GPhys::IO.open('DTempDtRadLBand2.nc', 'DTempDtRadLBand2').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,2] = GPhys::IO.open('DTempDtRadLBand3.nc', 'DTempDtRadLBand3').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,3] = GPhys::IO.open('DTempDtRadLBand4.nc', 'DTempDtRadLBand4').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,4] = GPhys::IO.open('DTempDtRadLBand5.nc', 'DTempDtRadLBand5').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,5] = GPhys::IO.open('DTempDtRadLBand6.nc', 'DTempDtRadLBand6').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,6] = GPhys::IO.open('DTempDtRadLBand7.nc', 'DTempDtRadLBand7').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,7] = GPhys::IO.open('DTempDtRadLBand8.nc', 'DTempDtRadLBand8').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,8] = GPhys::IO.open('DTempDtRadLBand9.nc', 'DTempDtRadLBand9').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val DTempDtRadL_a[true,9] = GPhys::IO.open('DTempDtRadLBand10.nc', 'DTempDtRadLBand10').cut('time'=>time_start..time_end,'lon'=>lon,'lat'=>0.0).mean('time').val