#!/usr/bin/ruby
# -*- coding: utf-8 -*-

require "numru/ggraph"
include NumRu
require "pp"

wsn = 4

# ベクトルの描画間隔 (単位は格子点数)
xintv = 6
yintv = 1

DATA_HOME = "."
#ITR=30   # 球
ITR=1

LON_CNT=125  # 球のときに使用
title=''
xfact1 = 1.0e-8
yfact1 = 1.0e-8
type = 1
cp = 1004.6

# 凡例の単位ベクトルの長さ (デフォルト値)
uxunit = 0.05
uyunit = 0.05

time_start = 725
time_end   = 1095

sig_start = 0.01
sig_end   = 0.8


sigm_start = 0.01
sigm_end   = 0.8


#time_start = 931
#time_end   = 1096

#omegas = [0.0, 0.15, 0.5, 1.0]
#omegas = [1.0]
#omegas = [0.0]

omega = 1.0

dir = nil


## オプション解析 new version ここから
## parser = GetoptLong.new
## parser.set_options(
##                    ###    global option   ###
##                    ['--lat',                      GetoptLong::REQUIRED_ARGUMENT],
##                    ['--time_start',                      GetoptLong::REQUIRED_ARGUMENT],
##                    ['--time_end',                      GetoptLong::REQUIRED_ARGUMENT],
##                    ['--range',                      GetoptLong::REQUIRED_ARGUMENT],
##                    ['--wsn',                      GetoptLong::REQUIRED_ARGUMENT]
##                    )
## parser.each_option do |name, arg|
##   eval "$OPT_#{name.sub(/^--/, '').gsub(/-/, '_')} = '#{arg}'"  # strage option value to $OPT_val
## end
# ここまで



# 簡易
while ARGV.size > 0
  case ARGV[0]
  when '--dir' then
    ARGV.shift
    dir = ARGV[0]
    ARGV.shift
  when '--wsn' then
    ARGV.shift
    wsn = ARGV[0]
    ARGV.shift
  when '--xintv' then
    ARGV.shift
    xintv = ARGV[0]
    ARGV.shift
  when '--yintv' then
    ARGV.shift
    yintv = ARGV[0]
    ARGV.shift
  when '--xfact1' then
    ARGV.shift
    xfact1 = ARGV[0]
    ARGV.shift
  when '--yfact1' then
    ARGV.shift
    yfact1 = ARGV[0]
    ARGV.shift
  when '--uxunit' then
    ARGV.shift
    uxunit = ARGV[0]
    ARGV.shift
  when '--uyunit' then
    ARGV.shift
    uyunit = ARGV[0]
    ARGV.shift
  when '--loncnt' then
    ARGV.shift
    lon_cnt = ARGV[0]
    ARGV.shift
  when '--range_time' then
    ARGV.shift
    tmp = ARGV[0].split(':')
    time_start = tmp[0].to_f
    time_end   = tmp[1].to_f
    ARGV.shift
  when '--range_val' then
    ARGV.shift
    tmp = ARGV[0].split(':')
    val_min = tmp[0].to_f
    val_max = tmp[1].to_f
    ARGV.shift
    val_interval = (val_max - val_min) / 30.0
  else
    raise "ERROR: invalid option: #{ARGV[0]}"
  end
end





if dir == nil then
  dir = File.join(DATA_HOME,"ExpOcean")
end

if type == 1 then
  puts "all"
elsif type == 2 then
  puts "disturbance"
else
  raise "Type #{type} is not supported."
end

#dir = File.join(DATA_HOME, 
#      "101220_diff-DryAir-Vapor_T21L16_sr_Omega" + omega.to_s + "_dt20m")

dir='.'


# geopotential 
temperature_file = dir + '/Temp.nc'
system("ruby /GFD_Dennou_Work11/masuda/script/geopotential_time.rb --var Temp.nc@Temp --range_time #{time_start}:#{time_end}")


# input
var = 'U'
path = var+'.nc'
gp_u = GPhys::IO.open(path, var).cut('time'=>time_start..time_end, 'lat'=>0, 'sig'=>sig_start..sig_end)

var = 'SigDot'
path = var+'.nc'
gp_sigdot = GPhys::IO.open(path, var).cut('time'=>time_start..time_end, 'lat'=>0, 'sigm'=>sigm_start..sigm_end)

var = 'QVap'
path = var+'.nc'
gp_qvap = GPhys::IO.open(path, var).cut('time'=>time_start..time_end, 'lat'=>0, 'sig'=>sig_start..sig_end)

var = 'Temp'
path = var+'.nc'
gp_h = GPhys::IO.open(path, var).cut('time'=>time_start..time_end, 'lat'=>0, 'sig'=>sig_start..sig_end)*cp

# Phi
gp_phi = GPhys::IO.open('Phi.nc', 'Phi').cut('time'=>time_start..time_end,'lat'=>0, 'sig'=>sig_start..sig_end)

# dry static energy
gp_q = gp_h + gp_phi

gp_u_q = gp_u * gp_q
gp_sigdot_q = gp_sigdot * gp_q

# 時間平均
gp_u_q = gp_u_q.mean('time')
gp_sigdot_q = gp_sigdot_q.mean('time')
if type == 2 then
  gp_u = gp_u.mean('time')
  gp_sigdot = gp_sigdot.mean('time')
  gp_q = gp_q.mean('time')
end

# 擾乱成分
if type == 2 then
  gp_ud_qd = gp_u_q - gp_u * gp_q
  gp_sigdotd_qd = gp_sigdot_q - gp_sigdot * gp_q
end

DCL.gropn(wsn)

DCL.sgpset('lcntl', false) ; DCL.uzfact(0.7)

GGraph.set_axes('xtickint'=>10, 'xlabelint'=>30)
GGraph.set_axes('ytickint'=>10, 'ylabelint'=>30)

# same as gpview
vx0,vx1,vy0,vy1 = 0.15,0.85,0.2,0.55
GGraph.set_fig 'itr'=>ITR, 'viewport'=>[vx0,vx1,vy0,vy1], 'map_axis'=>[lon_cnt, 0, 0]

DCL.uscset('cyspos', 'B' )              # move unit y axis

DCL.udlset('LMSG', false)

DCL.ugpset('XFACT1', xfact1)
DCL.ugpset('YFACT1', yfact1)
DCL.ugpset('UXUNIT', uxunit) # size of unit vector by dimentional value
DCL.ugpset('UYUNIT', uyunit)
DCL::ugrset('VXUNIT', 0.1)  # size of unit vector
DCL::ugrset('VYUNIT', 0.1)

DCL::uxmttl('T', '', 0.0)  # title (large character)

DCL.ugpset('LNRMAL', false)  # scale vector automatically
#DCL.ugpset('LMSG', false)   # no message
DCL.ugpset('LUMSG', false)   # no message

DCL.ugpset('INDEX', 5)
#DCL.ugpset('INDEX', 959)

opt_vector = {
  'xintv'=>xintv, 'yintv'=>yintv,
  'flow_vect'=>false, 'annotate'=>true, 'unit_vect'=>true
#  'flow_vect'=>true, 'annotate'=>false, 'unit_vect'=>true
}

opt_contour = {
  'min'=>val_min, 'max'=>val_max, 'interval'=>val_interval,
  'title'=>title, 'annotate'=>false
}


# 等値線とトーンのプロット
#DCL.sgscmn(9)   # colormap: white-blue-black
#DCL.sgscmn(8)   # colormap: white-yellow-red-black
#DCL.sgscmn(2)   # colormap: black-red-yellow-white

#GGraph.tone_and_contour( gp_qvap, true, opt_contour)

#GGraph.tone( gp_qvap, true)
#GGraph.color_bar('landscape'=>true)

#  GGraph.color_bar('landscape'=>true,'top'=>true, 'vlength'=>0.5, 'voff'=>0.03)

# ベクトルのプロット
#DCL.sgscmn(7)   # colormap: 
#DCL.sgscmn(77)   # colormap: 

if type == 1 then
  # mean
  GGraph.vector( gp_u_q, gp_sigdot_q, true, opt_vector)
elsif type == 2 then
  # disturbance
  GGraph.vector( gp_ud_qd, gp_sigdotd_qd, true, opt_vector)
else
  raise "Type #{type} is not supported."
end

DCL.grcls