# -*- coding: utf-8 -*-
# 運動エネルギー密度を計算し, 絵を描く

require "numru/ggraph"
include NumRu
include GGraph

velx_file = "./BS1998_VelX.nc"
velx_var ="VelX"
vely_file = "./BS1998_VelY.nc"
vely_var ="VelY"
velz_file = "./BS1998_VelZ.nc"
velz_var ="VelZ"
dens_file = "./BS1998_restart.nc"
dens_var = "DensBZ"

=begin
DCL.swpset('iwidth',700)      # window width                                    
DCL.swpset('iheight',700)     # window height                                   
###DCL.swpset('ldump',true)      # dump image files                             
DCL.swpset('lwait',false)     # do not wait mouse click to show the next page   
###DCL.swpset('lalt',true)       # background plot                              
DCL.sgscmn(10)                # change colomap (see below)                      
###DCL.sgscmn(5)                 # change colomap (see below)   
=end

#DCL.swlset( 'ldump', true )
DCL.gropn(2)

velx = GPhys::IO.open("#{velx_file}", "#{velx_var}")
vely = GPhys::IO.open("#{vely_file}", "#{vely_var}")
velz = GPhys::IO.open("#{velz_file}", "#{velz_var}")
dens = GPhys::IO.open("#{dens_file}", "#{dens_var}")
time = GPhys::IO.open("#{velx_file}", 't').val
x_dco = GPhys::IO.open("#{dens_file}", 'x').val
x_co = GPhys::IO.open("#{velx_file}", 'x').val
z_co = GPhys::IO.open("#{velz_file}", 'z').val
#kener_sum = Array.new(2)
#kener_tmp = Array.new(2)
#=begin

p dens

denscut = dens[5..104, 0..-1, 5..104]

puts denscut

kener = 0.50 * denscut * ( velx**2.0 + velz**2.0)
kener = kener.mean('x').mean('z').mean('y')



denscut = denscut.mean('x').mean('z').mean('y')
p kener
puts "denscut"
p denscut
kenerdens = kener 
#kenerdens = kener / denscut


#=begin
GGraph.set_axes('xtitle'=>'time','ytitle'=>'K')
GGraph.line( kenerdens, true, 'title'=>'Kinetic energy density')

#GGraph.line( kenerdens.cut( 'y'=>100 ), true)
#=end
DCL.grcls