# -*- coding: utf-8 -*-
require "numru/ggraph"
include NumRu

# ----- コマンドラインオプション -----
tro_file = ARGV[0].to_s         # 第1引数はスペクトルを取り出してくるファイル名(必須)
cli_file = ARGV[1].to_s         # 第1引数はスペクトルを取り出してくるファイル名(必須)
iws = ( ARGV[2] || 1 ).to_i     # 装置番号 (1,2,or 4)
tsleep = [ ( ARGV[3] || 0.0 ).to_f, 0.0 ].max  # 第5引数は描画間隔(秒)
wait = ( tsleep <= 0.0 )            #->true/false; 0以下(true)ならマウスクリックを待つ

# ----- 使用データ -----
tro_1deksp = GPhys::IO.open tro_file, 'tro_1deksp'
cli_1deksp = GPhys::IO.open cli_file, 'cli_1deksp'

tro_1deksp = tro_1deksp.copy
tro_1deksp.name = "barotropic"
cli_1deksp = cli_1deksp.copy
cli_1deksp.name = "baroclinic"
total = tro_1deksp + cli_1deksp
total.name = "total"

# k^{-5/3} の直線を描く
x1 = NArray.sfloat(150)
y1 = NArray.sfloat(150)
x1[0] = 20.0
y1[0] = 1.0*10**(-1)
for k in 0..148
  x1[k+1] = x1[k] + 1.0
  y1[k+1] = (x1[k+1]/x1[0])**(-5/3)*y1[0]
end
# k^{-3} の直線を描く
x2 = NArray.sfloat(150)
y2 = NArray.sfloat(150)
x2[0] = 20.0
y2[0] = 1.0*10.0**(-1)
for k in 0..148
  x2[k+1] = x2[k] + 1.0
  y2[k+1] = (x2[k+1]/x2[0])**(-3)*y2[0]
end

ytitle = 'E(k)'

# ----- DCL設定 -----
def prep_dcl(iws=1,wait=false) # iws : DCL出力デバイス．1,4:画面, 2:PS
  DCL.swpset("iwidth",700)     # 画面の幅
  DCL.swpset("iheight",600)    # 画面の高さ
  DCL.swpset("lwait",wait)     # 次の描画の前にマウスクリックを待つ
  DCL.swpset("lalt",true)      # 裏で描画（パラパラアニメ用）
  DCL.swpset("ifl",1)          # png で出力
#  DCL.swpset("lsysfnt",true)   # システムフォントを使用する
#  DCL.swcset("fontname","Hiragino Kaku Gothic")
  DCL.sgiset("ifont",2)
  DCL.sgscmn(39)               # カラーマップ番号の指定
  DCL.gropn(iws)
  DCL.sgpset('lfull',true)     # 全画面表示
  DCL.sglset("lcntl",true)     # 上付き下付き文字を有効
#  DCL.sgpset('isub', 96)       # 下付き添字を表す制御文字を '_' から '`' に
  DCL.glpset('lmiss',true)     # DCLの欠損値処理を on に
end

title = DCL.csgi(131) + DCL.csgi(164) + "= 0"
ytitle = 'E(k)_tot",' + 'E(k)_' + DCL.csgi(150) + '", E(k)_' + DCL.csgi(170) + '"'
prep_dcl(iws,wait)
t0=1800
t1=2000
# ----- エネルギースペクトルの図 -----
GGraph.set_fig("itr"=>4,"viewport"=>[0.15,0.72,0.15,0.77])
GGraph.set_axes("xtitle"=>"k","ytitle"=>ytitle,'xside'=>'b', 'yside'=>'l')
GGraph.line total.cut("t"=>t0..t1).mean("t"),true,"annotate"=>false,"title"=>"","max"=>5.0*10**0,"min"=>1.0*10**(-9),"legend"=>true,"type"=>1,"index"=>284
GGraph.line tro_1deksp.cut("t"=>t0..t1).mean("t"),false,"legend"=>true,"type"=>1,"index"=>862
GGraph.line cli_1deksp.cut("t"=>t0..t1).mean("t"),false,"legend"=>true,"type"=>1,"index"=>502
DCL.uulinz(x1,y1,3,2)
DCL.uulinz(x2,y2,4,2)

DCL.sgstxs(0.02)
DCL.sgtxu(260,1.5*10**(-3),'k|-5/3"')
DCL.sgtxu(240,1.5*10**(-4),'k|-3"')

DCL.grcls      # 描画終了処理