#!/usr/bin/env ruby

# ----------------------------------------------
# 必要なライブラリ, モジュールの読み込み
require "numru/ggraph"
require "numru/dcl"
include NumRu

# ----------------------------------------------
# ファイル名, 変数の読み込み

filename = "wavecisk.nc"
varname = "eigen_val"
eig = GPhys::NetCDF_IO.open(filename,varname)
varname = "growth"
growth = GPhys::NetCDF_IO.open(filename,varname)


# ----------------------------------------------
# 窓開く

iws = (ARGV[0] || (puts ' WORKSTATION ID (I)  ? ;'; DCL::sgpwsn; gets)).to_i
DCL::gropn(-iws)
GGraph.set_fig('viewport'=>[0.13,0.83,0.23,0.58],
	       "window" => [0.0, 4.0, -4.0, 4.0] ) 
DCL.slmgn(0.0, 0.0, 0.0, 0.0)  # 描画マージン指定
DCL.sgpset('lfull',true)       # 全画面表示
DCL.uzfact(0.6)
DCL.sgpset('lcorner',false)
DCL.uscset('cyspos', 'B' )      # y 軸単位を書く位置
DCL.sgpset('lfprop',true)      # プロポーショナルフォントを使う

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


END{

  plot_main(eig)  

  DCL::uusmki(2)
  markk(164,eig) #kelvin
  DCL::uusmki(22)

  markk(158,eig) #kelvin
#  markk(161,eig) #kelvin


  DCL::uxmttl("t", "2layer equatorial wave", -1.0 ) 

  plot_main(growth)

#  markk(113,growth)
#  DCL::uusmki(22)
#  markk(114,growth)


  DCL.grcls
  
}


def plot_main(eig)
  DCL.grfrm 

  k = eig.grid_copy.axis(0).pos.val
  eig = eig.val
  eig.shape[1].times{ |num|
    DCL.usspnt(k,eig[true,num])
  }
  
  NumRu::DCL.grswnd(0.0,4.0,-4.0,4.0)
  NumRu::DCL.grsvpt(0.23,0.58,0.13,0.83)
  NumRu::DCL.grstrn(1)
  
  DCL::grstrf
  DCL::ussttl( 'knum', '1', 'freq', '1' )
  DCL::usdaxs
  DCL::uusmks(0.015)
  DCL::uusmkt(1)
  
  (eig.shape[1]).times{ |num|
    DCL::uusmki(100+num*20+5)
    DCL.uumrk(k,eig[true,num])
  }
  
end

def markk(eigennum,eig)
#  DCL::uusmki(2)
  DCL::uusmkt(4)
  DCL::uusmks(0.01)
  puts eig.shape[1]
  a = [0.5,0.5]
  aa = [eig.val[4,eigennum],eig.val[4,eigennum]]
  DCL.uumrk(a,aa)
end


# lm=81
#DCL.uumrk(k[4..5],eig[4..5,327]) # ロスビー
# DCL.uumrk(k[4..5],eig[4..5,326]) # ロスビー
# DCL.uumrk(k[4..5],eig[4..5,325]) # 南北壁のケルビン1
# DCL.uumrk(k[4..5],eig[4..5,324]) # ロスビー1
# DCL.uumrk(k[4..5],eig[4..5,323]) # 赤道ケルビン1
# DCL.uumrk(k[4..5],eig[4..5,322]) # 混合1
# DCL.uumrk(k[4..5],eig[4..5,321]) # 南北壁のケルビン2
# DCL.uumrk(k[4..5],eig[4..5,320]) # n=0 東進重力2
# DCL.uumrk(k[4..5],eig[4..5,319])  # 赤道ケルビン2
# DCL.uumrk(k[4..5],eig[4..5,318])  # n=1 西進重力 2
# DCL.uumrk(k[4..5],eig[4..5,317])  # 南北壁のケルビン1
# DCL.uumrk(k[4..5],eig[4..5,316])  # 混合2
# DCL.uumrk(k[4..5],eig[4..5,315])  # n=1 東進重力1
# DCL.uumrk(k[4..5],eig[4..5,314])  # n=1 西進重力1
# DCL.uumrk(k[4..5],eig[4..5,313])  # n=2 東進重力1
# DCL.uumrk(k[4..5],eig[4..5,312])  # n=3 西進重力1
# DCL.uumrk(k[4..5],eig[4..5,311])  # n=0 東進重力2
# DCL.uumrk(k[4..5],eig[4..5,310])  # n=3 東進重力1


# 南北壁のケルビン1
# DCL.uumrk(k[4..5],eig[4..5,161]) # 南北壁のケルビン1
# DCL.uumrk(k[4..5],eig[4..5,162]) # 混合2
# DCL.uumrk(k[4..5],eig[4..5,163]) # 赤道ケルビン1
# DCL.uumrk(k[4..5],eig[4..5,164])  # 南北壁のケルビン2
# DCL.uumrk(k[4..5],eig[4..5,165]) # n=1 ロスビー
# DCL.uumrk(k[4..5],eig[4..5,159]) # 赤道ケルビン2
# DCL.uumrk(k[4..5],eig[4..5,158]) # n=0 東進重力2
# DCL.uumrk(k[4..5],eig[4..5,157]) # 混合 ? 壁ケルビン ? 
# DCL.uumrk(k[4..5],eig[4..5,156]) # n=1 西進重力
# DCL.uumrk(k[4..5],eig[4..5,155]) # n=1 東進重力2
# DCL.uumrk(k[4..5],eig[4..5,154]) # n=0 東進重力1

# lm=41
# DCL.uumrk(k[4..5],eig[4..5,160]) # 南北壁のケルビン1
# DCL.uumrk(k[4..5],eig[4..5,161]) # 南北壁のケルビン1
# DCL.uumrk(k[4..5],eig[4..5,162]) # 混合2
# DCL.uumrk(k[4..5],eig[4..5,163]) # 赤道ケルビン1
# DCL.uumrk(k[4..5],eig[4..5,164])  # 南北壁のケルビン2
# DCL.uumrk(k[4..5],eig[4..5,165]) # n=1 ロスビー
# DCL.uumrk(k[4..5],eig[4..5,159]) 
# 赤道ケルビン2
# DCL.uumrk(k[4..5],eig[4..5,158]) # n=0 東進重力2
# DCL.uumrk(k[4..5],eig[4..5,157]) # 混合 ? 壁ケルビン ? 
# DCL.uumrk(k[4..5],eig[4..5,156]) # n=1 西進重力
# DCL.uumrk(k[4..5],eig[4..5,155]) # n=1 東進重力2
# DCL.uumrk(k[4..5],eig[4..5,154]) # n=0 東進重力1