#!/usr/bin/env ruby

=begin

$  for file in *.xwd ; do xwdtopnm ${file} |  pnmcut 2 2 904 654 > test.pnm; ppmtogif test.pnm > ${file}.gif; done
$ gifsicle --delay=50 --no-loop dcl_*.gif > movie.gif

=end

# ----------------------------------------------
# プログラムメイン

END{

  num = 1
  set_filename(num)           # ファイル名等基本設定
  mradl_anim_main             # gif 作成
  print "gif convert start #{Time.now}\n"
  mkgifanim                   # gif アニメ作成シェルスクリプト

}

# ----------------------------------------------
# 必要なライブラリ, モジュールの読み込み

$local_path = '/home/yukiko/work/ape/yukiko/lib'
$: << $local_path

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

require "numru/ggraph"
require "colorbar"
include NumRu
include NMath


#----------------------------------------------
# 設定

def set_filename(num)

  numnum = 0
  num = ARGV[numnum+1].to_i  if numnum = ARGV.index("-num")

  $latitude = - 1.5
  $vector = true 
  $exp = "T39L48_eml"
  $timestep = "2hour"

  $timestep = "1hour"           if numnum = ARGV.index("-1hour")
  $latitude = ARGV[numnum+1].to_f if numnum = ARGV.index("-lat")
  $exp = ARGV[numnum+1]           if numnum = ARGV.index("-exp")

  $lattitude = -0.4  if $exp =~ /T159/

  $lat_label = 0  if $latitude == - 1.5
  $lat_label = 0  if $latitude == - 0.4

  $fig_path = "/home/yukiko/work/ape/GrADS/fig/"
  $dir_path = "/home/yukiko/work/ape/GrADS/#{$exp}_expID01/hourly/"

  ## へんすうじょうほう
  set_varinfo

  if num == 1
    varlist = "gph"
    $var_array.push($gph_hash) 
  elsif num == 2
    varlist = "t"
    $var_array.push($t_hash) 
  elsif num == 3
    varlist = "q"
    $var_array.push($q_hash) 
  elsif num == 4
    varlist = "u"
    $var_array.push($u_hash) 
  end


  if $timestep == "1hour"
    $gifname = 
      "#{$exp}_control_#{varlist}_y#{$lat_label}z0512-#{$timestep}"
  else
    $gifname = 
      "#{$exp}_control_#{varlist}_y#{$lat_label}z0512"
  end

  print "#{$gifname} name start #{Time.now}\n"

end

def set_varinfo

  rmiss = DCL.glpget('rmiss')
  $var_array = []

  # 凝結加熱 (cumulus + largescale)
  $qcnd_hash = {
    "variable" => "DTCON", 
    "diff"     => false, 
    "levels" => NArray[rmiss,-1.6e-4,-0.8e-4, 0, 
      0.8e-4, 1.6e-4, 2.4e-4, rmiss],
#    "levels" => NArray[0, 0.4e-4, 0.8e-4, 1.2e-4, 
#      1.6e-4, 2.0e-4, 2.4e-4, rmiss],
    "lost_axis" => " ",
    "title" => "condensation heating (K s-1)"
  }

  # 凝結加熱 (cumulus)
  $qcndc_hash = {
    "variable" => "DTCUM", 
    "diff"     => false, 
    "levels" => NArray[rmiss,-1.6e-4,-0.8e-4, 0, 
      0.8e-4, 1.6e-4, 2.4e-4, rmiss],
#    "levels" => NArray[0, 0.4e-4, 0.8e-4, 1.2e-4, 
#      1.6e-4, 2.0e-4, 2.4e-4, rmiss],
    "lost_axis" => " ",
    "title" => "cumulus condensation heating (K s-1)"
  }

  # 凝結加熱 (large scale)
  $qcndl_hash = {
    "variable" => "DTLSC", 
    "diff"     => false, 
    "levels" => NArray[rmiss,-1.6e-4,-0.8e-4, 0, 
      0.8e-4, 1.6e-4, 2.4e-4, rmiss],
#    "levels" => NArray[0, 0.4e-4, 0.8e-4, 1.2e-4, 
#      1.6e-4, 2.0e-4, 2.4e-4, rmiss],
    "lost_axis" => " ",
    "title" => "L.S. condensation heating (K s-1)"
  }

  # 温度偏差
  $t_hash = {
    "variable" => "T", 
    "diff"     => true, 
    "levels" => NArray[rmiss,-1,-0.5, 0, 0.5, 1, 1.5, rmiss], 
    "lost_axis" => "(diff) from (mean) zonal",
    "title" => "temperature (K)"
  }
  
  # 比湿偏差
  $q_hash = {
    "variable" => "Q", 
    "diff"     => true, 
    "levels" => NArray[rmiss,-2e-3,-1e-3, 0, 1e-3, 2e-3, 3e-3,rmiss], 
    "lost_axis" => "(diff) from (mean) zonal", 
    "title" => "specific humidity (1)"
  }

  # 鉛直流
  $sigdot_hash = {
    "variable" => "OMG", 
    "diff"     => false, 
    "levels" => NArray[rmiss,-30e-2,-15e-2, 0, 15e-2, 30e-2, 45e-2, rmiss], 
    "lost_axis" => " ",
    "title" => "p-vertical velocity (Pa s-1)"
  }

  # ジオポテンシャル
  $gph_hash = {
    "variable" => "Z", 
    "diff"     => true, 
    "levels" => NArray[rmiss,-10,-5, 0, 5, 10, 15, rmiss], 
    "lost_axis" => "(diff) from (mean) zonal",
    "title" => "geopotential height (m)"
  }

  # 東西風
  $u_hash = {
    "variable" => "U", 
    "diff"     => true, 
    "levels" => NArray[rmiss,-5,-2.5, 0, 2.5, 5, 7.5, rmiss], 
    "lost_axis" => "(diff) from (mean) zonal",
    "title" => "eastward wind (m s-1)"
  }

  # 南北風
  $v_hash = {
    "variable" => "V", 
    "diff"     => true, 
    "levels" => NArray[rmiss,-4,-2, 0, 2, 4, 6, rmiss], 
    "lost_axis" => "(diff) from (mean) zonal",
    "title" => "northward wind (m s-1)"
  }

  # 長波放射加熱
  $qradl_hash = {
    "variable" => "DTRADL", 
    "diff"     => false, 
    "levels" => NArray[rmiss,-1.5e-5,-1.2e-5,-0.9e-5,-0.6e-5,-0.3e-5, 0, rmiss], 
    "lost_axis" => " ",
    "title" => "longwave heating (K s-1)"
  }

  # 短波放射加熱
  $qrads_hash = {
    "variable" => "DTRADS", 
    "diff"     => false, 
    "levels" => NArray[0, 5e-6, 10e-6, 15e-6, 20e-6, 25e-6, 30e-6, rmiss], 
    "lost_axis" => " ",
    "title" => "shortwave heating (K s-1)"
  }

  # 鉛直拡散加熱
  $qvdf_hash = {
    "variable" => "DTVDF", 
    "diff"     => false, 
    "levels" => NArray[rmiss, 10e-6, 20e-6, 30e-6, 40e-6, 50e-6, 60e-6, rmiss], 
    "lost_axis" => " ",
    "title" => "diffusive heating (K s-1)"
  }


  # -----------------------
  $day_hash = {
    "T39L48_eml"=> 40,
    "T39L48_ias"=> 45,
    "T39L48_ksc"=> 40,
    "T39L48_kuo"=> 15,
    "T39L48_mca"=> 40,
    "T39L48_non"=> 20, 
    "T159L48_non"=> 60, 
    "T159L48_eml"=> 25
  }

end


# ----------------------------------------------
# お絵描きメイン

def mradl_anim_main

  # ------------------------------------------------------------
  # 温度

  var = $var_array[0]["variable"]
  gphys_1 = GPhys::IO.open("#{$dir_path}/0001/#{var}.ctl", var )
  gphys_11 = GPhys::IO.open("#{$dir_path}/0002/#{var}.ctl",var )
  gphys_111 = GPhys::IO.open("#{$dir_path}/0003/#{var}.ctl",var )

  gphys_u1 = GPhys::GrADS_IO.open("#{$dir_path}/0001/U.ctl", "U")
  gphys_u2 = GPhys::GrADS_IO.open("#{$dir_path}/0002/U.ctl", "U")
  gphys_u3 = GPhys::GrADS_IO.open("#{$dir_path}/0003/U.ctl", "U")

  gphys_v1 = GPhys::GrADS_IO.open("#{$dir_path}/0001/V.ctl", "V")
  gphys_v2 = GPhys::GrADS_IO.open("#{$dir_path}/0002/V.ctl", "V")
  gphys_v3 = GPhys::GrADS_IO.open("#{$dir_path}/0003/V.ctl", "V")
  
  gphys_sig1 = GPhys::GrADS_IO.open("#{$dir_path}/0001/OMG.ctl", "OMG")
  gphys_sig2 = GPhys::GrADS_IO.open("#{$dir_path}/0002/OMG.ctl", "OMG")
  gphys_sig3 = GPhys::GrADS_IO.open("#{$dir_path}/0003/OMG.ctl", "OMG")

 
  # ------------------------------------------------------------
  DCL.swpset('LDUMP', true) 
  DCL.swpset('LWAIT',false) ; DCL.swpset('LWAIT0',false) 
  DCL.swpset('LWAIT1',false)
  DCL.swpset('IPOSX', 50) ; DCL.swpset('IPOSY',50)  

  DCL.gropn(1)

  DCL.uzfact(0.6)                # 文字の大きさ倍
  DCL.sgpset('lcorner',false)    # コーナーマークはつけない
  DCL.slmgn(0.0, 0.0, 0.0, 0.0)  # 描画マージン指定
  DCL.sgpset('lfprop',true)      # プロポーショナルフォントを使う
  DCL.sgpset('lfull',true)       # 全画面表示
  DCL.sgpset('lcntl', false)     # アンダーバーは普通に表示
  DCL.uscset('cyspos', 'B' )      # y 軸単位を書く位置

#  ipara = NumRu::DCL.uzpget('RSIZEC1')
  DCL.uzpset('RSIZEC0',0.012)    # 小タイトル 文字大きさ
  DCL.uzpset('INDEXL0',1)       # 小タイトル 文字ラインインデックス
  rmiss = DCL.glpget('rmiss')
  colorbar_hash = {
    'vlength'=>0.13,
    "vwidth"=>0.01,
    'labelintv'=>1, 
    'voff'=>0.02,
    "constwidth" => true
  }
  patterns = NArray[30999, 35999, 45999, 55999, 70999, 75999, 85999]

  if $timestep == "1hour"
    timeall = 15*24
    timestep = 1
  elsif $timestep == "2hour"
    timeall = 15*12
    timestep = 2
  end

#  (15*12).times{ |num|
#     (1).times{ |num|
  timeall.times{ |num|

    $file_label = "day=#{format("%.1f",(num*1.0*timestep/24+ $day_hash[$exp]) )}" 
    print "  #{$file_label} #{Time.now}\n"

# 1 枚目: 凝結加熱
    DCL.uzpset('labelxb',false)
    GGraph.set_fig('viewport'=>[0.10,0.80,0.52,0.65],'new_frame'=>true) 
    GGraph.next_axes('xunits'=>'','yunits'=>'1','xtitle'=>'','ytitle'=>"sigma") 

    tone_hash = { 'patterns'=> patterns, 'levels'=>$var_array[0]["levels"] }

    if (120 > num*timestep)
      g = gphys_1.cut(true,$latitude,true,true)[true,true,num*timestep]
      g_tmp = gphys_1.cut(true,$latitude,true,true).mean(0)[true,num*timestep].
        set_lost_axes("")
    elsif (120 <= num*timestep) && (240 > num*timestep)
      g = gphys_11.cut(true,$latitude,true,true)[true,true,num*timestep-120]
      g_tmp = gphys_11.cut(true,$latitude,true,true).
        mean(0)[true,num*timestep-120].set_lost_axes("")
    else
      g = gphys_111.cut(true,$latitude,true,true)[true,true,num*timestep-240]
      g_tmp = gphys_111.cut(true,$latitude,true,true).
        mean(0)[true,num*timestep-240].set_lost_axes("")
    end

    g = g.set_att("long_name","").set_att("unit","").set_lost_axes("").
      lon_lotate

    if $var_array[0]["diff"] 
      GGraph.tone( g-g_tmp, true , tone_hash ) 
#      GGraph.tone( g-g_tmp ) 
    else 
      g = - g if $var_array[0]["variable"] == "OMG"
      GGraph.tone( g, true , tone_hash ) 
#      GGraph.tone( g ) 
    end
    DCL.uulin([-180,180],[0.25,0.25])
    DCL.uulin([-180,180],[0.85,0.85]) unless $surface

    if (120 > num*timestep)
      vect_vertical_plot(gphys_u1,gphys_sig1,num,timestep,0) if $vector 
    elsif (120 <= num*timestep) && (240 > num*timestep)
      vect_vertical_plot(gphys_u2,gphys_sig2,num,timestep,-120) if $vector 
    else
      vect_vertical_plot(gphys_u3,gphys_sig3,num,timestep,-240) if $vector 
    end

    DCL::uzinit 
    DCL::uxsttl("t", $var_array[0]["title"], -1.0 ) 
    DCL::uzinit 
    DCL.uxpttl("t", 0, "lat=#{$latitude}, #{$var_array[0]["lost_axis"]}", 1.0)
    DCL::uxsttl("t", "#{$file_label}" , 1.0 )
    GGraph.color_bar( colorbar_hash )
    DCL::uzinit 

# 2 枚目: 上層水平面 (z=12)

    level = 0.25

    DCL.uzpset('labelxb',false)
    GGraph.next_axes('xunits'=>'','xtitle'=>nil, 'yunits'=>"deg") 
    GGraph.set_fig('viewport'=>[0.10,0.80,0.34,0.47], 'new_frame'=>false )
    tone_hash = { 'patterns'=> patterns, 'levels'=>$var_array[0]["levels"] }

    if (120 > num*timestep)
      g = gphys_1.cut(true,-30..30,level,true)[true,true,num*timestep]
      g_tmp = gphys_1.cut(true,-30..30,level,true).mean(0)[true,num*timestep].
        set_lost_axes("")
    elsif (120 <= num*timestep) && (240 > num*timestep)
      g = gphys_11.cut(true,-30..30,level,true)[true,true,num*timestep-120]
      g_tmp = gphys_11.cut(true,-30..30,level,true).
        mean(0)[true,num*timestep-120].set_lost_axes("")
    else
      g = gphys_111.cut(true,-30..30,level,true)[true,true,num*timestep-240]
      g_tmp = gphys_111.cut(true,-30..30,level,true).
        mean(0)[true,num*timestep-240].set_lost_axes("")
    end

    g = g.set_att("long_name","").set_att("unit","").set_lost_axes("").
      lon_lotate

    if $var_array[0]["diff"] 
      GGraph.tone( g-g_tmp, true , tone_hash ) 
#      GGraph.tone( g-g_tmp ) 
    else 
      g = - g if $var_array[0]["variable"] =~ /OMG/
      GGraph.tone( g, true , tone_hash ) 
#      GGraph.tone( g ) 
    end

    if (120 > num*timestep)
      vect_horizontal_plot(gphys_u1,gphys_v1,num,timestep,level,0) if $vector 
    elsif (120 <= num*timestep) && (240 > num*timestep)
      vect_horizontal_plot(gphys_u2,gphys_v2,num,timestep,level,-120) if $vector 
    else
      vect_horizontal_plot(gphys_u3,gphys_v3,num,timestep,level,-240) if $vector 
    end

#    DCL::uxsttl("t", $var_array[0]["title"], -1.0 ) 
    DCL::uzinit 
    DCL.uxpttl("t", 0, "sigma=#{level}, #{$var_array[0]["lost_axis"]}", 1.0)
    GGraph.color_bar( colorbar_hash )
    DCL::uzinit 

# 3 枚目: 下層水平面 (z=1,5)

    if $surface
      level = 0.99499
    else
      level = 0.85
    end

    DCL.uzpset('labelxb',true)
    GGraph.set_fig('viewport'=>[0.10,0.80,0.16,0.29], 'new_frame'=>false )
    GGraph.next_axes('yunits'=>"deg",'xunits'=>"deg") 

    tone_hash = { 'patterns'=> patterns, 'levels'=>$var_array[0]["levels"] }

    if (120 > num*timestep)
      g = gphys_1.cut(true,-30..30,level,true)[true,true,num*timestep]
      g_tmp = gphys_1.cut(true,-30..30,level,true).mean(0)[true,num*timestep].
        set_lost_axes("")
    elsif (120 <= num*timestep) && (240 > num*timestep)
      g = gphys_11.cut(true,-30..30,level,true)[true,true,num*timestep-120]
      g_tmp = gphys_11.cut(true,-30..30,level,true).
        mean(0)[true,num*timestep-120].set_lost_axes("")
    else
      g = gphys_111.cut(true,-30..30,level,true)[true,true,num*timestep-240]
      g_tmp = gphys_111.cut(true,-30..30,level,true).
        mean(0)[true,num*timestep-240].set_lost_axes("")
    end

    g = g.set_att("long_name","").set_att("unit","").set_lost_axes("").
      lon_lotate

    if $var_array[0]["diff"] 
      GGraph.tone( g-g_tmp, true , tone_hash ) 
#      GGraph.tone( g-g_tmp ) 
    else 
      g = - g if $var_array[2]["variable"] =~ /sigdot/
      GGraph.tone( g, true , tone_hash ) 
#      GGraph.tone( g ) 
    end

    if (120 > num*timestep)
      vect_horizontal_plot(gphys_u1,gphys_v1,num,timestep,level,0) if $vector 
    elsif (120 <= num*timestep) && (240 > num*timestep)
      vect_horizontal_plot(gphys_u2,gphys_v2,num,timestep,level,-120) if $vector 
    else
      vect_horizontal_plot(gphys_u3,gphys_v3,num,timestep,level,-240) if $vector 
    end

    DCL::uzinit 
#    DCL::uxsttl("t", $var_array[0]["title"], -1.0 ) 
#    DCL::uzinit 
      DCL.uxpttl("t", 0, "sigma=#{level}, #{$var_array[0]["lost_axis"]}", 1.0) 
    GGraph.color_bar( colorbar_hash )

    
  }
  
  DCL::grcls 

end


# gif アニメ作成シェルスクリプト
def mkgifanim

  `for file in *.xwd ; do xwdtopnm ${file} |  pnmcut 2 2 904 654 > tmp.pnm; ppmtogif tmp.pnm > ${file}.gif; done >& /dev/null `
  `gifsicle --loopcount=100 dcl_*.gif > #{$fig_path}#{$gifname}_anim.gif`
  `mv dcl_001.xwd.gif #{$fig_path}#{$gifname}_t0.gif`
  `rm dcl_* tmp.pnm`

end



# ----------------------------------------------
# GPhys メソッド追加定義

def GGraph.annotate(str_ary) # GGraph の annotate メソッド再定義

  charsize = 0.7 * DCL.uzpget('rsizec1')
  dvx = 0.01
  dvy = charsize*1.5
  raise TypeError,"Array expected" if ! str_ary.is_a?(Array)
  vxmin,vxmax,vymin,vymax = DCL.sgqvpt
  vx = 0.66 
#  vy = 0.12 - charsize/2
#  vy = 0.10 - charsize/2
  vy = 0.11 - charsize/2
  str_ary.each{|str|
    DCL::sgtxzv(vx,vy,"[ #{str} ]",charsize,0,-1,1)
    vy -= dvy
  }
  nil
end

class GPhys

  # 座標, データ操作ログ記録関数の定義 (Grid.lost_axes 参照)
  def add_lost_axes( lost )
    GPhys.new( @grid.copy.add_lost_axes( lost.to_a ), @data)
  end
  def set_lost_axes( lost )
    GPhys.new( @grid.copy.set_lost_axes( lost.to_a ), @data)
  end

  # rename 再定義 (gphys class の値を返す)
  def rename(name)
    GPhys.new(@grid,@data.copy.rename(name))
  end

  # attribute 変更再定義 (gphys class の値を返す)
  def set_att(name,val)
    GPhys.new(@grid,@data.copy.set_att(name,val))
  end

  # 経度 180 度回転
  def lon_lotate

    # 軸の回転 (180 度引き算)
#    lon = VArray.new(@grid.axis("lon").pos.val-180).rename("lon")
    lon = VArray.new(@grid.axis(0).pos.val-180).rename("lon").
      put_att("units",@grid.axis(0).pos.get_att("units"))
    lon = Axis.new().set_pos(lon)
    @grid = @grid.change_axis(0, lon)

#    lat = VArray.new(@grid.axis(1).pos.val).rename("lat")
#    lat = Axis.new().set_pos(lat)
#    @grid = @grid.change_axis(1, lat)

    # データの回転
    dim = @data.val.shape[0]
    data_val =  @data.copy.val
    data_val[0..dim/2-1, false] = @data.val[dim/2..dim-1, false]
    data_val[dim/2..dim-1, false] = @data.val[0..dim/2-1, false]
    data = @data.copy
    data.val = (data_val)

    # gphys class の値を返す
    GPhys.new( @grid, data )

  end

end


# ----------------------------------------------
def vect_vertical_plot(gphys_u,gphys_sig,num,timestep,numdiff)

    g_u =  gphys_u.cut(true,$latitude,true,true)[true,true,num*timestep+numdiff].
      set_lost_axes($lost_axis)
    g_sig = gphys_sig.cut(true,$latitude,true,true)[true,true,num*timestep+numdiff].
      set_lost_axes($lost_axis)


      g_u = g_u.lon_lotate
      g_sig = g_sig.lon_lotate


    # ベクトルの間引
    u_data  = NArrayMiss.sfloat( g_u.axis(0).pos.val.size , 
			   g_u.axis(1).pos.val.size ).fill!(0.0)
 
    sig_data = NArrayMiss.sfloat( g_sig.axis(0).pos.val.size , 
			    g_sig.axis(1).pos.val.size ).fill!(0.0)

    num_l = 3
    num_l = 3 * 2*2 if $exp =~ /T159/
    num_z = 2

    g_u.axis(0).pos.val.size.times{ |num|
      if (num % num_l) == 0
	u_data[num,true] =  g_u.data.val[num,true]
	sig_data[num,true] =  g_sig.data.val[num,true]
      elsif
	u_data[num,true]  = 0.0
	sig_data[num,true]  =  0.0
      end
    }
    
    g_u.axis(1).pos.val.size.times{ |num|
      if (num % num_z) == 0
	u_data[true,num] = u_data[true,num]
	sig_data[true,num] = sig_data[true,num]
      elsif
	u_data[true,num] = 0.0
	sig_data[true,num] = 0.0
      end
    }

    #  (u,-sigdot) 
    DCL::uglset('LMSG', true)
    DCL::uglset('LNRMAL', false)

#    DCL::ugrset('XFACT1', 10E-6*100/2)
#    DCL::ugrset('YFACT1', 15.0*1000/2)
    DCL::ugrset('XFACT1', 10E-6*100)
    DCL::ugrset('YFACT1', 10.0/100)
    DCL::ugrset('UXUNIT', -999)
    DCL::ugrset('UYUNIT', -999)
    DCL::ugrset('VXULOC', -999)
    DCL::ugrset('VYULOC', -999)

    DCL::uglset('LUNIT', true)
    DCL::ugrset('VXUNIT', 0.027)
    DCL::ugrset('VYUNIT', 0.027)
    DCL::ugrset('VXUOFF', 0.044)
    DCL::ugrset('VYUOFF', -0.025)

  DCL::ugrset('RSIZET', 0.0115)
  DCL.uzpset('INDEXL2',1)        # 文字ラインインデックス
  DCL::ugpset('INDEX', 2)        # 矢印インデックス
  
  pad = DCL.uzpget('PAD1')
  DCL.uzpset('PAD1',-0.1)

  DCL::ugvect( u_data, -sig_data ) 
  DCL.uzpset('INDEXL2',3)        # 文字ラインインデックス
  DCL.uzpset('PAD1',pad)

end

def vect_horizontal_plot(gphys_u,gphys_v,num,timestep,level,numdiff)
  
  g_u =  gphys_u.cut(true,-30..30,level,true)[true,true,num*timestep+numdiff].
    set_lost_axes($lost_axis)
  g_v = gphys_v.cut(true,-30..30,level,true)[true,true,num*timestep+numdiff].
    set_lost_axes($lost_axis)

  g_u = g_u.lon_lotate
  g_v = g_v.lon_lotate


    # ベクトルの間引
    u_data  = NArrayMiss.sfloat( g_u.axis(0).pos.val.size , 
			   g_u.axis(1).pos.val.size ).fill!(0.0)
 
    v_data = NArrayMiss.sfloat( g_v.axis(0).pos.val.size , 
			    g_v.axis(1).pos.val.size ).fill!(0.0)

    num_l = 3
    num_m = 2

    num_l = 3 * 2*2 if $exp =~ /T159/
    num_m = 2 * 2*2 if $exp =~ /T159/

    g_u.axis(0).pos.val.size.times{ |num|
      if (num % num_l) == 0
	u_data[num,true] =  g_u.data.val[num,true]
	v_data[num,true] =  g_v.data.val[num,true]
      elsif
	u_data[num,true]  = 0.0
	v_data[num,true]  =  0.0
      end
    }
    
    g_u.axis(1).pos.val.size.times{ |num|
      if (num % num_m) == 0
	u_data[true,num] = u_data[true,num]
	v_data[true,num] = v_data[true,num]
      elsif
	u_data[true,num] = 0.0
	v_data[true,num] = 0.0
      end
    }

    #  (u,v) 
    DCL::uglset('LMSG', true)
    DCL::uglset('LNRMAL', false)

    DCL::ugrset('XFACT1', 0.0008)
    DCL::ugrset('YFACT1', 0.0008)
    DCL::ugrset('UXUNIT', -999)
    DCL::ugrset('UYUNIT', -999)
    DCL::ugrset('VXULOC', -999)
    DCL::ugrset('VYULOC', -999)

    DCL::uglset('LUNIT', true)
    DCL::ugrset('VXUNIT', 0.027)
    DCL::ugrset('VYUNIT', 0.027)
    DCL::ugrset('VXUOFF', 0.044)
    DCL::ugrset('VYUOFF', -0.025)

    DCL::ugrset('RSIZET', 0.0115)
    DCL.uzpset('INDEXL2',1)        # 文字ラインインデックス
    DCL::ugpset('INDEX', 2)        # 矢印インデックス

    pad = DCL.uzpget('PAD1')
    DCL.uzpset('PAD1',-0.1)

    DCL::ugvect( u_data, v_data ) 
    DCL.uzpset('INDEXL2',3)        # 文字ラインインデックス
    DCL.uzpset('PAD1',pad)

end