#!/usr/bin/env ruby =begin = mkfig.rb -- make figure for energy balance 南北エネルギー輸送関連図. 以下の物理量の経度鉛直積分値をプロット. * 乾燥静的エネルギー * CpvT + vh * 潜熱輸送 * LQv =end require "getopts" # for option_parse require "numru/netcdf_miss" require "numru/ggraph" require "colorbar" require "libgphys-n" include NumRu include Misc::EMath ## 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(var, 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, var) return gp end ###################################################### Rad = UNumeric.new(6.37E6, 'm') G = UNumeric.new(9.81 , 'm.s-2') year = [1979, 2003] # year = [1978, 1978] # year = [2000, 2000] winter = ["12", "01", "02"]; spring = ["03", "04", "05"] summer = ["06", "07", "08"]; fall = ["09", "10", "11"] annual = winter + spring + summer + fall djf_le = []; mam_le = []; jja_le = []; son_le = []; ann_le = [] djf_se = []; mam_se = []; jja_se = []; son_se = []; ann_se = [] djf_ta = []; mam_ta = []; jja_ta = []; son_ta = []; ann_ta = [] kikan = {# 各季節に対応する月名の配列を値にもつハッシュ 'winter'=>winter, 'spring'=>spring, 'summer'=>summer, 'fall' =>fall, 'annual'=>annual } gp_se = {# 各季節平均の GPhys オブジェクトを持つハッシュ 'winter'=>djf_se, 'spring'=>mam_se, 'summer'=>jja_se, 'fall' =>son_se, 'annual'=>ann_se } gp_le = {# 各季節平均の GPhys オブジェクトを持つハッシュ 'winter'=>djf_le, 'spring'=>mam_le, 'summer'=>jja_le, 'fall' =>son_le, 'annual'=>ann_le } gp_ta = {# 各季節平均の GPhys オブジェクトを持つハッシュ 'winter'=>djf_ta, 'spring'=>mam_ta, 'summer'=>jja_ta, 'fall' =>son_ta, 'annual'=>ann_ta } ['spring', 'summer', 'fall', 'winter', 'annual'].each {|ss| ## 乾燥静的エネルギーフラックス dsefl = make_mean_gphys('dsefl', 'DSEFL', 'NCEP', kikan[ss], year) dsefl_bar = make_mean_gphys('dsefl_bar', 'DSEFL', 'NCEP', kikan[ss], year) dsefl_dash = make_mean_gphys('dsefl_dash','DSEFL', 'NCEP', kikan[ss], year) ## 潜熱輸送量 lhefl = make_mean_gphys('lhefl', 'LHEFL', 'NCEP', kikan[ss], year) lhefl_bar = make_mean_gphys('lhefl_bar', 'LHEFL', 'NCEP', kikan[ss], year) lhefl_dash = make_mean_gphys('lhefl_dash','LHEFL', 'NCEP', kikan[ss], year) # get cos lat = lhefl.coord('lat') grid = Grid.new(lhefl.axis('lat')) cos_lat = GPhys.new(grid, cos(lat*2*PI/360.0)) dse = dsefl.integrate('level') /G * Rad * cos_lat * 2 * PI * 100 dse_bar = dsefl_bar.integrate('level') /G * Rad * cos_lat * 2 * PI * 100 dse_dash = dsefl_dash.integrate('level') /G * Rad * cos_lat * 2 * PI * 100 dse.data.set_att('units', 'W') dse_bar.data.set_att('units', 'W') dse_dash.data.set_att('units', 'W') lhe = lhefl.integrate('level') /G * Rad * cos_lat * 2 * PI * 100 lhe_bar = lhefl_bar.integrate('level') /G * Rad * cos_lat * 2 * PI * 100 lhe_dash = lhefl_dash.integrate('level') /G * Rad * cos_lat * 2 * PI * 100 lhe.data.set_att('units', 'W') lhe_bar.data.set_att('units', 'W') lhe_dash.data.set_att('units', 'W') gp_se[ss] << dse gp_se[ss] << dse_bar gp_se[ss] << dse_dash gp_le[ss] << lhe gp_le[ss] << lhe_bar gp_le[ss] << lhe_dash gp_ta[ss] << dse + lhe gp_ta[ss] << dse_bar + lhe_bar gp_ta[ss] << dse_dash + lhe_dash } ################################################################### # set User Path for dcldatabase DCL.glcset('DUPATH','/home/daktu32/.dcldir/') ## four in one(itr 2) itr = 1 min = 8e15 max = -8e15 ## annual mean DCL.swpset('lsep', true) # ページ毎別々のファイルに落す DCL.gropn(2) DCL.sgpset('lcntl', false) # DCL.sgpset('lfull',true) # フルスクリーン DCL.uzfact(0.35) # DCL.sgpset('lcorner',false) # コーナーを取っ払う DCL.sgpset('lfprop',true) # DCL.udpset('lmsg',false) # DCL.uscset('cyspos', 'B' ) # y 軸の単位の位置を下方へ ############################################################################## ##-------------------------------------- ## Sensible heat flux (seasonal) # ビューポート設定 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_se[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_se, '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_se, '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_se, '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_se, '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, 0.3, 0.1) ## -- extra information -- DCL::sgtxzv(0.5,vpt00[3]+0.03,"25 years mean Dry Static Energy Transport (1979-2003)", 1.15*DCL.uzpget('rsizec2'),0,0,3) DCL::sgtxzv(0.9,vpt00[3]+0.01,"(units:#{djf_se[0].data.units.to_s})", 0.9*DCL.uzpget('rsizec2'),0,0,3) ##--------------------------------------------------- ## Sensible Heat (annual) vpt = NArray[0.15,0.85,0.23,0.58] # dcl に欠損値情報を教える before = DCLExt.gl_set_params('lmiss'=>true,'rmiss'=>ann_se[0].data.get_att("missing_value")[0]) # 文字サイズ変更 DCL.uzfact(11.0/7) # annual mean DCL.uzpset('labelyl',true) GGraph.set_fig('viewport'=>vpt, 'itr'=>itr, 'new_frame'=>true) GGraph.set_axes('ytitle'=>'') names, idxs = plot_line_main(ann_se, 'annot'=>false, 'titl'=>'', 'min'=>min, 'max'=>max) __show_line_index(names,idxs) ## -- extra information -- DCL::sgtxzv(0.5,vpt[3]+0.028,"25 years annual mean Dry Static Energy Transport (1979-2003)", 1.15*DCL.uzpget('rsizec2'),0,0,3) DCL::sgtxzv(0.78,vpt[3]+0.01,"(units:#{ann_se[0].data.units.to_s})", 0.8*DCL.uzpget('rsizec2'),0,0,3) ############################################################################## ##-------------------------------------- ## Latent Heat (seasonal) # ビューポート設定 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_le[0].data.get_att("missing_value")[0]) # 文字サイズ変更 DCL.uzfact(7/11.0) # 冬平均(左上の絵) GGraph.set_fig('viewport'=>vpt00, 'itr'=>itr, 'new_frame'=>true) GGraph.set_axes('xunits'=>'','yunits'=>'','xtitle'=>'', 'ytitle'=>'') DCL.uzpset('labelxb',false); DCL.uzpset('labelyl',true) plot_line_main(mam_le, '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_le, '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_le, '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_le, '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, 0.3, 0.1) ## -- extra information -- DCL::sgtxzv(0.5,vpt00[3]+0.032,"25 years mean Latent Heat Transport (1979-2003)", 1.15*DCL.uzpget('rsizec2'),0,0,3) DCL::sgtxzv(0.89,vpt00[3]+0.01,"(units:#{djf_le[0].data.units.to_s})", 0.9*DCL.uzpget('rsizec2'),0,0,3) ##--------------------------------------------------- ## Latent Heat (annual) # annual mean vpt = NArray[0.15,0.85,0.23,0.58] # 文字サイズ変更 DCL.uzfact(11.0/7) DCL.uzpset('labelyl',true) GGraph.set_fig('viewport'=>vpt, 'itr'=>itr, 'new_frame'=>true) GGraph.set_axes('ytitle'=>'') names, idxs = plot_line_main(ann_le, 'annot'=>false, 'titl'=>'', 'min'=>min, 'max'=>max) __show_line_index(names,idxs) ## -- extra information -- DCL::sgtxzv(0.5,vpt[3]+0.04,"25 years annual mean Latent Heat Energy Transport (1979-2003)", 1.15*DCL.uzpget('rsizec2'),0,0,3) DCL::sgtxzv(0.82,vpt[3]+0.01,"(units:#{ann_le[0].data.units.to_s})", 0.8*DCL.uzpget('rsizec2'),0,0,3) ############################################################################## ##-------------------------------------- ## Total Energy Transport (seasonal) # ビューポート設定 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_ta[0].data.get_att("missing_value")[0]) # 文字サイズ変更 DCL.uzfact(7/11.0) # 冬平均(左上の絵) GGraph.set_fig('viewport'=>vpt00, 'itr'=>itr, 'new_frame'=>true) GGraph.set_axes('xunits'=>'','yunits'=>'','xtitle'=>'', 'ytitle'=>'') DCL.uzpset('labelxb',false); DCL.uzpset('labelyl',true) plot_line_main(mam_ta, '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_ta, '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_ta, '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_ta, '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, 0.3, 0.1) ## -- extra information -- DCL::sgtxzv(0.5,vpt00[3]+0.032,"25 years mean Atmospheric Energy Transport (1979-2003)", 1.15*DCL.uzpget('rsizec2'),0,0,3) DCL::sgtxzv(0.89,vpt00[3]+0.01,"(units:#{djf_ta[0].data.units.to_s})", 0.9*DCL.uzpget('rsizec2'),0,0,3) ##--------------------------------------------------- ## Latent Heat (annual) # annual mean vpt = NArray[0.15,0.85,0.23,0.58] # 文字サイズ変更 DCL.uzfact(11.0/7) DCL.uzpset('labelyl',true) GGraph.set_fig('viewport'=>vpt, 'itr'=>itr, 'new_frame'=>true) GGraph.set_axes('ytitle'=>'') names, idxs = plot_line_main(ann_ta, 'annot'=>false, 'titl'=>'', 'min'=>min, 'max'=>max) __show_line_index(names,idxs) ## -- extra information -- DCL::sgtxzv(0.5,vpt[3]+0.04,"25 years annual mean Atmospheric Energy Transport (1979-2003)", 1.15*DCL.uzpget('rsizec2'),0,0,3) DCL::sgtxzv(0.82,vpt[3]+0.01,"(units:#{ann_ta[0].data.units.to_s})", 0.8*DCL.uzpget('rsizec2'),0,0,3) DCL.grcls #------------------------------------------------------------------- ################################################## # 画像ファイル名を変更 File.rename('dcl_001.ps', 'ENERGY_TRANSPORT_1973-2003-SEASON_DS_NCEP.ps') File.rename('dcl_002.ps', 'ENERGY_TRANSPORT_1973-2003-ANNUAL_DS_NCEP.ps') File.rename('dcl_003.ps', 'ENERGY_TRANSPORT_1973-2003-SEASON_LE_NCEP.ps') File.rename('dcl_004.ps', 'ENERGY_TRANSPORT_1973-2003-ANNUAL_LE_NCEP.ps') File.rename('dcl_005.ps', 'ENERGY_TRANSPORT_1973-2003-SEASON_TA_NCEP.ps') File.rename('dcl_006.ps', 'ENERGY_TRANSPORT_1973-2003-ANNUAL_TA_NCEP.ps')