# -*- coding: euc-jp -*-
require "numru/ggraph"
include NumRu
include NMath

hitran_molname = [
  '', 'H_2"O', 'CO_2"', 'O_3"', 'N_2"O', 'CO', 'CH_4"', 'O_2"', 'NO', 'SO_2"',
  'N_2"O', 'NH_3"', '', '', 'HF', '', '', '', '', 'OCS',
  '', '', 'N_2"', '', '', '', '', '', '', '', '',
]
cldname = ["dust"]
molnum = [1,2,3]
ncfn  = 'out/Mars_Seiff1982+DustMMR_H2O-0.0_CO2-0.95_O3-0.0_DOD0.2.nc'

gpt   = GPhys::IO.open( ncfn, 'Temp' )
gpgas = GPhys::IO.open( ncfn, 'VMR' )
gpcld = GPhys::IO.open( ncfn, 'PtclMMR' )
gpmolnum  = GPhys::IO.open( ncfn, 'MolNum' )
gpcldnum  = GPhys::IO.open( ncfn, 'PtclNum' )
gpcldname = GPhys::IO.open( ncfn, 'PtclName' )

gpcld = gpcld.copy
gpcld.long_name = 'Mass mixing ratio'

#gpt = gpt.copy
#z = gpt.axis('Press').pos.convert_units( Units['atm'] )
#z.units = 'atm'
#gpt.axis('Press').set_pos(z)
#
#gpgas = gpgas.copy
#z = gpgas.axis('Press').pos.convert_units( Units['atm'] )
#z.units = 'atm'
#gpgas.axis('Press').set_pos(z)
#
#gpcld = gpcld.copy
#z = gpcld.axis('PressM').pos.convert_units( Units['atm'] )
#z.units = 'atm'
#gpcld.axis('PressM').set_pos(z)

gpp = gpt.coord('Press')


puts "# pressure (hPa), temperature (K), QH2O (1), QO3 (1)"
for k in 0..gpt.val.size-1
  puts gpp.val[k].to_s + ', ' + gpt.val[k].to_s #+ ', ' + gpq.val[k].to_s + ', ' + gpo3.val[k].to_s
end

iws = (ARGV[0] || (puts ' WORKSTATION ID (I)  ? ;'; DCL::sgpwsn; gets)).to_i
DCL.gropn(iws)

DCL.sldiv('y',3,1)           # 2x2に画面分割, 'y'=yoko: 左上→右上→左下...
DCL.sgpset('lcntl', true)   # 制御文字を解釈しない
DCL.sgpset('lfull', true)     # 全画面表示
DCL.uzfact(0.7)             # 座標軸の文字列サイズを 0.75 倍
DCL.sgpset('lfprop',true)    # プロポーショナルフォントを使う

DCL.sgpset('lclip', true)

#DCL.glpset('lmiss',true)
#DCL.glpset('rmiss',rmiss)

#< GGraph による 描画 >
vpx1 = 0.1; vpx2 = 0.4; vpy1 = 0.15; vpy2 = 0.8
y1 = 1050.0; y2 = 1.0e-4
#   Venus
x1 = 150  ; x2 = 750
y1 = 100.0; y2 = 1.0e-5
#   Mars
x1 = 120  ; x2 = 220
y1 = 1.0e3; y2 = 1.0e-3
#
GGraph.set_fig 'itr'=>2, 'viewport'=>[vpx1,vpx2,vpy1,vpy2], 'yrev'=>'units:Pa', 'window'=>[x1,x2,y1,y2]
#GGraph.set_axes('xlabelint'=>30)
GGraph.line( gpt, true, 'annotate'=>false, 'exchange'=>true )
#
GGraph.set_fig 'itr'=>4, 'viewport'=>[vpx1,vpx2,vpy1,vpy2], 'yrev'=>'units:Pa', 'window'=>[1e-9,1.0e0,y1,y2]
#GGraph.set_axes('xlabelint'=>30)
for m in 0..gpmolnum.val.size-1
  if m == 0 then
    flagfirst = true
  else
    flagfirst = false
  end
  p gpmolnum.val[m]
#  GGraph.line( gpgas.cut('molnum'=>molnum[m-1])+1e-10, flagfirst, 'annotate'=>false, 'exchange'=>true, 'type'=>m, 'legend'=>gpmolname.cut('molnum'=>molnum[m-1]).val.to_s, 'legend_vx'=>0.24 )
  GGraph.line( gpgas.cut('MolNum'=>gpmolnum.val[m])+1e-10, flagfirst, 'annotate'=>false, 'exchange'=>true, 'index'=>(m+1)*10+1, 'legend'=>hitran_molname[gpmolnum.val[m]], 'legend_vx'=>0.24 )
end
#
# cloud mass mixing ratio
GGraph.set_fig 'itr'=>4, 'viewport'=>[vpx1,vpx2,vpy1,vpy2], 'yrev'=>'units:Pa', 'window'=>[1e-8,1e-5,y1,y2]
# cloud number density
#GGraph.set_fig 'itr'=>4, 'viewport'=>[vpx1,vpx2,vpy1,vpy2], 'yrev'=>'units:Pa', 'window'=>[1e5,1e10,y1,y2]
#GGraph.set_axes('xlabelint'=>30)
for m in 0..gpcldnum.val.size-1
  if m == 0 then
    flagfirst = true
  else
    flagfirst = false
  end
#  GGraph.line( gpcld.cut('PtclNum'=>m)+1e-10, flagfirst, 'annotate'=>false, 'exchange'=>true, 'index'=>m*10+1, 'legend'=>gpcldname.cut('PtclNum'=>m).val.to_s, 'legend_vx'=>0.2 )
  GGraph.line( gpcld.cut('PtclNum'=>gpcldnum.val[m])+1e-10, flagfirst, 'annotate'=>false, 'exchange'=>true, 'index'=>(m+1)*10+1, 'legend'=>cldname[m], 'legend_vx'=>0.2 )
end

DCL.grcls