function zquad (args)

*******************************************************
****                 INPUT Variables               ****
*******************************************************

'numargs  'args
 numargs = result

DEBUG  = FALSE

        n   = 0
        num = 0
while ( num < numargs )
        num = num + 1

if( subwrd(args,num) = '-EXPID'  ) ; EXPID  = subwrd(args,num+1) ; endif
if( subwrd(args,num) = '-OUTPUT' ) ; OUTPUT = subwrd(args,num+1) ; endif
if( subwrd(args,num) = '-DEBUG'  ) ; DEBUG  = subwrd(args,num+1) ; endif

* Read SEASONS
* -----------
if( subwrd(args,num) = '-SEASON' )
     seasons = ''
           k = 1
    while( k > 0 )
           L = num + k
        season  = subwrd(args,L)
    if( season  = '' )
        k = -1
    else
        bit = substr(season,1,1)
        if( bit = '-' )
              k = -1
        else
              seasons = seasons % ' ' % season
              k = k+1
        endif
    endif
    endwhile
endif

endwhile
'uppercase 'seasons
            seasons = result
'numargs  ' seasons
 numseas =  result



* Initialize
* ----------
"reinit"
'set display color white'
'set csmooth on'
'c'


* Get Environment Variables
* -------------------------

'getenv "GEOSUTIL"'
         geosutil = result
'getenv "VERIFICATION"'
         verification = result


* Get Model Primary Fields
* ------------------------
'getvar PS DYN'
      pmname  = subwrd(result,1)
      pmfile  = subwrd(result,2)
      pmscale = subwrd(result,3)
  'd 'pmscale
      pmscale = subwrd(result,4)


'getvar U DYN'
      umname  = subwrd(result,1)
      umfile  = subwrd(result,2)
      umscale = subwrd(result,3)
       mdesc  = subwrd(result,4)
  'd 'umscale
      umscale = subwrd(result,4)

'getvar V DYN'
      vmname  = subwrd(result,1)
      vmfile  = subwrd(result,2)
      vmscale = subwrd(result,3)
  'd 'vmscale
      vmscale = subwrd(result,4)

'getvar T DYN'
      tmname  = subwrd(result,1)
      tmfile  = subwrd(result,2)
      tmscale = subwrd(result,3)
  'd 'tmscale
      tmscale = subwrd(result,4)

'getvar ZLE DYN'
      hmname  = subwrd(result,1)
      hmfile  = subwrd(result,2)
      hmscale = subwrd(result,3)
  'd 'hmscale
      hmscale = subwrd(result,4)

'getvar OMEGA DYN'
      wmname  = subwrd(result,1)
      wmfile  = subwrd(result,2)
      wmscale = subwrd(result,3)
  'd 'wmscale
      wmscale = subwrd(result,4)

'getvar Q MOIST'
      qmname  = subwrd(result,1)
      qmfile  = subwrd(result,2)
      qmscale = subwrd(result,3)
  'd 'qmscale
      qmscale = subwrd(result,4)

'define umscale = 'umscale
'define vmscale = 'vmscale
'define tmscale = 'tmscale
'define hmscale = 'hmscale
'define wmscale = 'wmscale
'define qmscale = 'qmscale

* Set Model Variance Names
* ------------------------
 result = setname (umname,umname,umfile,m)
          uumname  = subwrd(result,1)
          uumfile  = umfile
          uumscale = umscale*umscale

 result = setname (vmname,vmname,vmfile,m)
          vvmname  = subwrd(result,1)
          vvmfile  = vmfile
          vvmscale = vmscale*vmscale

 result = setname (tmname,tmname,tmfile,m)
          ttmname  = subwrd(result,1)
          ttmfile  = tmfile
          ttmscale = tmscale*tmscale

 result = setname (hmname,hmname,hmfile,m)
          hhmname  = subwrd(result,1)
          hhmfile  = hmfile
          hhmscale = hmscale*hmscale

 result = setname (wmname,wmname,wmfile,m)
          wwmname  = subwrd(result,1)
          wwmfile  = wmfile
          wwmscale = wmscale*wmscale

 result = setname (qmname,qmname,qmfile,m)
          qqmname  = subwrd(result,1)
          qqmfile  = qmfile
          qqmscale = qmscale*qmscale


* Set Model CoVariance Names
* --------------------------
 result = setname (umname,vmname,umfile,m)
          uvmname  = subwrd(result,1)
          uvmfile  = umfile
          uvmscale = umscale*vmscale

 result = setname (vmname,tmname,vmfile,m)
          vtmname  = subwrd(result,1)
          vtmfile  = vmfile
          vtmscale = vmscale*tmscale

 result = setname (vmname,qmname,vmfile,m)
          vqmname  = subwrd(result,1)
          vqmfile  = vmfile
          vqmscale = vmscale*qmscale

 result = setname (wmname,tmname,wmfile,m)
          wtmname  = subwrd(result,1)
          wtmfile  = wmfile
          wtmscale = wmscale*tmscale

 result = setname (wmname,qmname,wmfile,m)
          wqmname  = subwrd(result,1)
          wqmfile  = wmfile
          wqmscale = wmscale*qmscale

'define uvmscale = 'uvmscale
'define vtmscale = 'vtmscale
'define vqmscale = 'vqmscale
'define wtmscale = 'wtmscale
'define wqmscale = 'wqmscale

**************************************************
****            Compute Quadratics            ****
**************************************************

* Compute Zonal Means of PS Field
* -------------------------------
   'set dfile 'pmfile
   'setlons'
   'setdates'
   'set lat -90 90'
   'getdates'
   'set z 1'
   'run   setenv "ZDIM" 1'
   'makezf  'pmname' psmod z'
   'run unsetenv "ZDIM"'
   'set lon 0'
   'seasonal psmodz'
   'undefine psmodz'

* Compute Quadratics of Remaining Fields
* --------------------------------------
   'set dfile 'umfile
   'setlons'
   'setdates'
   'run getenv "BEGDATE"'
                begdate  = result
   'run getenv "ENDDATE"'
                enddate  = result
   'set lat -90 90'
   'getdates'
   'setz'

   'makequad -q1 'umname' -file1 'umfile' -scale1 'umscale' -q2 'umname' -file2 'umfile' -scale2 'umscale' -qq 'uumname' -tag m -zonal'
   'makequad -q1 'vmname' -file1 'vmfile' -scale1 'vmscale' -q2 'vmname' -file2 'vmfile' -scale2 'vmscale' -qq 'vvmname' -tag m -zonal'
   'makequad -q1 'tmname' -file1 'tmfile' -scale1 'tmscale' -q2 'tmname' -file2 'tmfile' -scale2 'tmscale' -qq 'ttmname' -tag m -zonal'
   'makequad -q1 'hmname' -file1 'hmfile' -scale1 'hmscale' -q2 'hmname' -file2 'hmfile' -scale2 'hmscale' -qq 'hhmname' -tag m -zonal'
   'makequad -q1 'wmname' -file1 'wmfile' -scale1 'wmscale' -q2 'wmname' -file2 'wmfile' -scale2 'wmscale' -qq 'wwmname' -tag m -zonal'
   'makequad -q1 'qmname' -file1 'qmfile' -scale1 'qmscale' -q2 'qmname' -file2 'qmfile' -scale2 'qmscale' -qq 'qqmname' -tag m -zonal'

   'makequad -q1 'umname' -file1 'umfile' -scale1 'umscale' -q2 'vmname' -file2 'vmfile' -scale2 'vmscale' -qq 'uvmname' -tag m -zonal'
   'makequad -q1 'vmname' -file1 'vmfile' -scale1 'vmscale' -q2 'tmname' -file2 'tmfile' -scale2 'tmscale' -qq 'vtmname' -tag m -zonal'
   'makequad -q1 'vmname' -file1 'vmfile' -scale1 'vmscale' -q2 'qmname' -file2 'qmfile' -scale2 'qmscale' -qq 'vqmname' -tag m -zonal'
   'makequad -q1 'wmname' -file1 'wmfile' -scale1 'wmscale' -q2 'tmname' -file2 'tmfile' -scale2 'tmscale' -qq 'wtmname' -tag m -zonal'
   'makequad -q1 'wmname' -file1 'wmfile' -scale1 'wmscale' -q2 'qmname' -file2 'qmfile' -scale2 'qmscale' -qq 'wqmname' -tag m -zonal'


*********************************************************************
****    Loop over Possible Experiment Datasets for Comparison    ****
*********************************************************************

'!/bin/mv HISTORY.T HISTORY.Tmp'
'getenv "CMPEXP"'
         cmpexp = result
         numexp = 1
            exp = subwrd( cmpexp,numexp )
while( exp != 'NULL' )
say ' '
say 'Comparing with: 'exp

'!chckfile 'exp'/.HOMDIR' 
 'getenv CHECKFILE'
         CHECKFILE  = result
     if( CHECKFILE != 'NULL' )
        '!/bin/cp `cat 'exp'/.HOMDIR`/HISTORY.rc .'
     else
        '!/bin/cp 'exp'/HISTORY.rc .'
     endif
'!remove CHECKFILE.txt'
'!cat HISTORY.rc | sed -e "s/,/ , /g" > HISTORY.T'

* Get CMPEXP Primary Fields
* -------------------------
'set dfile 1'
'set x 1'
'set y 1'
'set z 1'
'set t 1'

'getvar PS DYN 'exp
        poname = subwrd(result,1)
        pofile = subwrd(result,2)

'getvar U DYN 'exp
say 'RESULT: 'result
        uoname  = subwrd(result,1)
        uofile  = subwrd(result,2)
        uoscale = subwrd(result,3)
         odesc  = subwrd(result,4)
         obsid  = subwrd(result,5)
    'd 'uoscale
        uoscale = subwrd(result,4)

say 'ODESC: 'odesc

'getvar V DYN 'exp
        voname  = subwrd(result,1)
        vofile  = subwrd(result,2)
        voscale = subwrd(result,3)
    'd 'voscale
        voscale = subwrd(result,4)

'getvar T DYN 'exp
        toname  = subwrd(result,1)
        tofile  = subwrd(result,2)
        toscale = subwrd(result,3)
    'd 'toscale
        toscale = subwrd(result,4)

'getvar ZLE DYN 'exp
        honame  = subwrd(result,1)
        hofile  = subwrd(result,2)
        hoscale = subwrd(result,3)
    'd 'hoscale
        hoscale = subwrd(result,4)

'getvar OMEGA DYN 'exp
        woname  = subwrd(result,1)
        wofile  = subwrd(result,2)
        woscale = subwrd(result,3)
    'd 'woscale
        woscale = subwrd(result,4)

'getvar Q MOIST 'exp
        qoname  = subwrd(result,1)
        qofile  = subwrd(result,2)
        qoscale = subwrd(result,3)
    'd 'qoscale
        qoscale = subwrd(result,4)

'define uoscale = 'uoscale
'define voscale = 'voscale
'define toscale = 'toscale
'define hoscale = 'hoscale
'define woscale = 'woscale
'define qoscale = 'qoscale

* Set CMPEXP Variance Names
* -------------------------
 result = setname (uoname,uoname,uofile,o)
          uuoname  = subwrd(result,1)
          uuofile  = uofile
          uuoscale = uoscale*uoscale

 result = setname (voname,voname,vofile,o)
          vvoname  = subwrd(result,1)
          vvofile  = vofile
          vvoscale = voscale*voscale

 result = setname (toname,toname,tofile,o)
          ttoname  = subwrd(result,1)
          ttofile  = tofile
          ttoscale = toscale*toscale

 result = setname (honame,honame,hofile,o)
          hhoname  = subwrd(result,1)
          hhofile  = hofile
          hhoscale = hoscale*hoscale

 result = setname (woname,woname,wofile,o)
          wwoname  = subwrd(result,1)
          wwofile  = wofile
          wwoscale = woscale*woscale

 result = setname (qoname,qoname,qofile,o)
          qqoname  = subwrd(result,1)
          qqofile  = qofile
          qqoscale = qoscale*qoscale


* Set CMPEXP CoVariance Names
* ---------------------------
 result = setname (uoname,voname,uofile,o)
          uvoname  = subwrd(result,1)
          uvofile  = uofile
          uvoscale = uoscale*voscale

 result = setname (voname,toname,vofile,o)
          vtoname  = subwrd(result,1)
          vtofile  = vofile
          vtoscale = voscale*toscale

 result = setname (voname,qoname,vofile,o)
          vqoname  = subwrd(result,1)
          vqofile  = vofile
          vqoscale = voscale*qoscale

 result = setname (woname,toname,wofile,o)
          wtoname  = subwrd(result,1)
          wtofile  = wofile
          wtoscale = woscale*toscale

 result = setname (woname,qoname,wofile,o)
          wqoname  = subwrd(result,1)
          wqofile  = wofile
          wqoscale = woscale*qoscale

'define uvoscale = 'uvoscale
'define vtoscale = 'vtoscale
'define vqoscale = 'vqoscale
'define wtoscale = 'wtoscale
'define wqoscale = 'wqoscale


* Compute Zonal Means of PS Field
* -------------------------------
if( poname != 'NULL' )
   'set dfile 'pofile
   'setlons'
   'set lat -90 90'
   'getdates'
   'getinfo   tmin'
              tmin = result
   'getinfo   tmax'
              tmax = result
   'set t    'tmin
   'getinfo   date'
              begdateo = result
   'set t    'tmax
   'getinfo   date'
              enddateo = result
   'set time 'begdateo' 'enddateo
   'set z 1'
   'run   setenv "ZDIM" 1'
   'makezf  'poname' psobs z'
   'run unsetenv "ZDIM"'
   'set lon 0'
   'seasonal psobsz'
   'undefine psobsz'
else
   'set dfile 'uofile
   'setlons'
   'set lat -90 90'
   'getdates'
   'getinfo   tmin'
              tmin = result
   'getinfo   tmax'
              tmax = result
   'set t    'tmin
   'getinfo   date'
              begdateo = result
   'set t    'tmax
   'getinfo   date'
              enddateo = result
   'set time 'begdateo' 'enddateo
   'set z 1'
   'run   setenv "ZDIM" 1'
   'define psobs = 100000 + lon-lon + lat-lat'
   'makezf psobs psobs z'
   'run unsetenv "ZDIM"'
   'set lon 0'
   'seasonal psobsz'
   'undefine psobsz'
endif

   'run setenv   "BEGDATEO" 'begdateo
   'run setenv   "ENDDATEO" 'enddateo

* Compute Quadratics of Remaining Fields
* --------------------------------------
   'set dfile 'uofile
   'setlons'
   'set lat -90 90'
   'getdates'
   'getinfo   tmin'
              tmin = result
   'getinfo   tmax'
              tmax = result
   'set t    'tmin
   'getinfo   date'
              begdateo = result
   'set t    'tmax
   'getinfo   date'
              enddateo = result
   'set time 'begdateo' 'enddateo
   'setz'

   'makequad -q1 'uoname' -file1 'uofile' -scale1 'uoscale' -q2 'uoname' -file2 'uofile' -scale2 'uoscale' -qq 'uuoname' -tag o -zonal'
   'makequad -q1 'voname' -file1 'vofile' -scale1 'voscale' -q2 'voname' -file2 'vofile' -scale2 'voscale' -qq 'vvoname' -tag o -zonal'
   'makequad -q1 'toname' -file1 'tofile' -scale1 'toscale' -q2 'toname' -file2 'tofile' -scale2 'toscale' -qq 'ttoname' -tag o -zonal'
   'makequad -q1 'honame' -file1 'hofile' -scale1 'hoscale' -q2 'honame' -file2 'hofile' -scale2 'hoscale' -qq 'hhoname' -tag o -zonal'
   'makequad -q1 'woname' -file1 'wofile' -scale1 'woscale' -q2 'woname' -file2 'wofile' -scale2 'woscale' -qq 'wwoname' -tag o -zonal'
   'makequad -q1 'qoname' -file1 'qofile' -scale1 'qoscale' -q2 'qoname' -file2 'qofile' -scale2 'qoscale' -qq 'qqoname' -tag o -zonal'

   'makequad -q1 'uoname' -file1 'uofile' -scale1 'uoscale' -q2 'voname' -file2 'vofile' -scale2 'voscale' -qq 'uvoname' -tag o -zonal'
   'makequad -q1 'voname' -file1 'vofile' -scale1 'voscale' -q2 'toname' -file2 'tofile' -scale2 'toscale' -qq 'vtoname' -tag o -zonal'
   'makequad -q1 'voname' -file1 'vofile' -scale1 'voscale' -q2 'qoname' -file2 'qofile' -scale2 'qoscale' -qq 'vqoname' -tag o -zonal'
   'makequad -q1 'woname' -file1 'wofile' -scale1 'woscale' -q2 'toname' -file2 'tofile' -scale2 'toscale' -qq 'wtoname' -tag o -zonal'
   'makequad -q1 'woname' -file1 'wofile' -scale1 'woscale' -q2 'qoname' -file2 'qofile' -scale2 'qoscale' -qq 'wqoname' -tag o -zonal'

*******************************************************
****         Loop over Seasons to Process          ****
*******************************************************

ptop.1 = 100
ptop.2 = 1

       m = 1
while( m > 0 )
    season = subwrd(seasons,m)
if( season = '' )
         m = -1
else
         m = m+1
         say 'Processing Season: 'season

'set dfile 'umfile
'count "'season'" 'begdate' 'enddate
 nmod =  result

'set dfile 'uofile
'count "'season'"'
 nobs =  result

* Make Variance Plots
* -------------------
          VARSM = umname' 'vmname' 'tmname' 'qmname' 'wmname' 'hmname
          VARSO = uoname' 'voname' 'toname' 'qoname' 'woname' 'honame
'numargs 'VARSM
          nvars = result

n = 1
while( n<=nvars )
       varmname = subwrd( VARSM,n )
       varoname = subwrd( VARSO,n )

if( varmname != 'NULL' & varoname != 'NULL' )

'uppercase 'varmname
            varname = 'Var_'result

'set dfile 'umfile
'setz'
'define qsqsmz = sqrt(Var'varmname'sm'season')'
'define qpqpmz = sqrt(Var'varmname'pm'season')'
'set dfile 'uofile
'setz'
'define qsqsoz = sqrt(Var'varoname'so'season')'
'define qpqpoz = sqrt(Var'varoname'po'season')'

k = 1
while( k<=2 )
ptop = ptop.k

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/zplot -EXPID 'EXPID' -EXPORT 'varname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -MDESC 'mdesc' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

k = k+1
endwhile

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/lplot -EXPID 'EXPID' -EXPORT 'varname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

endif

n = n+1
endwhile

* Make KE Plots
* -------------
'set dfile 'umfile
'setz'
'define          KEsm'season' = 0.5*( Var'umname'sm'season' + Var'vmname'sm'season' )'
'define          KEpm'season' = 0.5*( Var'umname'pm'season' + Var'vmname'pm'season' )'
'define qsqsmz = KEsm'season
'define qpqpmz = KEpm'season
'set dfile 'uofile
'setz'
'define          KEso'season' = 0.5*( Var'uoname'so'season' + Var'voname'so'season' )'
'define          KEpo'season' = 0.5*( Var'uoname'po'season' + Var'voname'po'season' )'
'define qsqsoz = KEso'season
'define qpqpoz = KEpo'season

k = 1
while( k<=2 )
ptop = ptop.k

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/zplot -EXPID 'EXPID' -EXPORT Var_KE -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -MDESC 'mdesc' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

k = k+1
endwhile

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/lplot -EXPID 'EXPID' -EXPORT Var_KE -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile


* Make CoVariance Plots
* ---------------------
                  nvars = 1
          COVARSM.nvars = umname' 'vmname ; nvars = nvars + 1
          COVARSM.nvars = vmname' 'tmname ; nvars = nvars + 1
          COVARSM.nvars = vmname' 'qmname ; nvars = nvars + 1
          COVARSM.nvars = wmname' 'tmname ; nvars = nvars + 1
          COVARSM.nvars = wmname' 'qmname ; nvars = nvars + 1
                  nvars = nvars - 1

                  nvars = 1
          COVARSO.nvars = uoname' 'voname ; nvars = nvars + 1
          COVARSO.nvars = voname' 'toname ; nvars = nvars + 1
          COVARSO.nvars = voname' 'qoname ; nvars = nvars + 1
          COVARSO.nvars = woname' 'toname ; nvars = nvars + 1
          COVARSO.nvars = woname' 'qoname ; nvars = nvars + 1
                  nvars = nvars - 1

n = 1
while( n<=nvars )
       varmname1 = subwrd( COVARSM.n,1 )
       varmname2 = subwrd( COVARSM.n,2 )
       varoname1 = subwrd( COVARSO.n,1 )
       varoname2 = subwrd( COVARSO.n,2 )

'uppercase 'varmname1
             varname1 = result
'uppercase 'varmname2
             varname2 = result
               cname  = 'Cov_'varname1'_'varname2

'set dfile 'umfile
'setz'
'define qsqsmz = Cov'varmname1''varmname2'sm'season
'define qpqpmz = Cov'varmname1''varmname2'pm'season
'set dfile 'uofile
'setz'
'define qsqsoz = Cov'varoname1''varoname2'so'season
'define qpqpoz = Cov'varoname1''varoname2'po'season

k = 1
while( k<=2 )
ptop = ptop.k

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/zplot -EXPID 'EXPID' -EXPORT 'cname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -MDESC 'mdesc' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

k = k+1
endwhile

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/lplot -EXPID 'EXPID' -EXPORT 'cname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else 
            flag = "next"
        endif
       'c'
   endwhile

n = n+1
endwhile

endif
endwhile

* Check next Comparison Experiment Dataset
* ----------------------------------------
  numexp = numexp + 1
     exp = subwrd( cmpexp,numexp )
endwhile
'!/bin/mv HISTORY.Tmp HISTORY.T'


*********************************************************************
****   Loop over Possible Verification Datasets for Comparison   ****
*********************************************************************

' getnumrc 'geosutil'/plots/quad'
     rcinfo = result
     numrc  = subwrd( rcinfo,1 )
       num  = 1
       cnt  = 0

while( num <= numrc )
        loc = num + 1
     rcfile = subwrd( rcinfo,loc )


* Get OBS Primary Fields
* ----------------------
'set dfile 1'
'set x 1'
'set y 1'
'set z 1'
'set t 1'

'getobs PS DYN 'rcfile
        poname = subwrd(result,1)
        pofile = subwrd(result,2)

'getobs U DYN 'rcfile
        uoname  = subwrd(result,1)
        uofile  = subwrd(result,2)
        uoscale = subwrd(result,3)
         odesc  = subwrd(result,4)
         obsid  = subwrd(result,5)
    'd 'uoscale
        uoscale = subwrd(result,4)

'getobs V DYN 'rcfile
        voname  = subwrd(result,1)
        vofile  = subwrd(result,2)
        voscale = subwrd(result,3)
    'd 'voscale
        voscale = subwrd(result,4)

'getobs T DYN 'rcfile
        toname  = subwrd(result,1)
        tofile  = subwrd(result,2)
        toscale = subwrd(result,3)
    'd 'toscale
        toscale = subwrd(result,4)

'getobs ZLE DYN 'rcfile
        honame  = subwrd(result,1)
        hofile  = subwrd(result,2)
        hoscale = subwrd(result,3)
    'd 'hoscale
        hoscale = subwrd(result,4)

'getobs OMEGA DYN 'rcfile
        woname  = subwrd(result,1)
        wofile  = subwrd(result,2)
        woscale = subwrd(result,3)
    'd 'woscale
        woscale = subwrd(result,4)

'getobs Q MOIST 'rcfile
        qoname  = subwrd(result,1)
        qofile  = subwrd(result,2)
        qoscale = subwrd(result,3)
    'd 'qoscale
        qoscale = subwrd(result,4)

'define uoscale = 'uoscale
'define voscale = 'voscale
'define toscale = 'toscale
'define hoscale = 'hoscale
'define woscale = 'woscale
'define qoscale = 'qoscale

* Set OBS Variance Names
* ----------------------
 result = setname (uoname,uoname,uofile,o)
          uuoname  = subwrd(result,1)
          uuofile  = uofile
          uuoscale = uoscale*uoscale

 result = setname (voname,voname,vofile,o)
          vvoname  = subwrd(result,1)
          vvofile  = vofile
          vvoscale = voscale*voscale

 result = setname (toname,toname,tofile,o)
          ttoname  = subwrd(result,1)
          ttofile  = tofile
          ttoscale = toscale*toscale

 result = setname (honame,honame,hofile,o)
          hhoname  = subwrd(result,1)
          hhofile  = hofile
          hhoscale = hoscale*hoscale

 result = setname (woname,woname,wofile,o)
          wwoname  = subwrd(result,1)
          wwofile  = wofile
          wwoscale = woscale*woscale

 result = setname (qoname,qoname,qofile,o)
          qqoname  = subwrd(result,1)
          qqofile  = qofile
          qqoscale = qoscale*qoscale


* Set OBS CoVariance Names
* ------------------------
 result = setname (uoname,voname,uofile,o)
          uvoname  = subwrd(result,1)
          uvofile  = uofile
          uvoscale = uoscale*voscale

 result = setname (voname,toname,vofile,o)
          vtoname  = subwrd(result,1)
          vtofile  = vofile
          vtoscale = voscale*toscale

 result = setname (voname,qoname,vofile,o)
          vqoname  = subwrd(result,1)
          vqofile  = vofile
          vqoscale = voscale*qoscale

 result = setname (woname,toname,wofile,o)
          wtoname  = subwrd(result,1)
          wtofile  = wofile
          wtoscale = woscale*toscale

 result = setname (woname,qoname,wofile,o)
          wqoname  = subwrd(result,1)
          wqofile  = wofile
          wqoscale = woscale*qoscale

'define uvoscale = 'uvoscale
'define vtoscale = 'vtoscale
'define vqoscale = 'vqoscale
'define wtoscale = 'wtoscale
'define wqoscale = 'wqoscale


* Compute Zonal Means of PS Field
* -------------------------------
if( poname != 'NULL' )
   'set dfile 'pofile
   'setlons'
   'set lat -90 90'
   'getdates'
   'getinfo   tmin'
              tmin = result
   'getinfo   tmax'
              tmax = result
   'set t    'tmin
   'getinfo   date'
              begdateo = result
   'set t    'tmax
   'getinfo   date'
              enddateo = result
   'set time 'begdateo' 'enddateo
   'set z 1'
   'run   setenv "ZDIM" 1'
   'makezf  'poname' psobs z'
   'run unsetenv "ZDIM"'
   'set lon 0'
   'seasonal psobsz'
   'undefine psobsz'
else
   'set dfile 'uofile
   'setlons'
   'set lat -90 90'
   'getdates'
   'getinfo   tmin'
              tmin = result
   'getinfo   tmax'
              tmax = result
   'set t    'tmin
   'getinfo   date'
              begdateo = result
   'set t    'tmax
   'getinfo   date'
              enddateo = result
   'set time 'begdateo' 'enddateo
   'set z 1'
   'run   setenv "ZDIM" 1'
   'define psobs = 100000 + lon-lon + lat-lat'
   'makezf psobs psobs z'
   'run unsetenv "ZDIM"'
   'set lon 0'
   'seasonal psobsz'
   'undefine psobsz'
endif

   'run setenv   "BEGDATEO" 'begdateo
   'run setenv   "ENDDATEO" 'enddateo

* Compute Quadratics of Remaining Fields
* --------------------------------------
   'set dfile 'uofile
   'setlons'
   'set lat -90 90'
   'getdates'
   'getinfo   tmin'
              tmin = result
   'getinfo   tmax'
              tmax = result
   'set t    'tmin
   'getinfo   date'
              begdateo = result
   'set t    'tmax
   'getinfo   date'
              enddateo = result
   'set time 'begdateo' 'enddateo
   'setz'

   'makequad -q1 'uoname' -file1 'uofile' -scale1 'uoscale' -q2 'uoname' -file2 'uofile' -scale2 'uoscale' -qq 'uuoname' -tag o -zonal'
   'makequad -q1 'voname' -file1 'vofile' -scale1 'voscale' -q2 'voname' -file2 'vofile' -scale2 'voscale' -qq 'vvoname' -tag o -zonal'
   'makequad -q1 'toname' -file1 'tofile' -scale1 'toscale' -q2 'toname' -file2 'tofile' -scale2 'toscale' -qq 'ttoname' -tag o -zonal'
   'makequad -q1 'honame' -file1 'hofile' -scale1 'hoscale' -q2 'honame' -file2 'hofile' -scale2 'hoscale' -qq 'hhoname' -tag o -zonal'
   'makequad -q1 'woname' -file1 'wofile' -scale1 'woscale' -q2 'woname' -file2 'wofile' -scale2 'woscale' -qq 'wwoname' -tag o -zonal'
   'makequad -q1 'qoname' -file1 'qofile' -scale1 'qoscale' -q2 'qoname' -file2 'qofile' -scale2 'qoscale' -qq 'qqoname' -tag o -zonal'

   'makequad -q1 'uoname' -file1 'uofile' -scale1 'uoscale' -q2 'voname' -file2 'vofile' -scale2 'voscale' -qq 'uvoname' -tag o -zonal'
   'makequad -q1 'voname' -file1 'vofile' -scale1 'voscale' -q2 'toname' -file2 'tofile' -scale2 'toscale' -qq 'vtoname' -tag o -zonal'
   'makequad -q1 'voname' -file1 'vofile' -scale1 'voscale' -q2 'qoname' -file2 'qofile' -scale2 'qoscale' -qq 'vqoname' -tag o -zonal'
   'makequad -q1 'woname' -file1 'wofile' -scale1 'woscale' -q2 'toname' -file2 'tofile' -scale2 'toscale' -qq 'wtoname' -tag o -zonal'
   'makequad -q1 'woname' -file1 'wofile' -scale1 'woscale' -q2 'qoname' -file2 'qofile' -scale2 'qoscale' -qq 'wqoname' -tag o -zonal'

*******************************************************
****         Loop over Seasons to Process          ****
*******************************************************

ptop.1 = 100
ptop.2 = 1

       m = 1
while( m > 0 )
    season = subwrd(seasons,m)
if( season = '' )
         m = -1
else
         m = m+1
         say 'Processing Season: 'season

'set dfile 'umfile
'count "'season'" 'begdate' 'enddate
 nmod =  result

'set dfile 'uofile
'count "'season'"'
 nobs =  result

* Make Variance Plots
* -------------------
          VARSM = umname' 'vmname' 'tmname' 'qmname' 'wmname' 'hmname
          VARSO = uoname' 'voname' 'toname' 'qoname' 'woname' 'honame
'numargs 'VARSM
          nvars = result

n = 1
while( n<=nvars )
       varmname = subwrd( VARSM,n )
       varoname = subwrd( VARSO,n )

if( varmname != 'NULL' & varoname != 'NULL' )

'uppercase 'varmname
            varname = 'Var_'result

'set dfile 'umfile
'setz'
'define qsqsmz = sqrt(Var'varmname'sm'season')'
'define qpqpmz = sqrt(Var'varmname'pm'season')'
'set dfile 'uofile
'setz'
'define qsqsoz = sqrt(Var'varoname'so'season')'
'define qpqpoz = sqrt(Var'varoname'po'season')'

k = 1
while( k<=2 )
ptop = ptop.k

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/zplot -EXPID 'EXPID' -EXPORT 'varname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -MDESC 'mdesc' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

k = k+1
endwhile

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/lplot -EXPID 'EXPID' -EXPORT 'varname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

endif

n = n+1
endwhile

* Make KE Plots
* -------------
'set dfile 'umfile
'setz'
'define          KEsm'season' = 0.5*( Var'umname'sm'season' + Var'vmname'sm'season' )'
'define          KEpm'season' = 0.5*( Var'umname'pm'season' + Var'vmname'pm'season' )'
'define qsqsmz = KEsm'season
'define qpqpmz = KEpm'season
'set dfile 'uofile
'setz'
'define          KEso'season' = 0.5*( Var'uoname'so'season' + Var'voname'so'season' )'
'define          KEpo'season' = 0.5*( Var'uoname'po'season' + Var'voname'po'season' )'
'define qsqsoz = KEso'season
'define qpqpoz = KEpo'season

k = 1
while( k<=2 )
ptop = ptop.k

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/zplot -EXPID 'EXPID' -EXPORT Var_KE -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -MDESC 'mdesc' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

k = k+1
endwhile

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/lplot -EXPID 'EXPID' -EXPORT Var_KE -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile


* Make CoVariance Plots
* ---------------------
                  nvars = 1
          COVARSM.nvars = umname' 'vmname ; nvars = nvars + 1
          COVARSM.nvars = vmname' 'tmname ; nvars = nvars + 1
          COVARSM.nvars = vmname' 'qmname ; nvars = nvars + 1
          COVARSM.nvars = wmname' 'tmname ; nvars = nvars + 1
          COVARSM.nvars = wmname' 'qmname ; nvars = nvars + 1
                  nvars = nvars - 1

                  nvars = 1
          COVARSO.nvars = uoname' 'voname ; nvars = nvars + 1
          COVARSO.nvars = voname' 'toname ; nvars = nvars + 1
          COVARSO.nvars = voname' 'qoname ; nvars = nvars + 1
          COVARSO.nvars = woname' 'toname ; nvars = nvars + 1
          COVARSO.nvars = woname' 'qoname ; nvars = nvars + 1
                  nvars = nvars - 1

n = 1
while( n<=nvars )
       varmname1 = subwrd( COVARSM.n,1 )
       varmname2 = subwrd( COVARSM.n,2 )
       varoname1 = subwrd( COVARSO.n,1 )
       varoname2 = subwrd( COVARSO.n,2 )

'uppercase 'varmname1
             varname1 = result
'uppercase 'varmname2
             varname2 = result
               cname  = 'Cov_'varname1'_'varname2

'set dfile 'umfile
'setz'
'define qsqsmz = Cov'varmname1''varmname2'sm'season
'define qpqpmz = Cov'varmname1''varmname2'pm'season
'set dfile 'uofile
'setz'
'define qsqsoz = Cov'varoname1''varoname2'so'season
'define qpqpoz = Cov'varoname1''varoname2'po'season

k = 1
while( k<=2 )
ptop = ptop.k

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/zplot -EXPID 'EXPID' -EXPORT 'cname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -MDESC 'mdesc' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else
            flag = "next"
        endif
       'c'
   endwhile

k = k+1
endwhile

           flag = ""
   while ( flag = "" )
  'run 'geosutil'/plots/quad/lplot -EXPID 'EXPID' -EXPORT 'cname' -SEASON 'season' -NMOD 'nmod' -NOBS 'nobs' -OUTPUT 'OUTPUT' -PFILE 'umfile' -OFILE 'uofile' -OBDATE 'begdateo' -OEDATE 'enddateo' -ODESC 'odesc' -ONAME 'obsid' -PTOP 'ptop
        if( DEBUG = "debug" )
            say "Hit  ENTER  to repeat plot or non-blank to continue"
            pull flag
        else 
            flag = "next"
        endif
       'c'
   endwhile

n = n+1
endwhile

endif
endwhile

num = num + 1
endwhile

return

function setname (name1,name2,cfile1,tag)

    flag = false

        if( name1 != 'NULL' & name2 != 'NULL' )
            if( name1 != name2 )
                cname  = 'Cov_'name1'_'name2
               'lowercase 'cname
                          lcname = result
               'set dfile 'cfile1
               'query file'
                numvar = sublin(result,6)
                numvar = subwrd(numvar,5)
                m = 1
                while ( m<numvar+1 )
                field = sublin(result,6+m)
                field = subwrd(field,1)
                      if( lcname=field )
                      flag = true
                      endif
                m = m+1
                endwhile
            else
                cname  = 'Var_'name1
               'lowercase 'cname
                          lcname = result
               'set dfile 'cfile1
               'query file'
                numvar = sublin(result,6)
                numvar = subwrd(numvar,5)
                m = 1
                while ( m<numvar+1 )
                field = sublin(result,6+m)
                field = subwrd(field,1)
                      if( lcname=field )
                      flag = true
                      endif
                m = m+1
                endwhile
            endif
        endif

    if (flag=false)
        cname = 'NULL'
    endif

return cname