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

# 同期回転惑星計算用描画スクリプト
# == 説明
# * このスクリプトは水平平均した物理量の時系列をかくものである.
# * 下にある omegas に描画する Omega の値のリストを設定する.
#
# == オプション
# --time_start 開始時刻
# --time_end   終了時刻
# --region     平均領域の指定 (all: 全球, east: 東半球, west: 西半球)
# --range      プロットする値の範囲. min:max という形式で指定する.
# --expnum     アンサンブル実験の実験番号 (000 - 009)
# --wsn        出力先. ps だったら --wsn=2
#
# == USAGE
#   % OLR_timeseries.rb --time_start=0 --time_end=2000 --region=all --wsn=2
#
# == 履歴
# * 2013-6-19 石渡正樹 作成
# * 2014-04-22 石渡正樹 論文用の図の作成


require "getoptlong"        # for option_parse
require "numru/ggraph"
#require "../load_config.rb"

DATA_HOME = "/GFD_Dennou_Work10/momoko/SyncRotEarthRad-2"
DIR_HEADER = "SR_CLT1500_Omega"
DIR_FOOTER = "_T42L26"



TMPFILE = "dcl.ps"

## オプション解析
parser = GetoptLong.new
parser.set_options(
   ###    global option   ###
   ['--time_start',               GetoptLong::REQUIRED_ARGUMENT],
   ['--time_end',                 GetoptLong::REQUIRED_ARGUMENT],
   ['--region',                   GetoptLong::REQUIRED_ARGUMENT],
   ['--range',                    GetoptLong::REQUIRED_ARGUMENT],
   ['--expnum',                    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

# 変数の設定
vars = [
        {'name'=>'OLRA'}
       ]
#title = 'OLR'
title = ''
range = ($OPT_range||'320:360')

VAL_MIN =  range.split(':')[0].to_f
VAL_MAX =  range.split(':')[1].to_f

wsn = ($OPT_wsn||4)

flux_all = Hash.new
flux_east = Hash.new
flux_west = Hash.new

region = ($OPT_region||'all')

time_start = ($OPT_time_start||0).to_f
time_end = ($OPT_time_end||2000).to_f

#omegas = [0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.33, 0.5,
#          0.6, 0.67, 0.7, 0.75, 0.8, 0.85, 0.9, 1.0]
#omegas = [0.0, 0.15, 0.5, 1.0] 
omegas = [1.0] 

na_omegas = NArray.to_na(omegas)
nomega = omegas.size

# 線の色の設定
#index =[601,401,1,801]
index =[601]

# 各変数で回す
for v in vars

  var = v['name']
  epsfn = var + '_timeseries_' + region + '_' + '.eps'

  all  = NArray.sfloat(2001,nomega)
  west = NArray.sfloat(nomega)
  east = NArray.sfloat(nomega)

  ## 各 omega で回す.
  for i in 0..nomega-1

    omega = na_omegas[i]
    var.to_s + ' ' + omega.to_s

    ## データ読み込み, 加工
    dir = 'SR_CLT1500_Omega' + omega.to_s + '_T42L26'
    file = var + '.nc'
    path = File.join(DATA_HOME, dir, file)
    gphys = GPhys::IO.open(path, var)

    # 重みデータの読み込み
    na_lon = GPhys::IO.open(path, 'lon').val
    na_lat = GPhys::IO.open(path, 'lat').val
    na_lat_weight = GPhys::IO.open(path, 'lat_weight').val
    nlon = na_lon.size
    nlat = na_lat.size

    # 時間方向の切り出し
    gphys = gphys.cut('time'=>time_start..time_end)

    # 空間平均の計算
    gphys = gphys * na_lat_weight.reshape(1, nlat) / na_lat_weight.sum

    all = gphys.mean('lon').sum('lat')
    west = gphys.cut('lon'=>na_lon[0]..na_lon[nlon/2 -1]).mean('lon').sum('lat')
    east = gphys.cut('lon'=>na_lon[nlon/2]..na_lon[nlon-1]).mean('lon').sum('lat')

    # 軸のタイトル変える
    all.long_name = ('OLR')
    west.long_name = ('OLR')
    eastlong_name = ('OLR')

    all.units = ( 'W m|-2"')
    west.units = ( 'W m|-2"')
    east.units = ( 'W m|-2"')

## check
#    p all - (west + east)/2.0

    flux_all[i] = all
    flux_east[i] = east
    flux_west[i] = west

  end

  # 描画
  DCL.gropn(wsn)
    
  DCL.sgpset('lcntl', true) ; DCL.uzfact(0.7)
  GGraph.set_fig( 'itr'=> 1, 'viewport'=>[0.2,0.8,0.3,0.6] )
  DCL.sgpset("lclip", true)

  if region == 'all' then
    GGraph.line( flux_all[0] , true , 'index'=>index[0], 'max'=>VAL_MAX, 
                'min'=>VAL_MIN, 'annotate'=>false, 'title'=>title)
    for i in 1..nomega-1
      GGraph.line( flux_all[i] , false , 'index'=>index[i], 'max'=>VAL_MAX, 
                 'min'=>VAL_MIN, 'annotate'=>false, 'title'=>title)
    end

  elsif region == 'east' then
    GGraph.line( flux_east[0] , true , 'index'=>index[0], 'max'=>VAL_MAX, 
                 'min'=>VAL_MIN, 'annotate'=>false, 'title'=>title)
    for i in 1..nomega-1
      GGraph.line( flux_east[i] , false , 'index'=>index[i], 'max'=>VAL_MAX, 
                 'min'=>VAL_MIN, 'annotate'=>false, 'title'=>title)
    end

  elsif region == 'west' then
    GGraph.line( flux_west[0] , true , 'index'=>index[0], 'max'=>VAL_MAX, 
                 'min'=>VAL_MIN, 'annotate'=>false, 'title'=>title)
    for i in 1..nomega-1
      GGraph.line( flux_west[i] , false , 'index'=>index[i], 'max'=>VAL_MAX, 
                 'min'=>VAL_MIN, 'annotate'=>false, 'title'=>title)
    end

  else
    p 'region is invalid value!'
  end

  DCL.grcls

  if wsn == '2' then
    if File.exist?(TMPFILE) then
      puts("MESSAGE: after processing ... #{epsfn} will be generate.")
      `dclpsrmcm #{TMPFILE} | dclpsrot | dclpsmargin > #{epsfn}`
    end
  end

end