#!/usr/bin/env ruby

load "/work11/ape/yukiko/lib/ape-view.rb"

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

END{

$rezol = ["T39L48_eml",
  "T39L24_eml",
  "T39L96_eml",
  "T79L48_eml" ]

  $t= ape_new "/work11/ape/NetCDF/#{$rezol[0]}/AGUforAPE-03a_TR_control.nc"
  $t.gropn(2)
  'rm dcl*.ps'

  $rezol.each{ |rezol|
    
    $t= ape_new "/work11/ape/NetCDF/#{rezol}/AGUforAPE-03a_TR_control.nc"
    g = $t.go("tr_tppn")
    gs = g.cut(true,0,true)[true,-401..-1].stspct_fft("#{g.name}_spct").
      set_att("ape_name","precipitation flux S-T power spectrum").
      set_att("units","(kg m-2 s-1 day-1)2")
    dim1 = gs.coord(0).shape.to_s.to_i
    dim2 = gs.coord(1).shape.to_s.to_i
    $file_label= "#{rezol}_control"
    $t.plot_spct_bunsan(gs[((dim1+1)/2-31)..((dim1+1)/2+29),0..80])
    file_name = "#{$file_label}_#{gs.name}_wk1999.eps"
    `dclpsrmcm dcl*.ps | dclpsrot > tmp.ps`
    `eps2eps tmp.ps  #{$fig_path}#{file_name}` 
    `rm tmp.ps dcl*.ps`
  }


  $rezol = ["T159L48_eml","T319L48_eml"]
  
  $rezol.each{ |rezol|

    $gr2ncfile_path ="/work11/ape/yukiko/data/#{rezol}_expID01/"

    file_name = Array.new
    Dir.foreach($gr2ncfile_path) { |item|
      file_name.push("#{$gr2ncfile_path}/#{item}") if item =~ /.nc$/
    }
    $t = Ape.new(file_name)
    
    g = $t.go("tr_tppn")
    gs = g.cut(true,0,true)[true,-401..-1].stspct_fft("#{g.name}_spct").
      set_att("ape_name","precipitation flux S-T power spectrum").
      set_att("units","(kg m-2 s-1 day-1)2")
    dim1 = gs.coord(0).shape.to_s.to_i
    dim2 = gs.coord(1).shape.to_s.to_i
    $file_label= "#{rezol}_control"
    $t.plot_spct_bunsan(gs[((dim1+1)/2-31)..((dim1+1)/2+29),0..80])
    file_name = "#{$file_label}_#{gs.name}_wk1999.eps"
    `dclpsrmcm dcl*.ps | dclpsrot > tmp.ps`
    `eps2eps tmp.ps  #{$fig_path}#{file_name}` 
    `rm tmp.ps dcl*.ps`
  }
  
  $t.grcls
  
}



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

class Ape

  def plot_spct_bunsan(gphys)

    gropn(1) if $gropn == nil 
    plot_set(gphys)

#    $fig_set_hash = { "window" => 
#      [gphys.axis(0).pos.val.min, gphys.axis(0).pos.val.max, 0.0 ,0.5] }
#
#     GGraph.set_fig($fig_set_hash)
    
    # attribute の long_name を消す (タイトル描画位置の都合)
    if gphys.data.get_att("long_name")
      gphys = gphys.set_att("long_name","")
    end

      GGraph.tone( gphys, true , $tone_hash ) 
#      GGraph.contour( gphys, false, $cont_hash ) unless $cont_flag == false


    # 分散曲線
    x_size = gphys.grid_copy.axis(0).pos.val.size
    x = gphys.grid_copy.axis(0).pos.val[x_size/2..-1]
    h_array = NArray.to_na([ 12, 25, 50, 100, 200 ])

    omega = 2*3.1415/24/60/60
    beta =  2*(omega)/(4e+7) * ( 60*60*24*(6370e+3) )  # (2*Omega/a)*cos(0)


    for i in (0..4)
      
      h = h_array[i]
      c_g = ((h*9.8)**0.5)*60*60*24/(4e+7)    # 重力波速度
      
      kelvin = x*c_g
      rossby = -x/(x**2/beta +(2+1)/c_g )
      grav = -(c_g*( c_g*x**2 + beta*(2+1) ) )**0.5

      DCL::uuslnt(1)
      ii=(i)*20+4
      ii = 5*20+4 if ii == 2*20+4
      ii = 720+4 if ii == 4*20+4
      DCL::uuslni(ii)
      DCL::uulin(x, kelvin)
      DCL::uulin(-x, -rossby)
      DCL::uulin(-x, -grav)
      DCL::uulin(x, -grav)

    end

    name_ary = ["h = 12 ", "h = 25 ","h = 50 ","h =100", "h =200 "]
    line_index_ary = [4,24,104,64,724]

    __line_index(name_ary,line_index_ary) 

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

    # 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

    # トーンバー
    unless gphys.axnames.size == 1
      DCL::Util::color_bar($cbar_conf)  unless $tone_flag == false
    end

  end  
  
  def __line_index(str_ary,line_index_ary) 
    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.15
    vy = 0.12 - 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 == 2
        vx = 0.30
        vy = 0.12 - charsize/2
      end
    }
    nil
  end

end