#! /usr/bin/ruby
# -*- coding: euc-jp -*-
#
# 2013/05/26 石渡正樹
#
# = fig_heatflux-lon-dependence.rb
# * 南北平均熱フラックスの図
#
# == USAGE
#   % fig_heatflux-lon-dependence.rb --time_start=931 --time_end=1096 --ragne=0:400 lat=mean --wsn=2
#   % fig_heatflux-lon-dependence.rb --time_start=931 --time_end=1096 --ragne=0:400 lat=0 --wsn=2

require 'getoptlong'
require "numru/ggraph"
include NumRu

## 定数設定
latent_heat = 2.5e6

## オプション解析
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

wsn = ($OPT_wsn||4)

lat = ($OPT_lat||0)

time_start = ($OPT_time_start||1).to_f
time_end = ($OPT_time_end||9999).to_f

range = ($OPT_range||'0:400')
min = range.split(':')[0].to_f
max = range.split(':')[1].to_f

## データ読み込み, 加工
lat_weight = GPhys::IO.open('OLR.nc', 'lat_weight').val
nlat = lat_weight.size

olr = GPhys::IO.open('OLR.nc', 'OLR')
olr = olr.cut('time'=>time_start..time_end).mean('time')

slr = GPhys::IO.open('SLR.nc', 'SLR')
slr = slr.cut('time'=>time_start..time_end).mean('time')

rain = GPhys::IO.open('Rain.nc', 'Rain')
rain = rain.cut('time'=>time_start..time_end).mean('time')

evap = GPhys::IO.open('Evap.nc', 'Evap')
evap = evap.cut('time'=>time_start..time_end).mean('time')

sens = GPhys::IO.open('Sens.nc', 'Sens')
sens = sens.cut('time'=>time_start..time_end).mean('time')

if lat == 'mean' then
  olr = 0.5 * (olr*lat_weight.reshape(1, nlat)).sum('lat')
  slr = 0.5 * (slr*lat_weight.reshape(1, nlat)).sum('lat')
  rain = latent_heat * 0.5 * (rain*lat_weight.reshape(1, nlat)).sum('lat')
  evap = 0.5 * (evap*lat_weight.reshape(1, nlat)).sum('lat')
  sens = 0.5 * (sens*lat_weight.reshape(1, nlat)).sum('lat')
else
  lat = lat.to_f

  olr = olr.cut('lat'=>lat)
  slr = slr.cut('lat'=>lat)
  rain = latent_heat * rain.cut('lat'=>lat)
  evap = evap.cut('lat'=>lat)
  sens = sens.cut('lat'=>lat)
end

## 作図
DCL.gropn($OPT_wsn||4)

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

GGraph.set_fig( 'itr'=> 1, 'viewport'=>[0.2,0.8,0.3,0.6] )

GGraph.line( olr, true, "min"=>min, "max"=>max, "index"=>855, "annotate"=>false, "exchange"=>false, "title"=>"Heat flux")

GGraph.line( slr, false, "index"=>745 )
GGraph.line( rain, false, "index"=>255 )
GGraph.line( evap, false, "index"=>555 )
GGraph.line( sens, false, "index"=>405 )

DCL.grcls