function cgp (args)

cs=subwrd(args,1)
if (cs='-h')
  say ' '
  say ' Plots the Adimensional Group Velocity for the Horizontal Modes'
  say ' (Rossby + Mixed, Gravity West and Grevity East + Kelvin)'
  say ' in Function of the Zonal Wave Number'
  say ' '
  say ' Usage:'
  say ' cgp cs var [Mmax]'
  say ' '
  say ' Default: cs = br ; var = rb ; Mmax = 126'
  say ' '
  say ' Must be runned under grads -p (portrait)'
  say ' '
  say ' args meaning:'
  say ' '
  say ' cs : br (barotropic mode) or bc baroclinic mode closest to 250 m'
  say ' var : rb (Rossby + Mixed), gw (gravity west) or ge (gravity east + Kelvin)'
  say ' Mmax : Rhomb+1 (even), Rhomb is the zonal wave truncation (Rhomb odd number)'
  say '        is the total number of zonal waves plus the zonal mean (+1)'
  say ' '
  return
endif
if (cs='') ; cs='br' ; endif
if (cs!='br' & cs!='bc')
  say 'Wrong cs: 'cs' ; must be br or bc'
  return
endif
var=subwrd(args,2)
if (var='') ; var='rb' ; endif
if (var!='rb' & var!='gw' & var!='ge')
  say 'Wrong var: 'var' ; must be rb or gw or ge'
  return
endif
Mmax=subwrd(args,3)
if (Mmax='') ; Mmax=126 ; endif
mrb=Mmax-1
if (mrb < 100)
  Rhomb='RB00'mrb
endif
if (mrb >= 100 & mrb < 1000)
  Rhomb='RB0'mrb
endif
if (mrb >= 1000)
  Rhomb='RB'mrb
endif

rec=chkmdpt()
if (rec!=0)
  return
endif

* Getting Kmax (number of pressure levels) and csvm (reference case)
'run ../../Config/Case_Dates.gs 'MONAN' 1 'ZCAS' 'no
filecsdt='./Case_Dates.dgs'
rec=read(filecsdt)
lrec=sublin(rec,2)
Kmax=subwrd(lrec,1)
csvm=subwrd(lrec,5)
rec=close(filecsdt)

say ' '
say ' csvm = 'csvm'  Kmax = 'Kmax'  cs = 'cs'  var = 'var'  Mmax = 'Mmax
say ' '

if (cs='br')
  hh=01
  hl=00
endif
if (cs='bc')
  hh=08
  hl=07
endif
Trunc=Rhomb'L0'Kmax

if (var='rb') ; tit='Rossby - Mixed' ; endif
if (var='gw') ; tit='Gravity West' ; endif
if (var='ge') ; tit='Gravity East - Kelvin' ; endif

fcg='dir.cfg'
rec=read(fcg)
lrec=sublin(rec,2)
diri=subwrd(lrec,1)
rec=read(fcg)
lrec=sublin(rec,2)
dirfig=subwrd(lrec,1)%'png_'csvm'/'
rec=close(fcg)

filei='HF_'Trunc'_'csvm'_h0'hh'.ctl'

say ' '
say ' filei = 'filei
say ' '

'reinit'

dfilei=diri''filei
rec=chkopfl(dfilei)
if (rec!=0)
  return
endif

'open 'dfilei

if (cs='br')
  if (var='rb') ; cga=-0.11 ; cgb=0.21 ; endif
  if (var='gw') ; cga=-0.36 ; cgb=0.33 ; endif
  if (var='ge') ; cga=0.235 ; cgb=0.365 ; endif
endif
if (cs='bc')
  if (var='rb') ; cga=-0.02 ; cgb=0.028 ; endif
  if (var='gw') ; cga=-0.05 ; cgb=0.017 ; endif
  if (var='ge') ; cga=0.005 ; cgb=0.058 ; endif
endif

'query ctlinfo'
k=1
while (k<=100)
  leh=sublin(result,k)
  if (subwrd(leh,1)='title')
    eps=subwrd(leh,5)
    heq=subwrd(leh,7)
    break
  endif
  k=k+1
endwhile
if (heq >= 1000) ; fmt='%-6.1f' ; endif
if (heq >= 100 & heq < 1000) ; fmt='%-6.2f' ; endif
if (heq >= 10 & heq < 100) ; fmt='%-6.3f' ; endif
if (heq <  10) ; fmt='%-6.4f' ; endif
heqi=math_format(fmt,heq)
say ' 'eps' 'heqi
figname='cgp_'var'_'Trunc'_'csvm'_h0'hh

'set display grey white'
'clear'

'set grads off'
'set mproj off'

'set parea 1.25 8 1 10'
'set xlopts 1 8 0.14'
'set ylopts 1 8 0.14'

'set vrange 'cga' 'cgb
'set x 1 17'
'set xlopts 1 8 0.15'
'set ylopts 1 8 0.15'
'set clopts 1 8 0.125'
y=1
while (y<=11)
'set y 'y
'set cmark 0'
'set ccolor 1'
'set cthick 8'
'display cgp'var
y=y+1
endwhile
'draw title `1Group Velocity 'tit':\`2E'hl'`1 = 'eps'   `2H'hl'`1 = 'heqi' m - 'csvm
'draw xlab `1Zonal Wave Number'
'draw ylab `1Group Velocity'

'query gxinfo'
lx=sublin(result,3)
xa=subwrd(lx,4)
xb=subwrd(lx,6)
ly=sublin(result,4)
ya=subwrd(ly,4)
yb=subwrd(ly,6)

'set line 1 1 6'
'set string 1 c 6 0'

if (var='rb')
  xm=xa+0.15*(xb-xa)
  x1=xm-0.4
  x2=xm+0.4
  'set x 2'
  'set y 3'
  'display cgprb'
  cgrb=subwrd(result,4)
  ym=ya+(yb-ya)*(cgrb-cga)/(cgb-cga)
  y1=ym-0.15
  y2=ym+0.15
  'set line 0 1 8'
  'draw recf 'x1' 'y1' 'x2' 'y2
  'set line 1 1 8'
  'draw rec 'x1' 'y1' 'x2' 'y2
  'draw string 'xm' 'ym' `1Rossby'

  'set x 4'
  'set y 1'
  'display cgprb'
  cgmx=subwrd(result,4)
  xm=xa+0.25*(xb-xa)
  x1=xm-0.3
  x2=xm+0.3
  ym=ya+(yb-ya)*(cgmx-cga)/(cgb-cga)
  y1=ym-0.15
  y2=ym+0.15
  'set line 0 1 8'
  'draw recf 'x1' 'y1' 'x2' 'y2
  'set line 1 1 8'
  'draw rec 'x1' 'y1' 'x2' 'y2
  'draw string 'xm' 'ym' `1Mixed'
endif

if (var='gw')
xm=xa+0.15*(xb-xa)
x1=xm-0.75
x2=xm+0.75
'set x 2'
'set y 4'
'display cgpgw'
cggw=subwrd(result,4)
ym=ya+(yb-ya)*(cggw-cga)/(cgb-cga)
y1=ym-0.15
y2=ym+0.15
'set line 0 1 8'
'draw recf 'x1' 'y1' 'x2' 'y2
'set line 1 1 8'
'draw rec 'x1' 'y1' 'x2' 'y2
'draw string 'xm' 'ym' `1Gravity West'
endif

if (var='ge')
  xm=xa+0.15*(xb-xa)
  x1=xm-0.75
  x2=xm+0.75
  'set x 2'
  'set y 4'
  'display cgpge'
  cgge=subwrd(result,4)
  ym=ya+(yb-ya)*(cgge-cga)/(cgb-cga)
  y1=ym-0.15
  y2=ym+0.15
  'set line 0 1 8'
  'draw recf 'x1' 'y1' 'x2' 'y2
  'set line 1 1 8'
  'draw rec 'x1' 'y1' 'x2' 'y2
  'draw string 'xm' 'ym' `1Gravity East'

  'set x 4'
  'set y 1'
  'display cgpge'
  cgkv=subwrd(result,4)
  xm=xa+0.25*(xb-xa)
  x1=xm-0.3
  x2=xm+0.3
  ym=ya+(yb-ya)*(cgkv-cga)/(cgb-cga)
  y1=ym-0.15
  y2=ym+0.15
  'set line 0 1 8'
  'draw recf 'x1' 'y1' 'x2' 'y2
  'set line 1 1 8'
  'draw rec 'x1' 'y1' 'x2' 'y2
  'draw string 'xm' 'ym' `1Kelvin'
endif

xs=xb
ys=ya-0.6
xrb=xs+0.1
yrb=ys+0.25
xra=xrb-1.8
yra=yrb-0.35
'set line 0 1 6'
'draw recf 'xra' 'yra' 'xrb' 'yrb
'set line 1 1 6'
'draw rec 'xra' 'yra' 'xrb' 'yrb
'set strsiz 0.16'
'set string 1 br 6 0'
'draw string 'xs' 'ys' `1'Trunc

fig=dirfig''figname'.png'
say ' '
say ' 'figname
say ' '
'printim 'fig' png white'

return

function chkmdpt()
'query gxinfo'
lm=sublin(result,2)
xb=subwrd(lm,4)
yb=subwrd(lm,6)
rec=0
if (yb < xb)
  say ' '
  say ' Error: must use grads -p (portrait)'
  say ' '
  rec=1
endif
return rec

function chkopfl(file)
res=read(file)
rec=sublin(res,1)
if (rec=0)
   lrec=sublin(res,2)
   dst=subwrd(lrec,1)
   if (dst!='dset') ; rec=1 ; endif
endif
if (rec!=0)
  say ' '
  say ' Error : 'file' does not exist or is not a ctl file'
  say ' '
endif
if (rec=0)
  res=close(file)
endif
return rec