#!/usr/bin/env ruby =begin $ for file in *.xwd ; do xwdtopnm ${file} | pnmcut 2 2 904 654 > test.pnm; ppmtogif test.pnm > ${file}.gif; done $ gifsicle --delay=50 --no-loop dcl_*.gif > movie.gif =end # ---------------------------------------------- # プログラムメイン END{ num = 6 set_filename(num) # ファイル名等基本設定 mradl_anim_main # gif 作成 print "gif convert start #{Time.now}\n" mkgifanim # gif アニメ作成シェルスクリプト } # ---------------------------------------------- # 必要なライブラリ, モジュールの読み込み $local_path = '/home/yukiko/work/ape/yukiko/lib' $: << $local_path require "getopts" require "numru/netcdf" module NumRu class NetCDFVar alias put scaled_put alias get scaled_get end end require "numru/ggraph" require "colorbar" include NumRu include NMath #---------------------------------------------- # 設定 def set_filename(num) numnum = 0 num = ARGV[numnum+1].to_i if numnum = ARGV.index("-num") $exp = "T39L48_eml" $exp = ARGV[numnum+1] if numnum = ARGV.index("-exp") # $latitude = ARGV[numnum+1].to_f if numnum = ARGV.index("-lat") if $exp =~ /T159/ $latitude = - 0.4 $lat_label = 0 elsif $exp =~ /T39/ $latitude = - 1.5 $lat_label = 0 end $fig_path = "/home/yukiko/work/ape/GrADS/fig/" $dir_path = "/home/yukiko/work/ape/GrADS/#{$exp}_expID01/hourly/netcdf/" set_varinfo if num == 1 varlist = ["rain","q","t"] $sigma = [0.85,0.57] $var_array.push($rain_hash) $var_array.push($q_hash) $var_array.push($t_hash) $gifname = "#{$exp}_control_#{varlist[0]}-#{varlist[1]}-#{varlist[2]}_y#{$lat_label}" elsif num == 2 varlist = ["rain","q","t"] $sigma = [0.85,0.57] $var_array.push($rain_nolog_hash) $var_array.push($q_hash) $var_array.push($t_hash) $gifname = "#{$exp}_control_#{varlist[0]}-#{varlist[1]}-#{varlist[2]}_y#{$lat_label}_nologcon" elsif num == 3 varlist = ["rain","u","gph"] $sigma = [0.25,0.25] $var_array.push($rain_hash) $var_array.push($u_hash) $var_array.push($gph_hash) $gifname = "#{$exp}_control_#{varlist[0]}-#{varlist[1]}-#{varlist[2]}_y#{$lat_label}" end print "#{$gifname} name start #{Time.now}\n" end def set_varinfo rmiss = DCL.glpget('rmiss') $var_array = [] # こうすい $rain_hash = { "variable" => "tr_tppn", "diff" => true, "levels" => NArray[10.0**(-6), 10.0**(-5.5), 10.0**(-5), 10.0**(-4.5), 10.0**(-4), 10.0**(-3.5), DCL.glpget('rmiss')], "patterns" => NArray[35999, 45999, 55999, 70999, 75999, 85999], "lost_axis" => "(diff) from (mean) zonal", "title" => "precipitation flux (kg m-2 s-1)" } $rain_nolog_hash = { "variable" => "tr_tppn", "diff" => true, "levels" => NArray[0.0001, 0.0002, 0.0003, 0.0004, 0.0005, 0.0006, DCL.glpget('rmiss')], "patterns" => NArray[35999, 45999, 55999, 70999, 75999, 85999], "lost_axis" => "(diff) from (mean) zonal", "title" => "precipitation flux (kg m-2 s-1)" } # 温度偏差 $t_hash = { "variable" => "T", "diff" => true, "levels" => NArray[rmiss,-1,-0.5, 0, 0.5, 1, 1.5, rmiss], "lost_axis" => "(diff) from (mean) zonal", "title" => "temperature (K)" } # 比湿偏差 $q_hash = { "variable" => "Q", "diff" => true, "levels" => NArray[rmiss,-2e-3,-1e-3, 0, 1e-3, 2e-3, 3e-3,rmiss], "lost_axis" => "(diff) from (mean) zonal", "title" => "specific humidity (1)" } # 鉛直流 $sigdot_hash = { "variable" => "OMG", "diff" => false, "levels" => NArray[rmiss,-10e-7,-5e-7, 0, 5e-7, 10e-7, 15e-7, rmiss], "lost_axis" => " ", "title" => "sigma-vertical velocity (s-1)" } # ジオポテンシャル $gph_hash = { "variable" => "Z", "diff" => true, "levels" => NArray[rmiss,-10,-5, 0, 5, 10, 15, rmiss], "lost_axis" => "(diff) from (mean) zonal", "title" => "geopotential height (m)" } # 東西風 $u_hash = { "variable" => "U", "diff" => true, "levels" => NArray[rmiss,-5,-2.5, 0, 2.5, 5, 7.5, rmiss], "lost_axis" => "(diff) from (mean) zonal", "title" => "eastward wind (m s-1)" } # 南北風 $v_hash = { "variable" => "V", "diff" => true, "levels" => NArray[rmiss,-4,-2, 0, 2, 4, 6, rmiss], "lost_axis" => "(diff) from (mean) zonal", "title" => "northward wind (m s-1)" } # ----------------------- $day_hash = { "T39L48_eml"=> 40, "T39L48_ias"=> 45, "T39L48_ksc"=> 40, "T39L48_kuo"=> 15, "T39L48_mca"=> 40, "T39L48_non"=> 20, "T159L48_non"=> 60, "T159L48_eml"=> 25 } end # ---------------------------------------------- # お絵描きメイン def mradl_anim_main # ------------------------------------------------------------ # こうすい gphys_1 = GPhys::NetCDF_IO.open("#{$dir_path}#{$exp}_expID01_0001.nc", $var_array[0]["variable"]) gphys_2 = GPhys::NetCDF_IO.open("#{$dir_path}#{$exp}_expID01_0001.nc", $var_array[1]["variable"]) gphys_3 = GPhys::NetCDF_IO.open("#{$dir_path}#{$exp}_expID01_0001.nc", $var_array[2]["variable"]) # ------------------------------------------------------------ DCL.swpset('LDUMP', true) DCL.swpset('LWAIT',false) ; DCL.swpset('LWAIT0',false) DCL.swpset('LWAIT1',false) DCL.swpset('IPOSX', 50) ; DCL.swpset('IPOSY',50) DCL.gropn(1) DCL.uzfact(0.6) # 文字の大きさ倍 DCL.sgpset('lcorner',false) # コーナーマークはつけない DCL.slmgn(0.0, 0.0, 0.0, 0.0) # 描画マージン指定 DCL.sgpset('lfprop',true) # プロポーショナルフォントを使う DCL.sgpset('lfull',true) # 全画面表示 DCL.sgpset('lcntl', false) # アンダーバーは普通に表示 DCL.uscset('cyspos', 'B' ) # y 軸単位を書く位置 ipara = NumRu::DCL.uzpget('RSIZEC1') DCL.uzpset('RSIZEC0',0.012) # 小タイトル 文字大きさ DCL.uzpset('INDEXL0',1) # 小タイトル 文字ラインインデックス rmiss = DCL.glpget('rmiss') colorbar_hash = { 'vlength'=>0.13, "vwidth"=>0.01, 'labelintv'=>1, 'voff'=>0.02, "constwidth" => true } patterns = NArray[30999, 35999, 45999, 55999, 70999, 75999, 85999] # 1 枚目: 凝結加熱 $file_label = "lat=#{$latitude}" DCL.uzpset('labelxb',false) GGraph.set_fig('viewport'=>[0.10,0.80,0.52,0.65],'new_frame'=>true) GGraph.next_axes('xunits'=>'','yunits'=>'','xtitle'=>'') tone_hash = { 'patterns'=> $var_array[0]["patterns"], 'levels'=>$var_array[0]["levels"] } g = gphys_make(gphys_1,$latitude).lon_lotate g = g.set_att("long_name","").set_att("unit","").set_lost_axes("") GGraph.tone( g, true , tone_hash ) DCL::uxsttl("t", $var_array[0]["title"], -1.0 ) DCL::uzinit # DCL.uxpttl("t", 0, $var_array[0]["lost_axis"], 1.0) DCL::uxsttl("t", "#{$file_label}" , 1.0 ) GGraph.color_bar( colorbar_hash ) DCL::uzinit # 2 枚目: 温度偏差 sigma = $sigma[0] $file_label = "lat=#{$latitude}, sigma=#{sigma}" DCL.uzpset('labelxb',false) GGraph.next_axes('xunits'=>'','yunits'=>'','xtitle'=>nil) GGraph.set_fig('viewport'=>[0.10,0.80,0.34,0.47], 'new_frame'=>false ) tone_hash = { 'patterns'=> patterns, 'levels'=>$var_array[1]["levels"] } g = gphys_make(gphys_2,sigma).lon_lotate g_tmp = g. mean(0).set_lost_axes("") g = g.set_att("long_name","").set_att("unit","").set_lost_axes("") GGraph.tone( g-g_tmp, true , tone_hash ) DCL::uxsttl("t", $var_array[1]["title"], -1.0 ) DCL::uzinit DCL::uxsttl("t", "#{$file_label}" , 1.0 ) # DCL.uxpttl("t", 0, $var_array[1]["lost_axis"], 1.0) GGraph.color_bar( colorbar_hash ) DCL::uzinit # 3 枚目: 比湿偏差 sigma = $sigma[1] $file_label = "lat=#{$latitude}, sigma=#{sigma}" DCL.uzpset('labelxb',true) GGraph.set_fig('viewport'=>[0.10,0.80,0.16,0.29], 'new_frame'=>false ) tone_hash = { 'patterns'=> patterns, 'levels'=>$var_array[2]["levels"] } g = gphys_make(gphys_3,sigma).lon_lotate g_tmp = g. mean(0).set_lost_axes("") g = g.set_att("long_name","").set_att("unit","").set_lost_axes("") GGraph.tone( g-g_tmp, true , tone_hash ) DCL::uzinit DCL::uxsttl("t", $var_array[2]["title"], -1.0 ) DCL::uzinit DCL::uxsttl("t", "#{$file_label}" , 1.0 ) # DCL.uxpttl("t", 0, $var_array[2]["lost_axis"], 1.0) GGraph.color_bar( colorbar_hash ) DCL::grcls end # gif アニメ作成シェルスクリプト def mkgifanim `for file in *.xwd ; do xwdtopnm ${file} | pnmcut 2 2 904 654 > tmp.pnm; ppmtogif tmp.pnm > ${file}.gif; done >& /dev/null ` `mv dcl_001.xwd.gif #{$fig_path}#{$gifname}.gif` `rm dcl_* tmp.pnm` end # ---------------------------------------------- # GPhys メソッド追加定義 def GGraph.annotate(str_ary) # GGraph の annotate メソッド再定義 charsize = 0.7 * 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 = 0.66 # vy = 0.12 - charsize/2 # vy = 0.10 - charsize/2 vy = 0.11 - charsize/2 str_ary.each{|str| DCL::sgtxzv(vx,vy,"[ #{str} ]",charsize,0,-1,1) vy -= dvy } nil end class GPhys # 座標, データ操作ログ記録関数の定義 (Grid.lost_axes 参照) def add_lost_axes( lost ) GPhys.new( @grid.copy.add_lost_axes( lost.to_a ), @data) end def set_lost_axes( lost ) GPhys.new( @grid.copy.set_lost_axes( lost.to_a ), @data) end # rename 再定義 (gphys class の値を返す) def rename(name) GPhys.new(@grid,@data.copy.rename(name)) end # attribute 変更再定義 (gphys class の値を返す) def set_att(name,val) GPhys.new(@grid,@data.copy.set_att(name,val)) end # 経度 180 度回転 def lon_lotate # 軸の回転 (180 度引き算) # lon = VArray.new(@grid.axis("lon").pos.val-180).rename("lon") lon = VArray.new(@grid.axis(0).pos.val-180).rename("lon"). put_att("units",@grid.axis(0).pos.get_att("units")) lon = Axis.new().set_pos(lon) @grid = @grid.change_axis(0, lon) # lat = VArray.new(@grid.axis(1).pos.val).rename("lat") # lat = Axis.new().set_pos(lat) # @grid = @grid.change_axis(1, lat) # データの回転 dim = @data.val.shape[0] data_val = @data.copy.val data_val[0..dim/2-1, false] = @data.val[dim/2..dim-1, false] data_val[dim/2..dim-1, false] = @data.val[0..dim/2-1, false] data = @data.copy data.val = (data_val) # gphys class の値を返す GPhys.new( @grid, data ) end end def gphys_make(g,dimcut) grid0size = g.grid_copy.axis(0).pos.val.size grid_0 = g.grid_copy.axis(0) timesize = 5*24*3 +1 grid_time = Axis.new(). set_pos(VArray.new(NArray.sfloat(timesize).indgen!/24+ $day_hash[$exp]).rename("time"). put_att("units","days")) data = NArray.sfloat(grid0size,timesize).fill!(0.0) grid = Grid.new(grid_0, grid_time) data[true,1..120] = GPhys::NetCDF_IO.open( "#{$dir_path}#{$exp}_expID01_0001.nc", g.name).cut(true,dimcut,true).data.val.to_na data[true,121..240] = GPhys::NetCDF_IO.open( "#{$dir_path}#{$exp}_expID01_0002.nc", g.name).cut(true,dimcut,true).data.val.to_na data[true,241..360] = GPhys::NetCDF_IO.open( "#{$dir_path}#{$exp}_expID01_0003.nc", g.name).cut(true,dimcut,true).data.val.to_na data = VArray.new(data) #.rename("tr_tppn") gphys = GPhys.new(grid, data) return gphys end