#!/usr/bin/env ruby
=begin

= gtstspct.rb : 時空間(パワー)スペクトルを求めるスクリプト.
  
= USAGE :

  $ gtstrm.rb -[hHo] <--output=[outputname]> [ncfile]

= note

* 結果として得られる変数は "pws".

* パワースペクトル導出は以下の式を導入. (TeX の書式で書いてます.)

=end
## 保存する変数名
#$new_vname = "tp"
## 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"

#############  以下メインルーチン  ##############

require "getopts"
require "numru/netcdf"
module NumRu
  class NetCDFVar
#    alias put scaled_put
    alias get scaled_get
  end
end

require "numru/ggraph"
include NumRu

NetCDFVarDeferred.max_len_store = 250000000

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


if ($OPT_h) || ($OPT_H) || ($OPT_help) || ARGV.size == 0 then
  print <<EOF

EOF
exit 0
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


if file.size == 1 then  # 単体の nc ファイルの場合
  gphys = GPhys::NetCDF_IO.open(file[0].to_s, var)
elsif ($OPT_c) 
  gphys = GPhys::NetCDF_IO.open(file, var)
else         
  vgarray = Array.new
  file.each do |f|
    vgarray << GPhys::NetCDF_IO.open(f.to_s, var)
  end
  vgphys = vgarray[0]
  1.upto(vgarray.length-1) {|n|
    vgphys += vgarray[n]
  }
  gphys = vgphys/file.length
end



# range
if ($OPT_range)
  range = ($OPT_range).to_s.split(/\s*,\s*/)
  range.each{|rg|
    ran = rg.to_s.split(/\s*=\s*/)
    if ran[1] =~ /\s*\.\.\s*/
      rg1=ran[1].split(/\s*\.\.\s*/)[0].to_f
      rg2=ran[1].split(/\s*\.\.\s*/)[1].to_f
      gphys = gphys.cut(ran[0] => rg1..rg2)
    else
      gphys = gphys.cut(ran[0] => ran[1].to_f)
    end
  }
end


# title
$title = ($OPT_title) if ($OPT_title)

# mean
if ($OPT_mean)
  mean = ($OPT_mean).to_s.split(/\s*,\s*/)
  mean.each{|mn|
    gphys = gphys.mean(gphys.axnames.index(mn))
  }
end


# 保存ファイル名が与えられたら, その名前で保存. 
if ($OPT_o) then
  output = ($OPT_o).to_s
elsif ($OPT_output) then
  output = ($OPT_output).to_s
else 
  output = "gtstrm.nc"
end


p stspct = gphys.stspct_coh(var)
stspct.save(output, $global_attr, $new_vname, $var_attr, $history)  # nc ファイルに保存