require 'numru/ggraph' require 'numru/gphys/ep_flux' require 'numru/ggraph_on_merdional_section' require 'colorbar' require 'getopts' include NumRu ############## define methods module NumRu class NetCDFVar def get_with_miss_and_scaling2(*args) # make mask before scaling __interpret_missing_params if !defined?(@missval) packed_data = simple_get(*args) scaled_data = scaled_get(*args) sf = att('scale_factor') ao = att('add_offset') if @vmin || @vmax if sf && ao csf = sf.get cao = ao.get vmin = (@vmin-cao)/csf if @vmin vmax = (@vmax-cao)/csf if @vmax elsif vmin = @vmin; vmax = @vmax end if vmin mask = (packed_data >= vmin) mask = mask.and(packed_data <= vmax) if vmax else mask = (packed_data <= vmax) end data = NArrayMiss.to_nam(scaled_data, mask) elsif @missval # only missing_value is present. eps = 1e-6 missval = @missval[0].to_f vmin = missval - eps vmax = missval + eps mask = (packed_data <= vmin) mask = mask.or(packed_data >= vmax) data = NArrayMiss.to_nam(scaled_data, mask) else data = scaled_data end data end end module GGraph module_function def vector_for_dcchart(fx, fy, newframe=true, options=nil) if ! defined?(@@vector_for_dcchart_options) @@vector_for_dcchart_options = Misc::KeywordOptAutoHelp.new( ['newfig', true, 'if false, do not cleared before figure setting.'], ['title', nil, 'Title of the figure(if nil, internally determined)'], ['annotate', true, 'if false, do not put texts on the right margin even when newframe==true'], ['transpose', false, 'if true, exchange x and y axes'], ['xintv', 1, '(Effective only if flow_vect) interval of data sampling in x'], ['yintv', 1, '(Effective only if flow_vect) interval of data sampling in y'], ['xfact', 1, 'scale factor for x'], ['yfact', 1, 'scale factor for y'], ['itr', 1, 'log or linear'], ['unit_vect', false, 'Show the unit vector'], ['max_unit_vect', false, '(Effective only if flow_vect && unit_vect) If true, use the maximum arrows to scale the unit vector; otherwise, normalize in V coordinate.'] ) end opts = @@vector_for_dcchart_options.interpret(options) fx = fx.first2D.copy fy = fy.first2D.copy sh = fx.shape if sh != fy.shape raise ArgumentError, "shapes of fx and fy do not agree with each other" end fx = fx.transpose(1,0) if opts['transpose'] fy = fy.transpose(1,0) if opts['transpose'] if ((xi=opts['xintv']) >= 2) idx = NArray.int(sh[0]/xi).indgen!*xi # [0,xi,2*xi,..] fx = fx[idx, true] fy = fy[idx, true] end if ((yi=opts['xintv']) >= 2) idx = NArray.int(sh[1]/yi).indgen!*yi # [0,yi,2*yi,..] fx = fx[true, idx] fy = fy[true, idx] end xax = fx.coord(0) yax = fy.coord(1) if newframe nextfig = @@next_fig.dup if ( @@next_fig != nil ) # backup next_fig fig(xax, yax, {'itr'=>itr}) axes(xax, yax) if opts['title'] ttl = opts['title'] else fxnm = fx.data.get_att('long_name') || fx.name fynm = fy.data.get_att('long_name') || fy.name ttl = '('+fxnm+','+fynm+')' end title( ttl ) annotate(fx.lost_axes) if opts['annotate'] @@next_fig = nextfig if ( @@next_fig != nil ) end fig(xax, yax, {'new_frame'=>false, 'itr'=>1, 'yreverse'=>false}) \ if (opts['newfig']) DCL.uwsgxa(xax.val) DCL.uwsgya(yax.val) xfact = opts["xfact"]; yfact = opts["yfact"] before2=DCLExt.ug_set_params( {'LNRMAL'=>false, 'LMSG'=>false, 'XFACT1'=>xfact, 'YFACT1'=>yfact} ) before=DCLExt.ug_set_params({'lunit'=>true}) if opts['unit_vect'] DCL.ugvect(fx[true, -1..0].val, fy[true, -1..0].val) DCLExt.ug_set_params(before) if opts['unit_vect'] nil end end end def mean_gphys(file, var) gphys_file = Array.new() file.each {|f| temp = GPhys::NetCDF_IO.open(f, var) gphys_file << temp } sum = gphys_file[0] 1.upto(gphys_file.size-1) {|f| sum += gphys_file[f] } p gphys = sum/file.size.to_f return gphys end ############## main routine levels = NArray[-1.0e8, -5.0e1, -1.0e1, 0, 1.0e1, 5.0e1, 1.0e8]*1.0e-6 patterns = NArray[30999, 35999, 40999, 55999, 75999, 85999] #vect_ratio = NArray[1.0, 500.0]*5e-12 vect_ratio = NArray[1.0, 100.0]*8e-12 vect_fact = 1.0 unless getopts("", "title:" ,"vect_fact:") print "#{$0}:illegal options. please exec this program with -h/--help. \n" exit 1 end vect_fact = ($OPT_vect_fact) if ($OPT_vect_fact) files = ARGV epynm = 'epflx_phi'; epznm = 'epflx_p'; epdivnm = 'epflx_div' #epflx_y = mean_gphys(files, epynm).mean('time') #epflx_z = mean_gphys(files, epznm).mean('time') #epflx_div = mean_gphys(files, epdivnm).mean('time') epflx_y = mean_gphys(files, epynm) epflx_z = - mean_gphys(files, epznm) epflx_div = mean_gphys(files, epdivnm) filename = File::basename(files[0]) filename = ($OPT_title) if ($OPT_title) DCL.gropn(2) DCL::sglset('LFULL', true) # use full area DCL::slrat(1.0, 0.85) # set aspect ratio of drwable area DCL.sgpset('lcntl', false) # don't rede control character. DCL.sgpset('lfprop',true) # use prportional font DCL.uzfact(0.6) # set character size before = DCLExt.gl_set_params('lmiss'=>true,'rmiss'=>epflx_div.data.get_att('missing_value')[0]) #before = DCLExt.gl_set_params('lmiss'=>true,'rmiss'=>32766) GGraph::set_fig('view'=>[0.15,0.85,0.23,0.58],"itr"=>2) ## show tone div F GGraph.tone(epflx_div, true, "lev"=>levels, "patterns"=>patterns, "title"=>filename ) DCL::Util::color_bar({ "levels"=>levels, "colors"=>patterns, "nobound"=>0, "eqlev"=>true, "tick1"=>10,"tick2"=>1 }) GGraph.contour(epflx_div, false, "lev"=>levels, "title"=>filename ) ## show vector (Fy, Fz) GGraph::vector_for_dcchart(epflx_y, epflx_z, false, # 'xfact'=>5e-11, 'yfact'=>2.5e-8, 'unit'=>true, 'annot'=>false, 'itr'=>2 'xfact'=>vect_ratio[0]*vect_fact, 'yfact'=>vect_ratio[-1]*vect_fact, 'unit'=>true, 'annot'=>false, 'itr'=>2, 'yintv'=>2 ) #GGraph::vector_on_merdional_section(epflx_y, epflx_z, false, # ''unit'=>true, 'annot'=>false # ) DCL.grcls