=begin
= netcdf_cutter.rb
== 概要
:ORG_data に格納される ERA-40 オリジナルデータを 1 ヵ月毎に分割するスクリプト.

== USAGE
  $ netcdf_cutter.rb ファイル名 変数

現時点では以下のような難点がある.

* global attribute が消滅してしまう. 
=end

#############  以下メインルーチン  ##############
require "numru/netcdf"


module NumRu
  class NetCDFVar
    alias put scaled_put
    alias get scaled_get
  end
end

require "numru/ggraph"
include NumRu

require "getopts"

#======= ruby 1.6 で Date::parse をつかうためのおまじない
class Hash
  alias values_at indices # values_at は ruby-1.6 には無いメソッドなので
end
#=======

# 引数に与えた nc ファイルを格納する 配列.
file = Array.new

# 引数中の nc ファイルのみ格納
ARGV.each do |i| 
  if i =~ /.nc$/ then
    p i
    file << i
  end
end

# 変数 : デフォルトはヘッダ の一番最後の変数を自動的に選択
if ARGV[ARGV.size-1] =~ /.nc$/ then
  var = NetCDF.open(file[0], "r").var_names[-1]  
else # 陽に指定も可能. その場合, 一番最後に与える.
  var = ARGV[ARGV.size-1]
end

pname = var.to_s.upcase #=> 物理量を get.

#--------------------------------------
# title を得る. もしなければ, longname をそれとする.
nc = NetCDF.open(file[0])
p  longname = nc.var(var).att("long_name").get
p  title = nc.att("title")
nc.close

if title 
p  $title = title
else 
p  $title = "daily "+longname
end

#  gphys オープン
#
gphys = GPhys::NetCDF_IO.open(file[0], var)

#--------------------------------------
# ECMWF データから生成した gphys オブジェクトの開始日時を日付オブジェクトに格納するメソッド.
def get_start_date(gphys) 
  start_org = gphys.coord(gphys.axnames.index("date")).get_att('units')
  start = start_org[/(\d)+-(\d)+-(\d)+/]            # "2000-01-12" とか "99-1-3" で指定される部分を返す
  startd = Date.parse(start)
  return startd
end
# メソッドテスト
# p sd = get_start_date(gphys).to_s
#-------------------------------------

#------ function::get_mdays -------
# 引数に与えた日付オブジェクトの年の各月の日数を返す
#----------------------------------
def get_mdays(d)
  mdays = Array.new
  1.upto(12){|i|
    if i == 12 then 
      date = Date.new(d.year, 12, 31) 
      mdays << date.yday - (date << 1).yday
    else 
      fdate = Date.new(d.year, i) #forward_month_date
      adate = Date.new(d.year, i+1) #after_month_date
      mdays << adate.yday - fdate.yday
    end
  }
  return mdays
end
#メソッドテスト
# test = Date.new(2004)
#p get_mdays(test)
#--------------------------------------


#------ function::get_sum_mdays -------
# 引数に与えた日付オブジェクトの年の各月の累積日数
#--------------------------------------
def get_sum_mdays(d)
  ydays = Array.new
  1.upto(12){|i|
    date = Date.new(d.year, i)
    ydays << date.yday
  }
  return ydays
end
#メソッドテスト
#  test = Date.new(2004)
# p ydays = get_sum_mdays(test)
#----------------------------------

#------ function::get_drange -------
# 引数に与えた日付オブジェクトの年の各月の累積日数
#--------------------------------------
def get_drange(gphys, sd)
  mdays = get_mdays(sd) # [31, 28, 31, 30, ..., 31]
  sumd = get_sum_mdays(sd) # 1, 32, 
  na_sumd = NArray.to_na(sumd) + NArray.to_na(mdays) -2
end

#------------- main ------------------

# オプション解析 ----------------------------------------------------
unless getopts("hH", "month:")
  print "#{$0}:illegal options. please exec this program with -h/--help. \n"
  exit 1
end


p a = Time.now # 処理にかかった時間を計測するための布石.

p sd = get_start_date(gphys)
p ad = get_drange(gphys, sd)
ydays = NArray.sint(365).indgen!

# 設定部分
## RUBYLIB に追加
$: << '/GFD_Dennou_Club/dc-arch/dcchart/daktu32/ERA40'
$: << '/GFD_Dennou_Club/dc-arch/dcchart/daktu32/ERA40/bin'
$: << '/work31/daktu32/ERA40'
$: << '/work31/daktu32/ERA40/bin'
## 追加メソッド定義ファイル
require "libgphys-e.rb"
## 属性設定ファイル
require "ERA40-ncattr.conf.rb"

var_attr= {"units"=>"days since 2001-1-1"}

num_month = 5
num_month = ($OPT_month).to_i if ($OPT_month)

if num_month == 12 then
  sd.month.upto(sd.month + num_month-1){|i|
    if i == 1 then
      output = "daily_#{pname}_2001-0#{i}_ERA40.nc"
      p "now making "+output
      gphys.cut('date'=>0..ad[0]).save(output, $global_attr, nil, nil, $history)
      p "done"
    elsif (i < 10) then
      output = "daily_#{pname}_2001-0#{i}_ERA40.nc"
      p "now making "+output
      gphys.cut('date'=>ad[i-2]+1..ad[i-1]).save(output, $global_attr, nil, nil, $history)
      p "done"

    else 
      output = "daily_#{pname}_2001-#{i}_ERA40.nc"
      p "now making "+output
      gphys.cut('date'=>ad[i-2]+1..ad[i-1]).save(output, $global_attr, nil, nil, $history)
      p "done"
    end
  }    
else
  sd.month.upto(sd.month + num_month-1){|i|
    if i == 1 then
      output = "daily_#{pname}_2001-0#{i}_ERA40.nc"
      p "now making "+output
      gphys.cut('date'=>0..ad[0]).save(output, $global_attr, nil, nil, $history)
      p "done"
    elsif i < 7 then
      output = "daily_#{pname}_2001-0#{i}_ERA40.nc"
      p "now making "+output
      gphys.cut('date'=>ad[i-2]+1..ad[i-1]).save(output, $global_attr, nil, nil, $history)
      p "done"
    elsif (i >= 7) & (i < 10) then
      output = "daily_#{pname}_2001-0#{i}_ERA40.nc"
      p "now making "+output
      index = gphys.axnames.index("date")
      newgphys = gphys.cut('date'=>ad[i-2]-ad[5]..ad[i-1]-ad[5]-1)
      newgphys.save(output, $global_attr, "date", var_attr, $history)
      
      ng = NetCDF.open(output, "a")
      ng.var("date").put(ydays[ad[i-2]+1..ad[i-1]])
      ng.close
      
      p "done"
    else 
      output = "daily_#{pname}_2001-#{i}_ERA40.nc"
      p "now making "+output
      index = gphys.axnames.index("date")
      newgphys = gphys.cut('date'=>ad[i-2]-ad[5]..ad[i-1]-ad[5]-1)
      newgphys.save(output, $global_attr, "date", var_attr, $history)
      
      ng = NetCDF.open(output, "a")
      ng.var("date").put(ydays[ad[i-2]+1..ad[i-1]])
      ng.close
      
    end
  }
end

p Time.now - a