* +-======-+ 
*  Copyright (c) 2003-2007 United States Government as represented by 
*  the Admistrator of the National Aeronautics and Space Administration.  
*  All Rights Reserved.
*  
*  THIS OPEN  SOURCE  AGREEMENT  ("AGREEMENT") DEFINES  THE  RIGHTS  OF USE,
*  REPRODUCTION,  DISTRIBUTION,  MODIFICATION AND REDISTRIBUTION OF CERTAIN 
*  COMPUTER SOFTWARE ORIGINALLY RELEASED BY THE UNITED STATES GOVERNMENT AS 
*  REPRESENTED BY THE GOVERNMENT AGENCY LISTED BELOW ("GOVERNMENT AGENCY").  
*  THE UNITED STATES GOVERNMENT, AS REPRESENTED BY GOVERNMENT AGENCY, IS AN 
*  INTENDED  THIRD-PARTY  BENEFICIARY  OF  ALL  SUBSEQUENT DISTRIBUTIONS OR 
*  REDISTRIBUTIONS  OF THE  SUBJECT  SOFTWARE.  ANYONE WHO USES, REPRODUCES, 
*  DISTRIBUTES, MODIFIES  OR REDISTRIBUTES THE SUBJECT SOFTWARE, AS DEFINED 
*  HEREIN, OR ANY PART THEREOF,  IS,  BY THAT ACTION, ACCEPTING IN FULL THE 
*  RESPONSIBILITIES AND OBLIGATIONS CONTAINED IN THIS AGREEMENT.
*  
*  Government Agency: National Aeronautics and Space Administration
*  Government Agency Original Software Designation: GSC-15354-1
*  Government Agency Original Software Title:  GEOS-5 GCM Modeling Software
*  User Registration Requested.  Please Visit http://opensource.gsfc.nasa.gov
*  Government Agency Point of Contact for Original Software:  
*  			Dale Hithon, SRA Assistant, (301) 286-2691
*  
* +-======-+ 
function makequad (args)

'numargs  'args
 numargs = result

       num =  0
while( num < numargs )
       num = num + 1
if( subwrd(args,num)='-U'     ) ;  uname  = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-V'     ) ;  vname  = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-T'     ) ;  tname  = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-Q'     ) ;  qname  = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-H'     ) ;  hname  = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-W'     ) ;  wname  = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-UV'    ) ; uvname = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-VT'    ) ; vtname = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-VQ'    ) ; vqname = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-WT'    ) ; wtname = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-WQ'    ) ; wqname = subwrd(args,num+1) ; endif
if( subwrd(args,num)='-TAG'   ) ; tag    = subwrd(args,num+1) ; endif
endwhile

name = 'quad'tag

* Variable is in File
* -------------------
    say ''
    say 'Computing Quadratic Data:'
    say '-------------------------'

'getinfo file'
       infile = result

'q gxout'
   gxout = sublin(result,4)
   gxout = subwrd(gxout ,6)

'run getenv "GEOSUTIL"'
             geosutil = result

'getinfo lonmin'
         lonmin = result
'getinfo lonmax'
         lonmax = result

'getinfo latmin'
         latmin = result
'getinfo latmax'
         latmax = result

'getinfo tmin'
         tmin = result
'getinfo tmax'
         tmax = result

'set t 'tmin
'getinfo date'
      begdate = result

'getinfo zmin'
         zmin = result
'getinfo zmax'
         zmax = result
'set z  'zmin' 'zmax
         zdim = zmax-zmin+1

'sety'
'set t 'tmin' 'tmax
        tdim = tmax-tmin+1

'getinfo xdim'
         xdim = result
'getinfo ydim'
         ydim = result

'getinfo undef'
         undef = result
'getinfo dlat' 
         dlat  = result
'getinfo dlon' 
         dlon  = result
       lontot  = lonmax - lonmin + dlon/2

* Get INFO from CTL
* -----------------
'q ctlinfo'
         n = 1
      while( n >0 )
             line = sublin(result,n)
             word = subwrd(line,1)
         if( word = 'tdef' )
             dt = subwrd(line,5)
             n  = 0
          else
             n  = n + 1
          endif
      endwhile
         n = 1
      while( n >0 )
             line = sublin(result,n)
             word = subwrd(line,1)
         if( word = 'ydef' )
           lat0 = subwrd(line,4)
             n  = 0
          else
             n  = n + 1
          endif
      endwhile
         n = 1
      while( n >0 )
             line = sublin(result,n)
             word = subwrd(line,1)
         if( word = 'xdef' )
           lon0 = subwrd(line,4)
             n  = 0
          else
             n  = n + 1
          endif
      endwhile
* -----------------

      levs = ''
             z  = zmin
      while( z <= zmax )
            'set z 'z
            'getinfo level'
             levs = levs' 'result
             z = z + 1
      endwhile
       
    say ''
    say 'DFILE: 'infile
    say 'UNDEF: 'undef
    say ' LON0: 'lon0
    say ' DLON: 'dlon
    say ' YDEF: 'ydim
    say ' LAT0: 'lat0
    say ' DLAT: 'dlat
    say ' XDEF: 'xdim
    say ' ZDEF: 'zdim
    say ' TDEF: 'tdim
    say 'BDATE: 'begdate
    say '   DT: 'dt
    say ' LEVS: 'levs
    say ''
    say ' LONMIN:LONMAX = 'lonmin':'lonmax
    say ''

  '!remove sedfile'
  '!touch  sedfile'
  '!remove quad.ctl'
  '!remove 'name'.data'
   'set gxout fwrite'
   'set       fwrite 'name'.data'

* ----------------------------------------------------------------------------

    t = tmin
    while( t<=tmax )
    say '          Computing Quadratic Data for Tag: 'tag'  for t = 't
   'set t 't
   'sety'
   'set z 'zmin' 'zmax

   'set x 1 'xdim
       'define  qundef = lon-lon + lat-lat + 'undef

    if( uname = NULL )
       'define  locu   =  qundef'
    else
       'define  locu   = 'uname'*u'tag'scale'
    endif
    if( vname = NULL )
       'define  locv   =  qundef'
    else
       'define  locv   = 'vname'*v'tag'scale'
    endif
    if( tname = NULL )
       'define  loct   =  qundef'
    else
       'define  loct   = 'tname'*t'tag'scale'
    endif
    if( qname = NULL )
       'define  locq   =  qundef'
    else
       'define  locq   = 'qname'*q'tag'scale'
    endif
    if( hname = NULL )
       'define  loch   =  qundef'
    else
       'define  loch   = 'hname'*h'tag'scale'
    endif
    if( wname = NULL )
       'define  locw   =  qundef'
    else
       'define  locw   = 'wname'*w'tag'scale'
    endif

   'set x 1'
   'define locuz = ave(locu,x=1,x='xdim')'
   'define locvz = ave(locv,x=1,x='xdim')'
   'define loctz = ave(loct,x=1,x='xdim')'
   'define locqz = ave(locq,x=1,x='xdim')'
   'define lochz = ave(loch,x=1,x='xdim')'
   'define locwz = ave(locw,x=1,x='xdim')'

   'set x 1 'xdim
   'd (locu-locuz)*(locu-locuz)'
   'd (locv-locvz)*(locv-locvz)'
   'd (loct-loctz)*(loct-loctz)'
   'd (locq-locqz)*(locq-locqz)'
   'd (locw-locwz)*(locw-locwz)'
   'd (loch-lochz)*(loch-lochz)'

   'd (locu-locuz)*(locv-locvz)'
   'd (locv-locvz)*(loct-loctz)'
   'd (locv-locvz)*(locq-locqz)'
   'd (locw-locwz)*(loct-loctz)'
   'd (locw-locwz)*(locq-locqz)'

    if( uname = NULL )
       'd locu'
    else
       'd Var_'uname'*pow(u'tag'scale,2)'
    endif
    if( vname = NULL )
       'd locv'
    else
       'd Var_'vname'*pow(v'tag'scale,2)'
    endif
    if( tname = NULL )
       'd loct'
    else
       'd Var_'tname'*pow(t'tag'scale,2)'
    endif
    if( qname = NULL )
       'd locq'
    else
       'd Var_'qname'*pow(q'tag'scale,2)'
    endif
    if( wname = NULL )
       'd locw'
    else
       'd Var_'wname'*pow(w'tag'scale,2)'
    endif
    if( hname = NULL )
       'd loch'
    else
       'd Var_'hname'*pow(h'tag'scale,2)'
    endif

    if( uname = 'NULL' | vname = 'NULL' )
       'd  qundef'
    else
       'd 'uvname'*u'tag'scale*v'tag'scale'
    endif
    if( vname = 'NULL' | tname = 'NULL' )
       'd  qundef'
    else
       'd 'vtname'*v'tag'scale*t'tag'scale'
    endif
    if( vname = 'NULL' | qname = 'NULL' )
       'd  qundef'
    else
       'd 'vqname'*v'tag'scale*q'tag'scale'
    endif
    if( wname = 'NULL' | tname = 'NULL' )
       'd  qundef'
    else
       'd 'wtname'*w'tag'scale*t'tag'scale'
    endif
    if( wname = 'NULL' | qname = 'NULL' )
       'd  qundef'
    else
       'd 'wqname'*w'tag'scale*q'tag'scale'
    endif

    t = t+1
    endwhile

* ----------------------------------------------------------------------------

    say ''
   'disable fwrite'
   'undefine locu'
   'undefine locv'
   'undefine loct'
   'undefine locq'
   'undefine locw'
   'undefine loch'

   'undefine locuz'
   'undefine locvz'
   'undefine loctz'
   'undefine locqz'
   'undefine locwz'
   'undefine lochz'

   'set gxout 'gxout

'!echo "s?q.data?"'name'.data?g >> sedfile'
'!echo "s?UNDEF?"'undef'?g >> sedfile'
'!echo "s?DX?"'dlon'?g >> sedfile'
'!echo "s?DY?"'dlat'?g >> sedfile'
'!echo "s?DT?"'dt'?g >> sedfile'
'!echo "s?LON0?"'lon0'?g >> sedfile'
'!echo "s?LAT0?"'lat0'?g >> sedfile'
'!echo "s?XDIM?"'xdim'?g >> sedfile'
'!echo "s?YDIM?"'ydim'?g >> sedfile'
'!echo "s?ZDIM?"'zdim'?g >> sedfile'
'!echo "s?TDIM?"'tdim'?g >> sedfile'
'!echo "s?BDATE?"'begdate'?g >> sedfile'
'!echo "s?LEVS?"'levs'?g  >> sedfile'
'!/bin/cp 'geosutil'/plots/grads_util/makequad.tmpl .'
'!sed -f   sedfile makequad.tmpl > 'name'.ctl'

   'open 'name'.ctl'
   'getinfo    numfiles'
               newfile = result
   'set dfile 'newfile
   'set lon 'lonmin' 'lonmax
   'sety'
   'setz'
   'getinfo  tdim'
             tdim = result
   'set t 1 'tdim

   'define Var'uname's'tag' = varus'
   'define Var'vname's'tag' = varvs'
   'define Var'tname's'tag' = varts'
   'define Var'qname's'tag' = varqs'
   'define Var'wname's'tag' = varws'
   'define Var'hname's'tag' = varhs'

   'define Cov'uname''vname's'tag' = covuvs'
   'define Cov'vname''tname's'tag' = covvts'
   'define Cov'vname''qname's'tag' = covvqs'
   'define Cov'wname''tname's'tag' = covwts'
   'define Cov'wname''qname's'tag' = covwqs'

   'define Var'uname'p'tag' = varup'
   'define Var'vname'p'tag' = varvp'
   'define Var'tname'p'tag' = vartp'
   'define Var'qname'p'tag' = varqp'
   'define Var'wname'p'tag' = varwp'
   'define Var'hname'p'tag' = varhp'

   'define Cov'uname''vname'p'tag' = covuvp'
   'define Cov'vname''tname'p'tag' = covvtp'
   'define Cov'vname''qname'p'tag' = covvqp'
   'define Cov'wname''tname'p'tag' = covwtp'
   'define Cov'wname''qname'p'tag' = covwqp'

   'close 'newfile

   'set dfile 'infile
   'set lat 'latmin' 'latmax
   'set lon 'lonmin' 'lonmax
   'set t   '  tmin' '  tmax

return