#!/usr/bin/env ruby =begin = mkfig.rb -- make figure for energy balance エネルギーバランス関連図. 以下の物理量の緯度平均値をプロット. * NSWRS : 地表面ネット短波フラックス * NLWRS : 地表面ネット長波フラックス * USWRF : 大気上端外向き短波フラックス * ULWRF : 大気上端外向き長波フラックス * DSWRF : 大気上端内向き短波フラックス * PRATE : 降水強度 * 温度 0 ℃ における潜熱 2.5008 x 10e6 J/kg を用いた * LHTFL : 地表面における潜熱フラックス * SHTFL : 地表面における顕熱フラックス 描けたら, 正味のフラックスを描く. =end require "getopts" # for option_parse require "numru/netcdf_miss" require "numru/ggraph" require "colorbar" require "libgphys-n" include NumRu include Misc ## enable to the option min, max def GGraph::line(gphys, newframe=true, options=nil) if newframe!=true && newframe!=false raise ArgumentError, "2nd arg (newframe) must be true or false" end if ! defined?(@@line_options) @@line_options = Misc::KeywordOptAutoHelp.new( ['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'], ['exchange', false, 'whether to exchange x and y axes'], ['index', 1, 'line/mark index'], ['type', 1, 'line type'], ['label', nil, 'if a String is given, it is shown as the label'], ['max', nil, 'maximam value for draw line'], ['min', nil, 'minimam value for draw line'] ) end opts = @@line_options.interpret(options) gp = gphys.first1D if !opts['exchange'] x = gp.coord(0) y = gp.data xax = x if opts['min'] ymin = opts['min'].to_f else ymin = gp.data.val.min end if opts['max'] ymax = opts['max'].to_f else ymax = gp.data.val.max end yax = VArray.new(NArray[ymin, ymax], gp.data, gp.data.name) else y = gp.coord(0) x = gp.data yax = y if opts['min'] xmin = opts['min'] else xmin = gp.data.val.min end if opts['max'] xmax = opts['max'] if opts['max'] else xmax = gp.data.val.max end xax = VArray.new(NArray[xmin, xmax], gp.data, gp.data.name) end if newframe fig(xax, yax) axes(xax, yax) title( (opts['title'] || gp.data.get_att('long_name')) ) annotate(gp.lost_axes) if opts['annotate'] end if opts['label'] lcharbk = DCL.sgpget('lchar') DCL.sgpset('lchar',true) DCL.sgsplc(opts['label']) end DCL.uulinz(x.val, y.val, opts['type'], opts['index']) DCL.sgpset('lchar',lcharbk) if opts['label'] nil end ##----------------------- # ラインメソッド再定義? def plot_line_main(gphys_array, line_opts={}) line_index_ary = Array.new; name_ary = Array.new default_index = 12 line_hash = { "index"=> default_index }.update(line_opts) line_index_ary.push(default_index) name_ary.push(gphys_array[0].data.get_att("long_name")) GGraph.line( gphys_array[0], true, line_hash ) gphys_array.size.times{ |num| unless num == 0 line_hash = { "index"=> (num*10 + default_index)} line_index_ary.push(num*10 +default_index) name_ary.push(gphys_array[num].data.get_att("long_name")) GGraph.line( gphys_array[num], false, line_hash ) end } return name_ary, line_index_ary end ##----------------------- # ラインインデックスの種類を表示 def __show_line_index(str_ary,line_index_ary, x=0.15, y=0.12) charsize = 1.0 * DCL.uzpget('rsizec1') dvx = 0.01 dvy = charsize*1.5 raise TypeError,"Array expected" if ! str_ary.is_a?(Array) vxmin,vxmax,vymin,vymax = DCL.sgqvpt vx = x vy = y + charsize/2 str_ary.size.times{|num| DCL::sgtxzv(vx,vy,"--- #{str_ary[num]}", charsize, 0, -1, line_index_ary[num]) vy -= dvy if num == 6 vx = 0.30 vy = 0.12 - charsize/2 end } nil end def make_mean_gphys(phys, src, season,year) ary = [] y1 = year[0]; y2 = year[1] path = "../../#{phys}.#{src}/" season.each do |m| y1.upto(y2) do |y| fn = path + "#{phys}.#{y}.#{src}/#{phys}_#{y}-#{m}_#{src}.nc" ary << fn end end gp = mean_gphys(ary, phys.downcase) return gp.mean('lon') end ################################################################### # make gphys objects year = [1979, 2003] # year = [1979, 1979] L = UNumeric.new(2.5008e6, 'J.kg-1') winter = ["12", "01", "02"]; spring = ["03", "04", "05"] summer = ["06", "07", "08"]; fall = ["09", "10", "11"] djf = []; mam = []; jja = []; son = [] season_kikan = {# 各季節に対応する月名の配列を値にもつハッシュ 'winter'=>winter, 'spring'=>spring, 'summer'=>summer, 'fall' =>fall } season_gp = {# 各季節平均の GPhys オブジェクトを持つハッシュ 'winter'=>djf, 'spring'=>mam, 'summer'=>jja, 'fall' =>son } ['spring', 'summer', 'fall', 'winter'].each {|ss| ## 大気に吸収される正味の短波放射フラックス # nswrs = make_mean_gphys('NSWRS', 'NCEP', season_kikan[ss], year) uswrf = make_mean_gphys('USWRF', 'NCEP', season_kikan[ss], year) dswrf = make_mean_gphys('DSWRF', 'NCEP', season_kikan[ss], year) nswr = dswrf - uswrf# + nswrs nswr.data.set_att('long_name', 'net absorbed shortwave radiation flux at nominal top of atmosphere') ## 大気が射出する正味の惑星放射フラックス #nlwrs = make_mean_gphys('NLWRS', 'NCEP', season_kikan[ss], year) ulwrf = make_mean_gphys('ULWRF', 'NCEP', season_kikan[ss], year) nlwr = ulwrf# - nlwrs nlwr.data.set_att('long_name', 'net outward longwave radiation flux at nominal top of atmosphere') ## 潜熱フラックス lhtfl = make_mean_gphys('LHTFL', 'NCEP', season_kikan[ss], year) lhtfl.data.set_att('long_name', 'net latent heat flux at surface') ## 顕熱フラックス shtfl = make_mean_gphys('SHTFL', 'NCEP', season_kikan[ss], year) shtfl.data.set_att('long_name', 'net sensible heat flux at surface') ## 降水(潜熱に換算したもの) prate = make_mean_gphys('PRATE', 'NCEP', season_kikan[ss], year) lheat = prate*L lheat.data.set_att('long_name', 'latent heat with precipitation') season_gp[ss] << nswr season_gp[ss] << nlwr season_gp[ss] << lhtfl season_gp[ss] << shtfl season_gp[ss] << lheat } # set User Path for dcldatabase DCL.glcset('DUPATH','/home/daktu32/.dcldir/') ## four in one(itr 2) itr = 1 min = -300 max = +500 #min = -1e+10 #max = +2e+10 DCL.gropn(2) DCL.sgpset('lcntl', false) # DCL.sgpset('lfull',true) # DCL.sgpset('lcorner',false) # DCL.uzfact(0.35) # DCL.sgpset('lfprop',true) # DCL.udpset('lmsg',false) # # ビューポート設定 vpt = NArray[0.05, 0.45, 0.05, 0.25] # vpt00 = ( vpt + ([0.050]*2 + [0.32]*2) ).to_a # vpt01 = ( vpt + ([0.474]*2 + [0.32]*2) ).to_a # vpt10 = ( vpt + ([0.050]*2 + [0.10]*2) ).to_a # vpt11 = ( vpt + ([0.474]*2 + [0.10]*2) ).to_a # # dcl に欠損値情報を教える before = DCLExt.gl_set_params('lmiss'=>true,'rmiss'=>djf[0].data.get_att("missing_value")[0]) # 冬平均(左上の絵) GGraph.set_fig('viewport'=>vpt00, 'itr'=>itr) GGraph.set_axes('xunits'=>'','yunits'=>'','xtitle'=>'', 'ytitle'=>'') DCL.uzpset('labelxb',false) plot_line_main(mam, 'annot'=>false, 'titl'=>'', 'min'=>min, 'max'=>max) DCL.uzpset('pad1',0.2) ; DCL.uxsttl('t','Spring(MAM)', -1) ; DCL.uzpset('pad1',0.7) # 春平均(右上の絵) GGraph.set_fig('viewport'=>vpt01, 'new_frame'=>false, 'itr'=>itr) GGraph.set_axes('ytitle'=>'') DCL.uzpset('labelyl',false) plot_line_main(jja, 'annot'=>false, 'titl'=>'', 'min'=>min, 'max'=>max) DCL.uzpset('pad1',0.2) ; DCL.uxsttl('t','Summer(JJA)',-1) ; DCL.uzpset('pad1',0.7) # 夏平均(左下の絵) GGraph.set_fig('viewport'=>vpt10, 'new_frame'=>false, 'itr'=>itr) GGraph.set_axes('ytitle'=>'','xtitle'=>nil) DCL.uzpset('labelyl',true); DCL.uzpset('labelxb',true) plot_line_main(son, 'annot'=>false, 'titl'=>'', 'min'=>min, 'max'=>max) DCL.uzpset('pad1',0.2) ; DCL.uxsttl('t','Autumn(SON)',-1) ; DCL.uzpset('pad1',0.7) # 秋平均(右下の絵) GGraph.set_fig('viewport'=>vpt11, 'new_frame'=>false, 'itr'=>itr) GGraph.set_axes('ytitle'=>'') DCL.uzpset('labelyl',false) names, idxs = plot_line_main(djf, 'annot'=>false, 'titl'=>'', 'min'=>min, 'max'=>max) DCL.uzpset('pad1',0.2) ; DCL.uxsttl('t','Winter(DJF)',-1) ; DCL.uzpset('pad1',0.7) __show_line_index(names,idxs) ## -- extra information -- DCL::sgtxzv(0.5,vpt00[3]+0.028,"25 years mean Heat Balance (1979-2003)", 1.15*DCL.uzpget('rsizec2'),0,0,3) DCL::sgtxzv(0.88,vpt00[3]+0.01,"(units:#{djf[0].data.units.to_s})", 0.9*DCL.uzpget('rsizec2'),0,0,3) DCL.grcls