#!/usr/bin/env ruby

# ----------------------------------------------
# local load path

$local_path = '/work11/ape/yukiko/lib' 
$: << $local_path

# ----------------------------------------------
# 必要なライブラリ, モジュールの読み込み

load "#{$local_path}/ape-view.rb"
load "#{$local_path}/lib-ape-view.rb"
load "#{$local_path}/lib-ape-base.rb"
load "#{$local_path}/lib-ape-stspctrum.rb"

# --------------------------------------------------------------------


END{

  a = Ape_mkfig.new 3
  $expID = "control"
  $rezol = "T159L48_non"
  set_path
  a.gr_comp_point
  a.gr_comp_tuv
}

# ----------------------------------------------


def set_path
  $fig_path = "/work11/ape/yukiko/figs/tmp/"
  $gr2ncfile_path = "/work11/ape/yukiko/data/#{$rezol}_expID01"
  $file_label = "#{$rezol}_#{$expID}"
  $groupid = "AGUforAPE-03a"
end

# ----------------------------------------------


class Ape_mkfig

  def mkfig_plot_comp(gphys, gphys_comp)

    if @fig_num == 3 then
      @data.mkdump_comp(gphys, gphys_comp)
    else
      @data.plot_comp(gphys, gphys_comp)
    end

  end

  def gr_comp_point

    file_name = Array.new
    Dir.foreach($gr2ncfile_path) { |item|
      file_name.push("#{$gr2ncfile_path}/#{item}") if item =~ /.nc$/
    }
    @data = Ape.new(file_name)

    #  tr_tppn
    rain = @data.gphys_open("tr_tppn").cut(true,0,true).lon_lotate
    lost_axes = rain.lost_axes.to_s.sub("y=","lat=")
#    rain = rain[true,-400..-1].set_lost_axes(lost_axes)
    rain = rain.set_lost_axes(lost_axes)

    # t_conv composite
#    gphys = @data.gphys_open("T")[true,true,-400..-1].lon_lotate
    gphys = @data.gphys_open("T").lon_lotate
    comp = gphys.composite(rain, 0.0007, true)
    comp = comp.      
      set_att("ape_name",
	      "composite_point" ).
      set_att("units", "1" ).
      rename("comp_point").
      set_lost_axes(lost_axes)

    # 描画
#    mkfig_plot(rain)
#    Ape.plot_set(comp) 
#    mkfig_plot(comp)
    mkfig_plot_comp(rain, comp)
#    GGraph.contour( comp, false, $cont_hash )
   
  end

  def gr_comp_tuv

    file_name = Array.new
    Dir.foreach($gr2ncfile_path) { |item|
      file_name.push("#{$gr2ncfile_path}/#{item}") if item =~ /.nc$/
    }
    @data = Ape.new(file_name)

    #  tr_tppn
    rain = @data.gphys_open("tr_tppn").cut(true,0,true).lon_lotate
    lost_axes = rain.lost_axes.to_s.sub("y=","lat=")
#    rain = rain[true,-400..-1].set_lost_axes(lost_axes)
    rain = rain.set_lost_axes(lost_axes)

    # T composite
    #    gphys = @data.gphys_open("T")[true,true,-400..-1].lon_lotate
    gphys = @data.gphys_open("T").lon_lotate
    t_comp = gphys.composite(rain, 0.0007)
    t_comp = t_comp.
      set_att("ape_name", "T_(U,_-OMG)_composite").
      set_att("units", "K, (m s-1, Pa s-1)" )
    t_comp = (t_comp - t_comp.mean(0)).
      set_lost_axes(lost_axes).
      add_lost_axes("(diff) from (mean) zonal").
      rename("comp_tuv")

    # U composite
#    gphys = @data.gphys_open("U")[true,true,-400..-1].lon_lotate
    gphys = @data.gphys_open("U").lon_lotate
    u_comp = gphys.composite(rain, 0.0007)
    u_comp = u_comp - u_comp.mean(0)

    # OMG composite
#    gphys = @data.gphys_open("OMG")[true,true,-400..-1].lon_lotate
    gphys = @data.gphys_open("OMG").lon_lotate
    om_comp = gphys.composite(rain, 0.0007)    
    om_comp = om_comp - om_comp.mean(0)

    # ベクトルの間引
    u_data  = NArray.sfloat( u_comp.axis(0).pos.val.size , 
			    u_comp.axis(1).pos.val.size ).fill!(0.0)
    
    om_data = NArray.sfloat( om_comp.axis(0).pos.val.size , 
			    om_comp.axis(1).pos.val.size ).fill!(0.0)
    
    num_l = u_comp.axis(0).pos.val.size/120
    num_z = u_comp.axis(1).pos.val.size/48
    
    u_comp.axis(0).pos.val.size.times{ |num|
      if (num % num_l) == 0
	u_data[num,true]  =  u_comp.data.val[num,true]
	om_data[num,true] =  om_comp.data.val[num,true]
      elsif
	u_data[num,true]  = 0.0
	om_data[num,true] = 0.0
      end
    }
    
    num_z =1 if num_z < 1

    u_comp.axis(1).pos.val.size.times{ |num|
      if (num % num_z) == 0
	u_data[true,num]  =  u_data[true,num]
	om_data[true,num] =  om_data[true,num]
      elsif
	u_data[true,num]  = 0.0
	om_data[true,num] = 0.0
      end
    }


    # 軸の long_name 直し (pressure level => sigma)
    z = t_comp.grid_copy.axis(1).pos.
      set_att("long_name","sigma").
      set_att("units","1")
    z = Axis.new().set_pos(z)
    grid = t_comp.grid_copy.change_axis(1, z)
    t_comp = GPhys.new( grid, t_comp.data )

    # 描画
    dim = t_comp.axis(0).pos.val.size
    mkfig_plot(t_comp[(dim/3)..(dim*2/3),true],
	       u_data[(dim/3)..(dim*2/3),true], 
	       -om_data[(dim/3)..(dim*2/3),true])

  end

end


class Ape

  def plot_comp(gphys, gphys_u=nil, gphys_v=nil)
    gropn(1) if $gropn == nil
    
    if gphys.class == Array
      plot_set(gphys[0]) 
    else 
      plot_set(gphys) 
    end

    plot_main_comp(gphys, gphys_u)

    # タイトルを標準出力へ
    if gphys.class == Array then
      if  gphys[0].data.get_att("ape_name") then
	puts "#{gphys[0].name} (#{gphys[0].data.get_att("ape_name")})"
      end
    elsif gphys.data.get_att("ape_name") then
      puts "#{gphys.name} (#{gphys.data.get_att("ape_name")})"
    else
      puts "#{gphys.name}"
    end
  end

  def mkdump_comp(gphys, gphys_u=nil, gphys_v=nil)
    gropn(3) if $gropn == nil

    plot_comp(gphys, gphys_u)

    gphys = gphys[0]     if gphys.class == Array
    
    file_name = "#{$file_label}_#{gphys.name}.gif"
    `rm dcl_* tmp.pnm`  
    plot(gphys)
    `for file in *.xwd;do xwdtopnm ${file} | pnmcut 2 2 904 654 > tmp.pnm;done`
    `ppmtogif tmp.pnm > #{$fig_path}#{file_name}`
    `rm dcl_* tmp.pnm`  

    $file_name = $pre_file

  end

  def plot_main_comp(gphys, gphys_u=nil)

    # attribute の long_name を消す (タイトル描画位置の都合)
    if gphys.data.get_att("long_name")
      gphys = gphys.set_att("long_name","")
    end

    GGraph.set_fig($fig_set_hash)

    # 描画
    if gphys.axnames.size == 1 || gphys.name =~ "gt_"
      
      GGraph.line( gphys, true, $line_hash )
    elsif  $tone_flag == false
      GGraph.contour( gphys, true, $cont_hash )
    elsif $tone_hash['levels'][0] == $tone_hash['levels'][1] 
      GGraph.tone( gphys, true ) 
      GGraph.contour( gphys, false ) unless $cont_flag == false
    else
      GGraph.tone( gphys, true , $tone_hash ) 
      GGraph.contour( gphys, false, $cont_hash ) unless $cont_flag == false
    end

    # タイトル
    tani = "(#{gphys.data.get_att("units")})"
    DCL::uxsttl("t", tani , -1.0 )
    if  gphys.data.get_att("ape_name") then
      DCL::uxmttl("t", gphys.data.get_att("ape_name").gsub("_", " "), -1.0 ) 
    end

    $file_label = "#{@file}@#{gphys.name}"  if $file_label == "filename"

    # nc ファイル名
    if $gropn == 2
      charsize = 0.79 - $file_label.size * 0.0115 
      DCL::sgtxzv(charsize,0.62,$file_label,0.8 * DCL.uzpget('rsizec1'),0,-1,1)
    else
      charsize = 0.81 - $file_label.size * 0.0095 
      DCL::sgtxzv(charsize,0.6,$file_label,0.8 * DCL.uzpget('rsizec1'),0,-1,1) 
    end


    plot_set(gphys_u) 
    GGraph.contour( gphys_u, false, $cont_hash )

    # トーンバー
    unless gphys.axnames.size == 1
      unless $tone_hash['levels'][0] == $tone_hash['levels'][1] 
	DCL::Util::color_bar($cbar_conf)  unless $tone_flag == false
      end
    end

  end

end