#!/bin/bash
# L. Fita, June 2016. Generic script to plot model outputs

main='model_graphics.bash'

function uploadvars() {
# Function to upload variables to the system from an ASCII file as:
#   [varname] = [value]
  fileval=$1
  errormsg='ERROR -- error -- ERROR -- error'

  if test ! -f ${fileval}; then
    echo ${errormsg}
    echo "  "${fname}": file '"${fileval}"' does not exist!!"
    exit
  fi

  Nlin=`wc -l ${fileval} | awk '{print $1}'`

  ilin=1
  while test ${ilin} -le ${Nlin}; do
    line=`head -n ${ilin} ${fileval} | tail -n 1`
    varname=`echo ${line} | tr '=' ' ' | awk '{print $1}'`
    value=`echo ${line} | tr '=' ' ' | awk '{print $2}'`
    Lvarname=`expr length ${varname}'0'`

    if test ${Lvarname} -gt 1 && test ! ${varname:0:1} = '#'; then
      export ${varname}=${value}
    fi
    ilin=`expr ${ilin} + 1`
  done
}

function isInlist() {
# Function to check whether a value is in a list
  list=$1
  value=$2
  
  is=`echo ${list} | tr ':' '\n' | awk '{print "@"$1"@"}' | grep '@'${value}'@' | wc -w`
  if test ${is} -eq 1
  then
    true
  else
    false
  fi
}

function ferrmsg() {
# Function to exit and write an error message
#  comdexit: code number exit from the last execution
#  ref: script/reference from which error occurs
#  msg: message to write ('!' for spaces, '#' for end of line)
  comdexit=$1
  ref=$2
  msg=$3

  if test ${comdexit} -ne 0; then
    echo "ERROR -- error -- ERROR -- error"
    echo "  "${ref}": "$(echo ${msg} | tr '!' ' ' | tr '#' '\n')" !!"
    exit
  fi
}

function ferrmsgF() {
# Function to exit and write an error message and remove a file
#  comdexit: code number exit from the last execution
#  ref: script/reference from which error occurs
#  msg: message to write ('!' for spaces, '#' for end of line)
#  FileName: name of the file to remove in case of error
  comdexit=$1
  ref=$2
  msg=$3
  FileName=$4

  if test ${comdexit} -ne 0; then
    echo "ERROR -- error -- ERROR -- error"
    echo "  "${ref}": "$(echo ${msg} | tr '!' ' ' | tr '#' '\n')" !!"
    rm ${FileName} >& /dev/null
    exit
  fi
}

function WRF_toCF() {
# Function to pass a WRF original file to CF-conventions
#   wrff= WRF file
#   wrfvl= WRF longitude var name (it might come from WPS!)
#   wrfvL= WRF latitude var name (it might come from WPS!)

  wrff=$1
  wrfvl=$2
  wrfvL=$3

  fdims=`${pyHOME}/nc_var.py -o idims -f ${wrff} | grep alldims | tr '=' ' ' |       \
    awk '{print $3}'`

  CFdims='west_east@lon:south_north@lat:Time@time'
  CFds=`echo ${CFdims} | tr ':' ' '`
  
  pyout=`${pyHOME}/nc_var.py -o WRF_CFtime_creation -S 19491201000000,minutes        \
    -f ${wrff} -v time`
  pyn=$?
  Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
  ferrmsgF ${pyn} ${fname} "python!'WRF_CFtime_creation'!failed"${Spyout} ${wrff}
  
  pyout=`${pyHOME}/nc_var.py -o WRF_CFlonlat_creation -v ${wrfvl},${wrfvL}           \
    -f ${wrff} -S lon,lat`
  pyn=$?
  Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
  ferrmsgF ${pyn} ${fname} "python!'WRF_CFlonlat_creation'!failed"${Spyout} ${wrff}

  for CFd in ${CFds}; do
    cd=`echo ${CFd} | tr '@' ' ' | awk '{print $2}'`
    wd=`echo ${CFd} | tr '@' ' ' | awk '{print $1}'`
    if $(isInlist ${fdims} ${wd}); then
      pyout=`${pyHOME}/nc_var.py -o chdimname -f ${wrff} -S ${wd}':'${cd}`
      pyn=$?
      Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
      ferrmsgF ${pyn} ${fname} "python!'chdimname'!'"${wd}"'!failed#"${Spyout} ${wrff}
    fi
  done

}

function compute_variable() {
# Function to compute a variable
#   idir: directory with the input files
#   usefiles: ',' list of files as [headerf]@[file1]|...|[fileN]
#   odir: directory to write the output files
#   cvar: [CFvarn]|[varkind]|[headerf]|[varmod]|[vardiag]
#     [CFvarn]: CF name of the variable
#     [vark]: kind of variable: 
#        acc: temporal accumulated values
#        diff: differences between models
#        direct: no statistics
#        Lsec: latitudinal section (latitudinal value must be given, [var]@[lat]) 
#        lsec: longitudinal section (longitudinal value must be given, [var]@[lat])
#        lmean: longitudinal mean values
#        tmean: temporal mean values
#        zsum: vertical aggregated values
#     [headerf]: header of the files to use
#     [modvar]: variable to use from the file
#     [diagvar]: variable computed from the diagnostics (diagnostic.py)
#   mdim: name of the dimensions in the model
#   mvdim: name of the vaariable-dimensions in the model
#   scratch: whether is has to be done from the scratch or not
  fname='compute_variable'
  idir=$1
  usefiles=$2
  odir=$3
  cvar=$4
  mdim=$5
  mvdim=$6
  scratch=$7

  CFvarn=`echo ${cvar} | tr '|' ' ' | awk '{print $1}'`
  vark=`echo ${cvar} | tr '|' ' ' | awk '{print $2}'`
  headerf=`echo ${cvar} | tr '|' ' ' | awk '{print $3}'`
  modvar=`echo ${cvar} | tr '|' ' ' | awk '{print $4}'`
  diagvar=`echo ${cvar} | tr '|' ' ' | awk '{print $5}'`

  cfiles=`echo ${usefiles} | tr ',' '\n' | grep ${headerf} | tr '|' ' ' | \
    awk '{print $2}' | tr '@' ' '`

# dimensions
  dnx=`echo ${mdim} | tr ',' ' ' | awk '{print $1}'`
  dny=`echo ${mdim} | tr ',' ' ' | awk '{print $2}'`
# var-dimensions
  vdnx=`echo ${mvdim} | tr ',' ' ' | awk '{print $1}'`
  vdny=`echo ${mvdim} | tr ',' ' ' | awk '{print $2}'`

  cd ${odir}

  # Computing separately and then joinging for all files
  Ntotfiles=`echo ${cfiles} | wc -w | awk '{print $1}'`
  ifile=1

  for cf in ${cfiles}; do
    ifS=`printf "%0${Ntotfiles}d" ${ifile}`

    # Computing variable
    filen=${odir}/${CFvarn}_${headerf}_${ifile}-${Ntotfiles}.nc
    if ${scratch}; then 
      rm ${filen}
      rm ${odir}/${CFvarn}_${headerf}.nc
    fi

    if test ! -f ${filen} && test ! -f ${odir}/${CFvarn}_${headerf}.nc; then
      if test ! ${modvar} = 'None'; then
        # model variable
        values=${modvar}',0,-1,-1'
        vs=${modvar},${vdnx},${vdny},${vdnz},${vdnt}
        pyout=`${pyHOME}/nc_var.py -f ${cf} -o DataSetSection_multivars -v ${vs}     \
          -S ${values} | grep succesfull | awk '{print $6}' | tr '"' ' '`
        pyn=$?
        Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
        ferrmsg ${pyn} ${fname} "python!'DataSetSection_multivars'!failed"${Spyout}
        mv ${pyout} ${filen}

        pyout=`${pyHOME}/nc_var.py -f ${filen} -o chvarname -v ${modvar} -S ${CFvarn}`
        pyn=$?
        Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
        ferrmsgF ${pyn} ${fname} "python!'chvarname'!failed"${Spyout} ${filen}
      else
        # diagnostic variable
        dims=${dnt}@${vdnt},${dnz}@${vdnz},${dny}@${vdny},${dnx}@${vdnx}
        diagn=`echo ${diagvar} | tr ':' '\n' | head -n 1`
        Ndiagvars=`echo ${diagvar} | tr ':' ' ' | wc -w | awk '{print $1}'`
        idiagv=2
        diagc=''
        while test ${idiagv} -le ${Ndiagvars}; do
          diag1=`echo ${diagvar} | tr ':' '\n' | head -n ${idiagv} | tail -n 1`
          if test ${idiagv} -eq 2; then
            diagc=${diag1}
          else
            diagc=${diagc}'@'${diag1}
          fi
          idiagv=`expr ${idiagv} + 1`
        done

        pyout=`python ${pyHOME}/diagnostics.py -f ${cf} -d ${dims} -v ${diagn}'|'${diagc}`
        pyn=$?
        Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
        ferrmsg ${pyn} ${fname} "python!'diagnostics'!failed#"${Spyout}
        mv diagnostics.nc ${filen}
      fi

      # adding CF lon,lat,time in WRF files
      if test ${headerf:0:3} = 'wrf'; then
        WRF_toCF ${filen} ${vdnx} ${vdny}
      fi

LLUIS TO ADD ALL CF attributes!!!!!!!

    fi

    ifile=`expr ${ifile} + 1`
# End of files
  done

  # Joining variable files
  filen=${odir}/${CFvarn}_${headerf}.nc
  if ${scratch}; then rm ${filen}; fi

  if test ! -f ${filen}; then
    pyout=`python ${pyHOME}/nc_var.py -f ${CFvarn}_${headerf}_                       \
      -o netcdf_fold_concatenation -S ./,time -v all`
    pyn=$?
    Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
    ferrmsg ${pyn} ${fname} "python!'netcdf_fold_concatenation'!failed#"${Spyout}

    mv netcdf_fold_concatenated.nc ${filen}
    if test -f ${filen}; then
      rm ${CFvarn}_${headerf}_*
    fi
  fi
}

function compute_statistics(){
# Function to compute different statistics
#   idir: directory with the input files
#   usefiles: ',' list of files to use [file1],...,[fileN]
#   odir: directory to write the output files
#   cvar: [CFvarn]|[varkind]|[headerf]|[varmod]|[vardiag]
#     [CFvarn]: CF name of the variable
#     [vark]: kind of variable: 
#        acc: temporal accumulated values
#        diff: differences between models
#        direct: no statistics
#        Lmean: latitudinal mean values
#        Lsec: latitudinal section (latitudinal value must be given, [var]@[lat]) 
#        lmean: longitudinal mean values
#        lsec: longitudinal section (longitudinal value must be given, [var]@[lat])
#        tmean: temporal mean values
#        xmean: x-axis mean values
#        ymean: y-axis mean values
#        zsum: vertical aggregated values
#     [headerf]: header of the files to use
#     [modvar]: variable to use from the file
#     [diagvar]: variable computed from the diagnostics (diagnostic.py)
#   mdim: name of the dimensions in the model
#   mvdim: name of the vaariable-dimensions in the model
#   scratch: whether is has to be done from the scratch or not
  fname='compute_statistics'

  idir=$1
  usefiles=$2
  odir=$3
  cvar=$4
  mdim=$5
  mvdim=$6
  scratch=$7

  CFvarn=`echo ${cvar} | tr '|' ' ' | awk '{print $1}'`
  vark=`echo ${cvar} | tr '|' ' ' | awk '{print $2}'`
  headerf=`echo ${cvar} | tr '|' ' ' | awk '{print $3}'`
  modvar=`echo ${cvar} | tr '|' ' ' | awk '{print $4}'`
  diagvar=`echo ${cvar} | tr '|' ' ' | awk '{print $5}'`

  cfiles=`echo ${usefiles} | tr ',' ' '`

# dimensions
  dnx=`echo ${mdim} | tr ',' ' ' | awk '{print $1}'`
  dny=`echo ${mdim} | tr ',' ' ' | awk '{print $2}'`
# var-dimensions
  vdnx=`echo ${mvdim} | tr ',' ' ' | awk '{print $1}'`
  vdny=`echo ${mvdim} | tr ',' ' ' | awk '{print $2}'`

  cd ${odir}

  filen=${CFvarn}_${headerf}_${vark}.nc
  if ${scratch}; then rm ${filen} >& /dev/null; fi

  if test ! -f ${filen}; then
    # Computing stats
    case ${vark} in
       # temporal accumulated values
      'acc') 
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
      # differences between models
      'diff')
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
      # no statistics
      'direct')
        cp ${cfiles} ${filen}
      ;;
      # latitudinal mean values
      'Lmean')
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
      # latitudinal section (latitudinal value must be given, [var]@[lat]) 
      'Lsec')
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
      # longitudinal section (longitudinal value must be given, [var]@[lat])
      'lsec')
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
      # longitudinal mean values
      'lmean')
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
      # temporal mean values
      'tmean')
        vals='time:-1,time,mean,lon:lat:'${vdnz}':time'
        dims='time@time,'${dnz}'@'${vdnz}',lat@lat,lon@lon'

        pyout=`python ${pyHOME}/nc_var.py -o file_oper_alongdims -S ${vals}          \
          -f ${cfiles} -v ${CFvarn}`
        pyn=$?
        Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
        ferrmsg ${pyn} ${fname} "python!'file_oper_alongdims'!'tmean'!failed#"${Spyout}
        mv file_oper_alongdims_mean.nc ${filen}
      ;;
      # x-axis mean values
      'xmean')
        vals='lon:-1,lon,mean,lon:lat:'${vdnz}':time'
        dims='time@time,'${dnz}'@'${vdnz}',lat@lat,lon@lon'
        echo "python ${pyHOME}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${cfiles} -v ${CFvarn}"

        pyout=`python ${pyHOME}/nc_var.py -o file_oper_alongdims -S ${vals}         \
          -f ${cfiles} -v ${CFvarn}`
        pyn=$?
        Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
        ferrmsg ${pyn} ${fname} "python!'file_oper_alongdims'!'tmean'!failed#"${Spyout}
        mv file_oper_alongdims_mean.nc ${filen}
      ;;
      # y-axis mean values
      'ymean')
        vals='lat:-1,lat,mean,lon:lat:'${vdnz}':time'
        dims='time@time,'${dnz}'@'${vdnz}',lat@lat,lon@lon'
        echo "python ${pyHOME}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${cfiles} -v ${CFvarn}"

        pyout=`python ${pyHOME}/nc_var.py -o file_oper_alongdims -S ${vals}         \
          -f ${cfiles} -v ${CFvarn}`
        pyn=$?
        Spyout=`echo ${pyout} | tr '\n' '#' | tr ' ' '!'`
        ferrmsg ${pyn} ${fname} "python!'file_oper_alongdims'!'tmean'!failed#"${Spyout}
        mv file_oper_alongdims_mean.nc ${filen}
      ;;
      # vertical aggregated values
      'zsum')
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
      *)
        echo ${errmsg}
        echo "  "${fname}": kind '"${vark}"' not ready !!"
        exit
      ;;
    esac
  fi
#  exit
# LLUIS
}

#######    #######
## MAIN
    #######
rootsh=`pwd`

# Uploading environment
uploadvars model_graphics.dat

if test scratch = 'true'; then
  scratch=true
  echo ${warnmsg}
  echo "  "${main}": starting from the SCRATCH !!"
else
  scratch=false
fi

timeval='tstep'
zval='null'
dimval='lon,lat,pressure,time'

compute=0
computetlon=false
plot=false
single2Dplot=false
lonZsec=false
differences=true

combosfile=${HOMEpy}/diagnostics.inf

####### ###### ##### #### ### ## #

mods=`echo ${models} | tr ':' ' '`
varks=`echo ${varkinds} | tr ':' ' '`

# Models loop
##
for mod in ${mods}; do
  case ${mod} in
    'WRF')
      exps=`echo ${WRFexps} | tr ':' ' '`
      fheaders=`echo ${WRFheaders} | tr ':' ' '`
    ;;
    'LMDZ')
      exps=`echo ${LMDZexps} | tr ':' ' '`
      fheaders=`echo ${LMDZheaders} | tr ':' ' '`
    ;;
    'WRF_LMDZ')
      exps=`echo ${WRF_LMDZexps} | tr ':' ' '`
      fheaders=`echo ${WRF_LMDZheaders} | tr ':' ' '`
    ;;
    '*')
      echo ${errmsg}
      echo "  "${main}": model '"${mod}"' not ready!!"
      exit
    ;;
  esac  

  modinf=`$pyHOME/nc_var.py -o model_characteristics -f None -S ${mod} | \
    grep singleline | awk '{print $2}'`
  ferrmsg $? $main "python!'model_characteristics'!failed"
  dnx=`echo ${modinf} | tr ';' '\n' | grep dimxn | tr '=' ' ' | awk '{print $2}'`
  dny=`echo ${modinf} | tr ';' '\n' | grep dimyn | tr '=' ' ' | awk '{print $2}'`
  dnz=`echo ${modinf} | tr ';' '\n' | grep dimzn | tr '=' ' ' | awk '{print $2}'`
  dnt=`echo ${modinf} | tr ';' '\n' | grep dimtn | tr '=' ' ' | awk '{print $2}'`
  vdnx=`echo ${modinf} | tr ';' '\n' | grep vardxn | tr '=' ' ' | awk '{print $2}'`
  vdny=`echo ${modinf} | tr ';' '\n' | grep vardyn | tr '=' ' ' | awk '{print $2}'`
  vdnz=`echo ${modinf} | tr ';' '\n' | grep vardzn | tr '=' ' ' | awk '{print $2}'`
  vdnt=`echo ${modinf} | tr ';' '\n' | grep vardtn | tr '=' ' ' | awk '{print $2}'`
  echo ${mod}
  echo "  dims: "${dnx}", "${dny}", "${dnz}", "${dnt}
  echo "  var dims: "${vdnx}", "${vdny}", "${vdnz}", "${vdnt}
  moddims=${dnx}','${dny}','${dnz}','${dnt}
  modvdims=${vdnx}','${vdny}','${vdnz}','${vdnt}

# Experiments loop
##
  for exp in ${exps}; do
    echo "  "${exp}"..."
    iwdir=${ifold}/${mod}/${exp}

    # Does input folder exist?
    if test ! -d ${iwdir}; then
      echo ${errmsg}
      echo "  "${main}": folder '"${iwdir}"' does not exist !!"
      exit
    fi

    owdir=${ofold}/${mod}/${exp}
    mkdir -p ${owdir}

    # Need to pass to analyze all the data?
    if ${scratch}; then rm ${owdir}/varcompute.inf >& /dev/null; fi
    if test ! -f ${owdir}/varcompute.inf; then

      # Does input folder has header files?
      cd ${iwdir}
      files=''
      testfiles=''
      ih=1
      for fh in ${fheaders}; do
        files1h=`ls -1 ${fh}* | tr '\n' '@'`
        Lfiles1h=`expr length ${files1h}'0'`
        if test ${Lfiles1h} -lt 2; then
          ferrmsg 1 $main "folder!:!"${iwdir}"!does!not!contain!files!"${fh}"*"
        fi
        testfiles1h=`ls -1 ${fh}* | head -n 1`
        if test ${ih} -eq 1; then
          files=${fh}'|'${files1h}
          testfiles=${fh}'|'${testfiles1h}
        else
          files=${files}','${fh}'|'${files1h}
          testfiles=${testfiles}'@'${fh}'|'${testfiles1h}
        fi
        ih=`expr ${ih} + 1`
      done
      testfs=`echo ${testfiles} | tr '@' ' '`
      cd ${rootsh}

# Kind variables loop
##
      ik=1
      itotv=1
      for vark in ${varks}; do
        echo "    "${vark}" ..."
        case ${vark} in
          'acc')
            varvks=${varacc}
          ;;
          'diff')
            varvks=${vardiff}
          ;;
          'direct')
            varvks=${vardirect}
          ;;
          'Lmean')
            varvks=${varLmean}
          ;;
          'Lsec')
            varvks=${varLsec}
          ;;
          'lmean')
            varvks=${varlmean}
          ;;
          'lsec')
            varvks=${varlsec}
          ;;
          'tmean')
            varvks=${vartmean}
          ;;
          'xmean')
            varvks=${varxmean}
          ;;
          'ymean')
            varvks=${varymean}
          ;;
          'zsum')
            varvks=${varzsum}
          ;;
          '*')
            echo ${errmsg}
            echo "  "${main}": variable kind '"${vark}"' not ready!!"
            exit
          ;;
        esac

# Do we have variables for this kind?
        Lvarvks=`expr length ${varvks}'0'`
        if test ${Lvarvks} -lt 2; then
          ferrmsg 1 ${main} "variable!kind!"${vark}"!without!variables"
        fi
        if test ${ik} -eq 1; then
          allvars=${varvks}
        else
          allvars=${allvars}':'${varvks}
        fi
        ik=`expr ${ik} + 1`
        vars=`echo ${varvks} | tr ':' ' '`

# Variables loop
##
        iv=1
        for var in ${vars}; do
          echo "      "${var}
          # How to copmute it?
          cancompute=true
          for ftest in ${testfs}; do
            headerf=`echo ${ftest} | tr '|' ' ' | awk '{print $1}'`
            filen=`echo ${ftest} | tr '|' ' ' | awk '{print $2}'`
            compute=`${pyHOME}/nc_var.py -o computevar_model -f ${iwdir}/${filen}      \
              -S ${var}`
            ferrmsg $? ${main} "python!'computevar_model'!failed!#"$(echo ${compute} | \
               tr ' ' '!')
            varmod=`echo ${compute} | tr ' ' '#' | tr '|' ' ' | awk '{print $2}' |   \
               tr '@' ' ' | awk '{print $2}' | tr '=' ' ' | awk '{print $2}'`
            vardiag=`echo ${compute} | tr ' ' '#' | tr '|' ' ' | awk '{print $2}' |  \
               tr '@' ' ' | awk '{print $3}' | tr '=' ' ' | awk '{print $2}'`
            echo "  "${var}" mod:"${varmod}" diag: "${vardiag}
            if test ${varmod} = 'None' && test ${vardiag} = 'None'; then 
              cancompute=false
            else
              cancompute=true
              # Should be considered that variable can also be computed by both ways?
              break
            fi
          done
          if ! ${cancompute}; then
            msg="there!is!no!way!to!compute!'"${var}"'!for!model!"${mod}
# Too extrict!
#            ferrmsg 1 ${main} ${msg}
            echo ${warnmsg}
            echo $(echo $msg | tr '!' ' ')
          fi

          # A ';' list 'varcompute' it is created for each variable giving:
          #   [var]|[vark]|[headerf][varmod]|[vardiag]
          # This list will be used to compute a new file for each variable
          varcomp=${var}'|'${vark}'|'${headerf}'|'${varmod}'|'${vardiag}
          if test ${itotv} -eq 1; then
            varcompute=${varcomp}
          else
            varcompute=${varcompute}';'${varcomp}
          fi

          iv=`expr ${iv} + 1`
          itotv=`expr ${itotv} + 1`

# end loop of variables
        done
# end of kind of variables
      done

      # Outwritting the varcompute to avoid next time (if it is not scratch!)
      cat << EOF > ${owdir}/varcompute.inf
files: ${files}
varcompute: ${varcompute}
itotv: ${itotv}
EOF
    else
      echo $warnmsg
      echo "  "${main}": getting already data information from the experiment!"
      files=`cat ${owdir}/varcompute.inf | grep files | awk '{print $2}'`
      varcompute=`cat ${owdir}/varcompute.inf | grep varcompute | awk '{print $2}'`
      itotv=`cat ${owdir}/varcompute.inf | grep itotv | awk '{print $2}'`
# End of avoiding to repeat all the experiment search
    fi

    echo "  For experiment '"${exp}"' is required to compute: "${itotv}" variables"
# Computing files for each variable
## 
    echo "      Computing all variables ..."
    cd $owdir
    ic=1
    isc=1
    cvars=`echo ${varcompute} | tr ';' ' '`
    for cvar in ${cvars}; do
      CFv=`echo ${cvar} | tr '|' ' ' | awk '{print $1}'`
      vark=`echo ${cvar} | tr '|' ' ' | awk '{print $2}'`
      fileh=`echo ${cvar} | tr '|' ' ' | awk '{print $3}'`
      modv=`echo ${cvar} | tr '|' ' ' | awk '{print $4}'`
      diagv=`echo ${cvar} | tr '|' ' ' | awk '{print $5}'`
      echo "        "${CFv}"; "${modv}" "${diagv}

      if test ! ${modv} = 'None' || test ! ${diagv} = 'None'; then
        compute_variable ${iwdir} ${files} ${owdir} ${cvar} ${moddims} ${modvdims}   \
          ${scratch}
        ic=`expr ${ic} + 1`

        if test ! ${vark} = 'diff'; then
          compute_statistics ${iwdir} ${CFv}_${fileh}.nc ${owdir} ${cvar}            \
            ${moddims} ${modvdims} ${scratch}
          isc=`expr ${isc} + 1`
        else
          echo "  "${main}": differences will be calculated when all the " \
            "model/experiments will be done !"
        fi
      else
        echo "        not computing '"${CFv}"' for model '"${mod}"' !!"
      fi

      # exit
# end of computing vars
    done
    echo "  "${main}": "${ic}" variables has been computed"
    echo "  "${main}": "${isc}" statistics has been computed"

    cd ${rootsh}
    exit
# end of experiments
  done
# end of models
done

exit

istep=0

if test ${compute} -eq 1; then
  for exp in ${exps}; do
    echo "exp: "${exp}
    if test ${exp} == 'AR40'; then
      var2Dtmeans=${var2Dtmeans_AR40}
    else
      var2Dtmeans=${var2Dtmeans_NPv31}
    fi
    for vark in ${varks}; do
      echo "  vark: "${vark}

      if test ${vark} = 'mean'; then
        fileorig=${ifoldmean}/${exp}/${filenmean}
        if test ! -f ${fileorig}; then
          fileworig=${filensfc}
          dvals='T:Time,Z:bottom_top,Y:south_north,X:west_east'
          dimnames='T:Time,Z:bottom_top,Y:south_north,X:west_east'
          python ${HOMEpy}/vertical_interpolation.py                            \
            -f ${ofold}/${exp}/${filensfc} -o WRFp -i ${plevels} -k 'lin'            \
            -v WRFt,U,V,WRFrh,WRFght -d ${dimnames} -D T:Times,Z:ZNU,Y:XLAT,X:XLONG
          if test $? -ne 0; then
            echo ${errormsg}
            echo "  python failed!"
            echo "  python ${HOMEpy}/vertical_interpolation.py"                            \
            "-f ${ofold}/${exp}/${filensfc} -o WRFp -i ${plevels} -k 'lin'"            \
            "-v WRFt,U,V,WRFrh,WRFght -d ${dimnames} -D T:Times,Z:ZNU,Y:XLAT,X:XLONG"
            rm ${ofold}/${exp}/vertical_interpolation_WRFp.nc >& /dev/null
            exit
          fi
          mv vertical_interpolation_WRFp.nc ${ofold}/${exp}
        fi
        vars=`echo ${mean_variables} | tr ':' ' '`
        ofile='mean_variables'
      elif test ${vark} = 'sfc'; then
        fileorig=${ifoldsfc}/${exp}/${filensfc}
        if test ! -f ${fileorig}; then
          fileworig=${filensfc}
          cp ${ifoldsfc}/${exp}/${fileworig} ${ifoldsfc}/${exp}/${filensfc}
          python ${HOMEpy}/nc_var.py -o WRF_CFxtime_creation -f ${ifoldsfc}/${exp}/${filensfc} \
            -S 19491201000000,hours -v time
          python ${HOMEpy}/nc_var.py -o WRF_CFlonlat_creation -f ${ifoldsfc}/${exp}/${filensfc}\
            -S longitude,latitude -v XLONG,XLAT
        fi
        if test ${exp} = 'AR40'; then
          vars=`echo ${sfc_variables_AR40} | tr ':' ' '`
        else
          vars=`echo ${sfc_variables_NPv31} | tr ':' ' '`
        fi
        ofile='sfc_variables'
      elif test ${vark} = 'diag'; then
        fileorig=${ifolddiag}/${exp}/${filendiag}
        vars=`echo ${diag_variables} | tr ':' ' '`
        if test ! -f ${fileorig}; then
          values='Time@XTIME,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
          ivar=0
          varcombos=''
          for var in ${vars}; do
            echo "    var:"${var}
            varn=`echo ${var} | tr '_' ' ' | awk '{print $2}'`
            combo=`python ${HOMEpy}/diagnostics.py -f ${combosfile} -d variable_combo -v ${varn}\
              | grep COMBO | awk '{print $2}'`
            if test ${combo} = 'ERROR'; then
              echo ${errormsg}
              echo "  No variable '"${varn}"' in '"${combosfile}"' !!"
              exit
            fi
            if test ${ivar} -eq 0; then
              varcombos=${varn}'|'${combo}
            else
              varcombos=${varcombos}','${varn}'|'${combo}
            fi
            ivar=`expr ${ivar} + 1`
          done
          python ${HOMEpy}/diagnostics.py -d ${values} -f ${ifoldsfc}/${exp}/${filensfc} -v \
            ${varcombos}
          if test $? -ne 0; then
            echo ${errormsg}
            echo "  python failed!!"
            echo python ${HOMEpy}/diagnostics.py -d ${values} -f ${ifoldsfc}/${exp}/${filensfc}\
              -v ${varcombos}
            exit
          fi
          python ${HOMEpy}/nc_var.py -o WRF_CFxtime_creation -f diagnostics.nc                 \
            -S 19491201000000,hours -v time
           python ${HOMEpy}/nc_var.py -o WRF_CFlonlat_creation -f diagnostics.nc               \
            -S longitude,latitude -v XLONG,XLAT
          mv diagnostics.nc ${ifolddiag}/${exp}
        fi
        ofile='diag_variables'
      elif test ${vark} = 'z'; then
        fileorig=${ifoldmean}/${exp}/${filenmean}
        vars=`echo ${z_variables} | tr ':' ' '`
        ofile='z_variables'
      fi

# Averaging
##
      echo "  averaging..."
      ivar=0
      for varn in ${vars}; do
        file=${fileorig}
        echo "    var: "${varn}
        var=`echo ${varn} | tr '|' ' ' | awk '{print $2}'`
        varval=`python ${HOMEpy}/drawing.py -o variable_values -S ${var} | grep all_values | awk '{print $3}'`
        varname=`echo ${varval} | tr ',' ' ' | awk '{print $1}'`

        echo "      std name: "${varname}
        varhead=`echo ${varname} | tr '_' ' ' | awk '{print $1}'`
        vartail=`echo ${varname} | tr '_' ' ' | awk '{print $2}'`

        if test ${vark} = 'mean'; then
          vals='time:-1|z:-1|y:-1|x:-1,time:x,mean,pressure:lat'
          python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${var}
          pyexec=$?
          if test ${pyexec} -ne 0; then
            echo ${errormsg}
            echo "  "${main}": python fai1ls!"
            echo "python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${var}"
            exit
          else
            oper=`echo ${vals} | tr ',' ' ' | awk '{print $3}'`
            mv file_oper_alongdims_mean.nc ${ofold}/${exp}/${varname}${oper}.nc
          fi
        elif test ${vark} = 'z'; then
          vals='time:-1|z:-1|y:-1|x:-1,time:x,mean,pressure:lat'
          echo "NOT FINISHED!!!!"
          exit
        else
          vals='Time:-1|south_north:-1|west_east:-1,west_east,mean,time:XLAT'
          if test ${var} = 'ACRAINTOT'; then
            files='add|'${file}'|RAINC,add|'${file}'|RAINNC'
            python ${HOMEpy}/nc_var.py -S 'time|XLAT|XLONG@'${files} -o compute_opersvarsfiles \
              -v ${var}
            mv opersvarsfiles_${var}.nc ${ofold}/${exp}
            file=${ofold}/${exp}/opersvarsfiles_${var}.nc
          elif test ${var} = 'RAINTOT'; then
            dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
            python ${HOMEpy}/diagnostics.py -d ${dims} -v 'RAINTOT|RAINC@RAINNC@time' -f ${file}
            mv diagnostics.nc ${ofold}/${exp}/diagnostics_${varname}.nc
            file=${ofold}/${exp}/diagnostics_${varname}.nc
            var='pr'
          elif test ${var} = 'cllmh'; then
            dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
            clouds='cll:clm:clh'
            cls=`echo ${clouds} | tr ':' ' '`
            for cl in ${cls}; do
              file=${ofold}/${exp}/diagnostics.nc
              var=${cl}
              echo "    var: "${var}
              python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${cl}
              pyexec=$?
              if test ${pyexec} -ne 0; then
                echo ${errormsg}
                echo "  "${main}": python fails!"
                echo python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${cl}
                exit
              fi
              mv file_oper_alongdims_mean.nc ${ofold}/${exp}/${cl}${oper}.nc
              vals='Time:-1|south_north:-1|west_east:-1,Time,mean,XLONG:XLAT'
              dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
              python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${cl}
              pyexec=$?
              if test ${pyexec} -ne 0; then
                echo ${errormsg}
                echo "  "${main}": python fails!"
                exit
              else
                oper=`echo ${vals} | tr ',' ' ' | awk '{print $3}'`
                mv file_oper_alongdims_mean.nc ${ofold}/${exp}/${cl}t${oper}.nc
              fi
            done
#            break
          elif test ${var} = 'WRFbils' && test ! -f ${ofold}/${exp}/diagnostics_${varname}.nc; then
            dims='Time@time,south_north@XLAT,west_east@XLONG'
            python ${HOMEpy}/diagnostics.py -d ${dims} -v 'WRFbils|HFX@LH' -f ${file}
            mv diagnostics.nc ${ofold}/${exp}/diagnostics_${varname}.nc
            file=${ofold}/${exp}/diagnostics_${varname}.nc
            var='bils'
          elif test ${var} = 'WRFprw' && test ! -f ${ofold}/${exp}/diagnostics_${varname}.nc; then
            dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
            python ${HOMEpy}/diagnostics.py -d ${dims} -v 'WRFprw|WRFdens@QVAPOR' -f ${file}
            mv diagnostics.nc ${ofold}/${exp}/diagnostics_${varname}.nc
            file=${ofold}/${exp}/diagnostics_${varname}.nc
            var='prw'
          elif test ${var} = 'WRFrhs' && test ! -f ${ofold}/${exp}/diagnostics_${varname}.nc; then
            dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
            python ${HOMEpy}/diagnostics.py -d ${dims} -v 'WRFrhs|PSFC@T2@Q2' -f ${file}
            mv diagnostics.nc ${ofold}/${exp}/diagnostics_${varname}.nc
            file=${ofold}/${exp}/diagnostics_${varname}.nc
            var='huss'
          elif test ${var} = 'WRFprc'; then
            vals='Time:-1|south_north:-1|west_east:-1,west_east,mean,time:XLAT'
            dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
            var='prc'
          elif test ${var} = 'WRFprls'; then
            vals='Time:-1|south_north:-1|west_east:-1,west_east,mean,time:XLAT'
            dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
            var='prls'
          fi
          if test ! -f ${ofold}/${exp}/${varname}${oper}.nc; then
            python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${var}
            pyexec=$?
            if test ${pyexec} -ne 0; then
              echo ${errormsg}
              echo "  "${main}": python fails!"
              exit
            else
              mv file_oper_alongdims_mean.nc ${ofold}/${exp}/${varname}${oper}.nc
            fi
          fi
        fi
        echo "HereLluis"${var}
# Time means
        if $(isin_list ${var2Dtmeans} ${var}); then 
          echo "Computing time means!!"
          if test ${var} = 'cllmh'; then
            clouds='cll:clm:clh'
            clds=`echo ${clouds} | tr ':' ' '`
            for cld in ${clds}; do
              vals='Time:-1|south_north:-1|west_east:-1,Time,mean,XLONG:XLAT'
              dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
              python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${cld}
              pyexec=$?
              if test ${pyexec} -ne 0; then
                echo ${errormsg}
                echo "  "${main}": python fails!"
                exit
              else
                oper=`echo ${vals} | tr ',' ' ' | awk '{print $3}'`
                mv file_oper_alongdims_mean.nc ${ofold}/${exp}/${cld}t${oper}.nc
              fi
            done
          else
            vals='Time:-1|south_north:-1|west_east:-1,Time,mean,XLONG:XLAT'
            dims='Time@time,bottom_top@ZNU,south_north@XLAT,west_east@XLONG'
            python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${file} -v ${var}
            pyexec=$?
            if test ${pyexec} -ne 0; then
              echo ${errormsg}
              echo "  "${main}": python fails!"
              exit
            else
              oper=`echo ${vals} | tr ',' ' ' | awk '{print $3}'`
              mv file_oper_alongdims_mean.nc ${ofold}/${exp}/${varname}t${oper}.nc
            fi
          fi
        fi
        if test ${var} = 'cllmh'; then exit; fi
        ivar=`expr ${ivar} + 1`
#        exit
# end of vars
      done

# end of kind vars
    done

  done
# end of compute
fi
# Longitudinal means of tmeans
##
if ${computetlon}; then
  echo "Computing Longitudinal means of temporal means..."
  for exp in ${exps}; do
    echo "  exp:"${exp}
    for tmean in ${ofold}/${exp}/*tmean.nc; do
      fname=`basename ${tmean}`
      var=`ncdump -h ${tmean} | grep float | grep mean | tr '(' ' ' | awk '{print $2}'`
      vals='south_north:-1|west_east:-1,west_east,mean,XLAT'
      dims='bottom_top@bottom_top,south_north@south_north,west_east@west_east'
      if test ! -f ${ofold}/${exp}/${var}t-lonmean.nc; then
        python ${HOMEpy}/nc_var.py -o file_oper_alongdims -S ${vals} -f ${tmean} -v ${var}
        pyexec=$?
        if test ${pyexec} -ne 0; then
          echo ${errormsg}
          echo "  "${main}": python fails!"
          exit
        else
          mv file_oper_alongdims_mean.nc ${ofold}/${exp}/${var}t-lonmean.nc
        fi
#      exit
      fi
    done

  done
fi

if ${plot}; then
# Plots
##

for exp in ${exps}; do
  echo "exp: "${exp}
  for vark in ${varks}; do
    echo "  vark: "${vark}
    if test ${vark} = 'mean'; then
      gcoups=`echo ${graph_couples_mean} | tr ':' ' '`
    elif test ${vark} = 'sfc'; then
      if test ${exp} = 'AR40'; then
        gcoups=`echo ${graph_couples_sfc_AR40} | tr ':' ' '`
      else
        gcoups=`echo ${graph_couples_sfc_NPv31} | tr ':' ' '`
      fi
    elif test ${vark} = 'diag'; then
      gcoups=`echo ${graph_couples_diag} | tr ':' ' '`
    fi

    for gpair in ${gcoups}; do
      shdvar=`echo ${gpair} | tr '@' ' ' | awk '{print $1}'`
      cntvar=`echo ${gpair} | tr '@' ' ' | awk '{print $2}'`
      echo "  couple: "${shdvar}'-'${cntvar}
      shdvals=`python ${HOMEpy}/drawing.py -o variable_values -S ${shdvar} | grep all_values | awk '{print $3}'`
      cntvals=`python ${HOMEpy}/drawing.py -o variable_values -S ${cntvar} | grep all_values | awk '{print $3}'`
      files=${ofold}'/'${exp}'/'${shdvar}'mean.nc,'${ofold}'/'${exp}'/'${cntvar}'mean.nc'

      Lshdvals=`expr length ${shdvals}0`
      if test ${Lshdvals} -lt 2; then
        echo ${errormsg}
        echo "  Error in drawing_tools.py 'variables_values'!!"
        echo "    variable '"${shdvar}"' NOT found!"
        exit
      fi

      Lcntvals=`expr length ${cntvals}0`
      if test ${Lcntvals} -lt 2; then
        echo ${errormsg}
        echo "  Error in drawing_tools.py 'variables_values'!!"
        echo "    variable '"${cntvar}"' NOT found!"
        exit
      fi

      shdstdn=`echo ${shdvals} | tr ',' ' ' | awk '{print $1}'`
      shdcbar=`echo ${shdvals} | tr ',' ' ' | awk '{print $7}'`
      shdmin=`echo ${shdvals} | tr ',' ' ' | awk '{print $3}'`
      shdmax=`echo ${shdvals} | tr ',' ' ' | awk '{print $4}'`

      cntstdn=`echo ${cntvals} | tr ',' ' ' | awk '{print $1}'`
      cntmin=`echo ${cntvals} | tr ',' ' ' | awk '{print $3}'`
      cntmax=`echo ${cntvals} | tr ',' ' ' | awk '{print $4}'`

      if test ${shdstdn} = 'ERROR' || test ${cntstdn} = 'ERROR'; then
        echo ${errmsg}
        echo "  "${main}": wrong variable names!!!"
        echo "    shdvals: "${shdvals}
        echo "    cntvals: "${cntvals}
        echo "    shdvar: "${shdvar}" cntvar: "${cntvar}
        exit
      fi

      cntcolor='black'
      cntfmt='%g'

      if test ${gpair} = 'va@ua'; then
        shdmin=`echo ${shdmin} | awk '{print $1/4.}'`
        shdmax=`echo ${shdmax} | awk '{print $1/4.}'`
        cntmin=`echo ${cntmin} | awk '{print $1/3.}'`
        cntmax=`echo ${cntmax} | awk '{print $1 + 20.}'`
      elif test ${gpair} = 'tas@ps'; then
        shdmin=`echo ${shdmin} | awk '{print $1 + 50.}'`
        shdmax=`echo ${shdmax} | awk '{print $1 - 15.}'`
        cntmin=`echo ${cntmin} | awk '{print $1 + 15000.}'`
        cntmax=`echo ${cntmax} | awk '{print $1 - 2000.}'`
      elif test ${gpair} = 'uas@vas'; then
        cntmin=`echo ${cntmin} | awk '{print $1 * 0.1}'`
        cntmax=`echo ${cntmax} | awk '{print $1 * 0.1}'`
      elif test ${gpair} = 'pr@rsds'; then
#        shdmax=`echo ${shdmax} | awk '{print $1 / 20.}'`
        cntmax=`echo ${cntmax} | awk '{print $1 / 3.}'`
      elif test ${gpair} = 'clt@cll' || test ${gpair} = 'clh@clm'; then
        cntcolor='red'
        cntfmt='%.1g'
      elif test ${gpair} = 'clh@clm'; then
        cntmax=`echo ${cntmax} | awk '{print $1}'`
      elif test ${gpair} = 'prw@huss'; then
        shdmax=`echo ${shdmax} | awk '{print $1*4}'`
      elif test ${gpair} = 'prls@prc'; then
        shdmax=`echo ${shdmax} | awk '{print $1/1.}'`
        cntmax=`echo ${cntmax} | awk '{print $1 * 0.5}'`
      elif test ${gpair} = 'hfls@hfss'; then
        cntmin='-50.'
        cntmax='50.'
      fi
  
      if test ${vark} = 'mean'; then
        graphvals=${shdstdn}','${cntstdn}':z|-1,y|-1:z|-1,y|-1:lat:pressure:'
        graphvals=${graphvals}${shdcbar}':fixsigc,'${cntcolor}':'${cntfmt}':'${shdmin}','
        graphvals=${graphvals}${shdmax}':'${cntmin}','${cntmax}',9:LMDZ+WRF|'${exp}
        graphvals=${graphvals}'|meridional|monthly|average|of|'${shdstdn}'|&|'
        graphvals=${graphvals}${cntstdn}':pdf:flip@y:None'
      else
        graphvals=${shdstdn}','${cntstdn}':y|-1,time|-1:y|-1,time|-1:time:XLAT:'
        graphvals=${graphvals}${shdcbar}':fixsigc,'${cntcolor}':'${cntfmt}':'${shdmin}','
        graphvals=${graphvals}${shdmax}':'${cntmin}','${cntmax}',9:WRF+LMDZ|'${exp}
        graphvals=${graphvals}'|mean|meridional|evolution|of|'${shdstdn}'|&|'
        graphvals=${graphvals}${cntstdn}':pdf:None:time|hours!since!1949-12-01|'
        graphvals=${graphvals}'exct,5,d|%d|date!([DD]):None'
      fi

      echo ${files}
      echo ${graphvals}
      echo ${shdstdn}'mean,'${cntstdn}'mean'

      if test ${vark} = 'sfc' || test ${vark} = 'diag' && 
       ! $(isin_list ${coup2D} ${gpair}); then
        python ${HOMEpy}/drawing.py -f ${files} -o draw_2D_shad_cont_time -S ${graphvals} -v   \
          ${shdstdn}'mean,'${cntstdn}'mean'
        pyexc=$?
      else
        python ${HOMEpy}/drawing.py -f ${files} -o draw_2D_shad_cont -S ${graphvals} -v        \
          ${shdstdn}'mean,'${cntstdn}'mean'
        pyexc=$?
      fi
      if test ${pyexc} -ne 0; then
        echo ${errormsg}
        echo "  "${main}": drawing.py fails!"
        exit
      else
        mv 2Dfields_shadow-contour.pdf ${ofold}/${exp}/${shdvar}mean_${cntvar}mean.pdf
        evince ${ofold}/${exp}/${shdvar}mean_${cntvar}mean.pdf &
      fi

      if $(isin_list ${coup2D} ${gpair}); then
        graphvals=${shdstdn}','${cntstdn}':south_north|-1,west_east|-1:'
        graphvals=${graphvals}'south_north|-1,west_east|-1:longitude:latitude:'
        graphvals=${graphvals}${shdcbar}':fixsigc,'${cntcolor}':'${cntfmt}':'${shdmin}','
        graphvals=${graphvals}${shdmax}':'${cntmin}','${cntmax}',9:LMDZ+WRF|'${exp}
        graphvals=${graphvals}'|monthly|average|of|'${shdstdn}'|&|'
        graphvals=${graphvals}${cntstdn}':pdf:flip@y:None'
        echo '  '${shdstdn}'mean,'${cntstdn}'mean'

        python ${HOMEpy}/drawing.py -f ${files} -o draw_2D_shad_cont -S ${graphvals} -v        \
          ${shdstdn}'mean,'${cntstdn}'mean'
        pyexc=$?
        if test ${pyexc} -ne 0; then
          echo ${errormsg}
          echo "  "${main}": drawing.py fails!"
          exit
        else
          mv 2Dfields_shadow-contour.pdf ${ofold}/${exp}/${shdvar}tmean_${cntvar}tmean.pdf
          evince ${ofold}/${exp}/${shdvar}tmean_${cntvar}tmean.pdf &
        fi
      fi
#      exit

# end of couples
    done
# end of kinds
  done

done
fi

# 2D single plots
##
if ${single2Dplot}; then
  echo "2D single plots"
  for exp in ${exps}; do
    echo "  exp:"${exp}
    for tmean in ${ofold}/${exp}/*tmean.nc ;do
      fname=`basename ${tmean}`
      var=`ncdump -h ${tmean} | grep float | grep mean | tr '(' ' ' | awk '{print $2}'`
      if test ! ${var} = 'cltmean'; then
        shdvals=`python ${HOMEpy}/drawing.py -o variable_values -S ${var} |          \
          grep all_values | awk '{print $3}'`
        cbar=`echo ${shdvals} | tr ',' ' ' | awk '{print $7}'`
        min=`echo ${shdvals} | tr ',' ' ' | awk '{print $3}'`
        max=`echo ${shdvals} | tr ',' ' ' | awk '{print $4}'`

        if test ${var} = 'prlsmean'; then xtrms='0,0.00015';
        elif test ${var} = 'prwmean'; then xtrms='0,40';
        elif test ${var} = 'psmean'; then xtrms='99000,102500';
        elif test ${var} = 'r2mean'; then xtrms='0,0.025';
        elif test ${var} = 'tasmean'; then xtrms='275,310';
        elif test ${var} = 'tmlamean'; then xtrms='260,310';
        elif test ${var} = 'uasmean'; then xtrms='-20,20';
        elif test ${var} = 'vasmean'; then xtrms='-20,20';
        else xtrms=${min}','${max}; fi

#        vals=${var}':XLONG|-1,XLAT|-1:XLONG:XLAT:'${cbar}':'${xtrms}':monthly|mean:pdf:'
#        vals=${vals}'None:None:true'
#        dims='south_north@XLAT,west_east@XLONG'
        python ${HOMEpy}/nc_var.py -o WRF_CFlonlat_creation -S longitude,latitude    \
          -f ${tmean} -v XLONG,XLAT
        vals=${var}':longitude|-1,latitude|-1:longitude:latitude:'${cbar}':'${xtrms}':monthly|mean:pdf:'
        vals=${vals}'None:None:true'
        dims='south_north@latitude,west_east@longitude'
        python ${HOMEpy}/drawing.py -o draw_2D_shad -S ${vals} -f ${tmean} -v ${var}
        pyexec=$?
        if test ${pyexec} -ne 0; then
          echo ${errormsg}
          echo "  "${main}": python fails!"
          exit
        else
          mv 2Dfields_shadow.pdf ${ofold}/${exp}/${var}tmean.pdf
#          evince ${ofold}/${exp}/${var}tmean.pdf &
#          exit
        fi
      fi
    done
  done
fi

# lon 2 vertical section
##
echo "lon 2 vertical section ..."
vars=`echo ${lon2DZsecvars} | tr ':' ' '`
pts=`echo ${lon2DZsecpts} | tr ':' ' '`

if ${lonZsec}; then
  for exp in ${exps}; do
    file=${ofold}/${exp}/vertical_interpolation_WRFp.nc
#    python ${HOMEpy}/nc_var.py -o WRF_CFlonlat_creation -S longitude,latitude        \
#     -f ${file} -v lon,lat
## ALREADY done!
#    python ${HOMEpy}/nc_var.py -o valmod -S lowthres@oper,0.,sumc,360. -f ${file}    \
#       -v lon
    if test ${exp} = 'AR40'; then labexp='wlmdza'
    else labexp='wlmdzb'; fi

    for pt in ${pts}; do
      echo "  pt: "${pt}
      vals='x:15|y:'${pt}'|time:23'
      lonval=`python ${HOMEpy}/nc_var.py -o varout -f ${file} -S ${vals} -v lon |    \
        awk '{print $2}'`
      latval=`python ${HOMEpy}/nc_var.py -o varout -f ${file} -S ${vals} -v lat |    \
        awk '{print $2}'`
      tval=`python ${HOMEpy}/nc_var.py -o varout -f ${file} -S ${vals} -v time |     \
        awk '{print $2}'`
      for var in ${vars}; do
        echo "    var: "${var}

        shdvals=`python ${HOMEpy}/drawing.py -o variable_values -S ${var} |          \
          grep all_values | awk '{print $3}'`
        cbar=`echo ${shdvals} | tr ',' ' ' | awk '{print $7}'`
        min=`echo ${shdvals} | tr ',' ' ' | awk '{print $3}'`
        max=`echo ${shdvals} | tr ',' ' ' | awk '{print $4}'`

# WRFt,U,V,WRFrh,WRFght
        if test ${var} = 'WRFght'; then xtrms='0,40000';
        elif test ${var} = 'WRFt'; then xtrms='200,300';
        elif test ${var} = 'U'; then xtrms='-40,40';
        elif test ${var} = 'V'; then xtrms='-20,20';
        else xtrms=${min}','${max}; fi

# Plotting
        vals=${var}':x|-1,y|'${pt}',z|-1,time|24:lon:pressure:'${cbar}':'
        vals=${vals}${xtrms}':'${labexp}'|vertical|longitudinal|section|('${latval}
        vals=${vals}'$^{\circ}$):pdf:flip@y:None:true'
        python ${HOMEpy}/drawing.py -o draw_2D_shad -S ${vals} -f ${file} -v ${var}
        pyexec=$?
        if test ${pyexec} -ne 0; then
          echo ${errormsg}
          echo "  "${main}": python fails!"
          exit
        else
          mv 2Dfields_shadow.pdf ${ofold}/${exp}/${var}_lonZsec_${pt}pt.pdf
          evince ${ofold}/${exp}/${var}_lonZsec_${pt}pt.pdf &
        fi

        cntcolor='black'
        cntfmt='%g'

        if test ${var} = 'WRFght'; then cxtrms='0,40000';
        elif test ${var} = 'WRFt'; then cxtrms='200,300';
        elif test ${var} = 'U'; then cxtrms='-40,40';
        elif test ${var} = 'V'; then cxtrms='-20,20';
        else cxtrms=${min}','${max}; fi

        graphvals=${var}','${var}':x|-1,y|'${pt}',z|-1,time|24:'
        graphvals=${graphvals}'x|-1,y|'${pt}',z|-1,time|24:lon:pressure:'
        graphvals=${graphvals}${cbar}':fixsigc,'${cntcolor}':'${cntfmt}':'${xtrms}':'
        graphvals=${graphvals}${cxtrms}',15:'${labexp}
        graphvals=${graphvals}'|vertical|zonal|section|('${latval}'$^{\circ}$):pdf:'
        graphvals=${graphvals}'flip@y:None'

        python ${HOMEpy}/drawing.py -o draw_2D_shad_cont -S ${graphvals}             \
          -f ${file},${file} -v ${var},${var}
        pyexec=$?
        if test ${pyexec} -ne 0; then
          echo ${errormsg}
          echo "  "${main}": python fails!"
          exit
        else
          mv 2Dfields_shadow-contour.pdf ${ofold}/${exp}/${var}_lonZsec-cnt_${pt}pt.pdf
          evince ${ofold}/${exp}/${var}_lonZsec-cnt_${pt}pt.pdf &
        fi

#        exit
      done
    done
  done
fi

echo "Computing differences"
exp1=`echo ${experiments} | tr ':' ' ' | awk '{print $1}'`
exp2=`echo ${experiments} | tr ':' ' ' | awk '{print $2}'`

diffks=`echo ${diffkinds} | tr ':' ' '`

if test ${differences}; then
  for kind in ${diffks}; do
    echo ${kind}
    case ${kind} in
      'diffZ')
        vars=`echo ${diffZvars} | tr ':' ' '`
        fileorigd='mean.nc'
      ;;
      'diffH')
        vars=`echo ${diffHvars} | tr ':' ' '`
        fileorigd='tmean.nc'
      ;;
    esac

    echo "  "${var}
    for var in ${vars}; do
      if test ! -f ${ofold}/${var}_diff.nc; then
#        cdo sub ${ofold}/${exp1}/${var}${fileorigd} ${ofold}/${exp2}/${var}${fileorigd} ${ofold}/${var}_diff.nc
        if test ${kind} = 'diffZ'; then
          values='pressure|lat@add|'${ofold}'/'${exp1}'/'${var}${fileorigd}'|'${var}'mean,'
          values=${values}'sub|'${ofold}'/'${exp2}'/'${var}${fileorigd}'|'${var}'mean'
          python ${HOMEpy}/nc_var.py -o compute_opersvarsfiles -S ${values} -v ${var}
          pyexec=$?
          if test ${pyexec} -ne 0; then
            echo ${errormsg}
            echo "  "${main}": python fails!"
            exit
          fi
          mv opersvarsfiles_${var}.nc ${ofold}/${var}_diff.nc
        else
          values='longitude|latitude@add|'${ofold}'/'${exp1}'/'${var}${fileorigd}'|'${var}'mean,'
          values=${values}'sub|'${ofold}'/'${exp2}'/'${var}${fileorigd}'|'${var}'mean'
          python ${HOMEpy}/nc_var.py -o compute_opersvarsfiles -S ${values} -v ${var}
          pyexec=$?
          if test ${pyexec} -ne 0; then
            echo ${errormsg}
            echo "  "${main}": python fails!"
            exit
          fi
          mv opersvarsfiles_${var}.nc ${ofold}/${var}_diff.nc
        fi
      fi
    done

    if test ${kind} = 'diffZ'; then
      coups=`echo ${graph_couples_mean} | tr ':' ' '`
      for coup in ${coups}; do
        shdvar=`echo ${coup} | tr '@' ' ' | awk '{print $1}'`
        cntvar=`echo ${coup} | tr '@' ' ' | awk '{print $2}'`

        echo "  couple: "${shdvar}'-'${cntvar}

        shdvals=`python ${HOMEpy}/drawing.py -o variable_values -S ${shdvar} | grep all_values | awk '{print $3}'`
        cntvals=`python ${HOMEpy}/drawing.py -o variable_values -S ${cntvar} | grep all_values | awk '{print $3}'`
        files=${ofold}'/'${shdvar}'_diff.nc,'${ofold}'/'${cntvar}'_diff.nc'

        shdstdn=`echo ${shdvals} | tr ',' ' ' | awk '{print $1}'`
        shdcbar=`echo ${shdvals} | tr ',' ' ' | awk '{print $7}'`
        shdmin=`echo ${shdvals} | tr ',' ' ' | awk '{print $3}'`
        shdmax=`echo ${shdvals} | tr ',' ' ' | awk '{print $4}'`

        cntstdn=`echo ${cntvals} | tr ',' ' ' | awk '{print $1}'`
        cntmin=`echo ${cntvals} | tr ',' ' ' | awk '{print $3}'`
        cntmax=`echo ${cntvals} | tr ',' ' ' | awk '{print $4}'`

        shdcbar='seismic'
        if test ${coup} = 'va@ua'; then
          shdmin=-4.
          shdmax=4.
          cntmin=-10.
          cntmax=10.
        elif test ${coup} = 'hus@ta'; then 
          shdmin=-0.2
          shdmax=0.2
          cntmin=-2.
          cntmax=2.
        fi

        cntcolor='black'
        cntfmt='%g'

        graphvals=${shdstdn}','${cntstdn}':x|-1,y|-1,x_2|-1:x|-1,y|-1,x_2|-1:lat:'
        graphvals=${graphvals}'pressure:'${shdcbar}':fixsigc,'${cntcolor}':'${cntfmt}':'
        graphvals=${graphvals}${shdmin}','${shdmax}':'${cntmin}','${cntmax}',9:'
        graphvals=${graphvals}${exp1}'-'${exp2}'|meridional|monthly|average|differences|of|'
        graphvals=${graphvals}${shdstdn}'|&|'${cntstdn}':pdf:flip@y:None'

        python ${HOMEpy}/drawing.py -o draw_2D_shad_cont -S ${graphvals}             \
          -f ${files} -v ${shdvar},${cntvar}
        pyexec=$?
        if test ${pyexec} -ne 0; then
          echo ${errormsg}
          echo "  "${main}": python fails!"
          exit
        else
          mv 2Dfields_shadow-contour.pdf ${ofold}/${shdvar}-${cntvar}_diff.pdf
          evince ${ofold}/${shdvar}-${cntvar}_diff.pdf &
        fi
      done
#      exit
    else
      for var in ${vars}; do
        echo "  "${var}
        vals=`python ${HOMEpy}/drawing.py -o variable_values -S ${var} | grep all_values | awk '{print $3}'`
        file=${ofold}/${var}_diff.nc

        stdn=`echo ${vals} | tr ',' ' ' | awk '{print $1}'`
        cbar='seismic'

        if test ${var} = 'uas'; then xtrms='-10.,10.';
        elif test ${var} = 'vas'; then xtrms='-5,5';
        elif test ${var} = 'ps'; then xtrms='-500,500';
        elif test ${var} = 'pr'; then xtrms='-0.001,0.001';
        else xtrms=${min}','${max}; fi

        vals=${var}':x|-1,y|-1:longitude:latitude:'${cbar}':'
        vals=${vals}${xtrms}':'${exp1}'-'${exp2}'|meridional|monthly|average|differences:'
        vals=${vals}'pdf:None:None:true'
        python ${HOMEpy}/drawing.py -o draw_2D_shad -S ${vals} -f ${file} -v ${var}
        pyexec=$?
        if test ${pyexec} -ne 0; then
          echo ${errormsg}
          echo "  "${main}": python fails!"
          exit
        else
          mv 2Dfields_shadow.pdf ${ofold}/${var}_diff.pdf
          evince ${ofold}/${var}_diff.pdf &
        fi
#        exit
      done
    fi
  done
fi
