#!/bin/bash
## e.g. # concatenate.bash hydrostatic -1:14:2:25
if test $1 = '-h'; then
  echo "*******************************"
  echo "*** Files WRF concatenation ***"  
  echo "*******************************"
  echo "concatenate.bash [fold](fodler to use) [point](as Time:z:y:x, -1 last value)"
else
function isin_list() {
# 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
}

####### ###### ##### #### ### ## #
  allvars='XLONG,XLAT,Times,LDQAJS,LDTAJS,LDQCON,LDTCON,LDUCON,LDVCON,LDQDYN,LDTDYN,LDUDYN,LDVDYN,LDQEVA,LDTEVA,LDQLSCST,LDTLSCST,LDQLSCTH,LDTLSCTH,LDQLSC,LDTLSC,LDQPHY,LDTPHY,LDQTHE,LDTTHE,LDQVDF,LDTVDF,LDUVDF,LDVVDF,LDQWAK,LDTWAK,QVAPOR,QCLOUD,LGEOP,LTEMP,LPRES,LVITU,LVITV,LVITW,LVPRECIP,LRHUM,LH2O,LPR_CON_I,LPR_CON_L,LIWCON,LPTCONV,LSNOW,LEVAP,LPRECIP,LPLUC,LPLUL,LRH2M,T2,U10,V10,LCLDL,LCLDM,LCLDH,PSFC'

# Variables to be printed centered at zero
  ZEROcenteredvars='LDQAJS:LDTAJS:LDQCON:LDTCON:LDUCON:LDVCON:LDQDYN:LDTDYN:LDUDYN:LDVDYN:LDQEVA:LDTEVA:LDQLSCST:LDTLSCST:LDQLSCTH:LDTLSCTH:LDQLSC:LDTLSC:LDQPHY:LDTPHY:LDQTHE:LDTTHE:LDQVDF:LDTVDF:LDUVDF:LDVVDF:LDQWAK:LDTWAK:LVITU:LVITV:LVITW:U10:V10'

# Variables printed following min/max or quantiles
  MINMAXvars='QVAPOR:QCLOUD:LGEOP:LTEMP:LPRES:LVPRECIP:LRHUM:LH2O:LPR_CON_I:LPR_CON_L:LIWCON:LPTCONV:LSNOW:LEVAP:LPRECIP:LPLUC:LPLUL:LRH2M:T2:LCLDL:LCLDM:LCLDH:PSFC'

# 1D variables
  OneDvars='LSNOW:LEVAP:LPRECIP:LPLUC:LPLUL:LRH2M:T2:U10:V10:LCLDL:LCLDM:LCLDH:PSFC'

  NOdraw='XLONG:XLAT:Times:longitude:latitude:time:intTime:varDIMbottom_top'

# File creation
  file=0

# Cutting file
  cut=0

# Drawing
  draw=0

# Graphic visualization
  vis=0

  hovzvarns='dqajs:dtajs:dqcon:dtcon:ducon:dvcon:dqdyn:dtdyn:dudyn:dvdyn:dqeva:dteva:dqlscst:dtlscst:dqlscth:dtlscth:dqlsc:dtlsc:dqphy:dtphy:dqthe:dtthe:dqvdf:dtvdf:duvdf:dvvdf:dqwak:dtwak:r:c:zg:ta:pres:ua:va:wap:prprof:hus:h2o:prci:prcl:ri:cpt'

  dh2varns='dqajs:dtajs:dqcon:dtcon:ducon:dvcon:dqdyn:dtdyn:dudyn:dvdyn:dqeva:dteva:dqlscst:dtlscst:dqlscth:dtlscth:dqlsc:dtlsc:dqphy:dtphy:dqthe:dtthe:dqvdf:dtvdf:duvdf:dvvdf:dqwak:dtwak:r:c:zg:ta:pres:ua:va:wap:prprof:hus:h2o:prci:prcl:ri:cpt:prsn:evspsbl:pr:prc:prls:huss:tas:uas:vas:cll:clm:clh:ps'

  tevolvarns='dqajs:dtajs:dqcon:dtcon:ducon:dvcon:dqdyn:dtdyn:dudyn:dvdyn:dqeva:dteva:dqlscst:dtlscst:dqlscth:dtlscth:dqlsc:dtlsc:dqphy:dtphy:dqthe:dtthe:dqvdf:dtvdf:duvdf:dvvdf:dqwak:dtwak:r:c:zg:ta:pres:ua:va:wap:prprof:hus:h2o:prci:prcl:ri:cpt:prsn:evspsbl:pr:prc:prls:huss:tas:uas:vas:cll:clm:clh:ps'

# Tex
  tex=1
  Ncoldraw=5
 
  graphkinds='hor:Hovmoller:TimeSeries' 

  kvalues='2dh:hovz:tevol'

#######    #######
## MAIN
    #######
  rootsh=`pwd`
  PYhome=${HOME}/PY
  main='concatenation.bash'
  errormsg='ERROR -- error -- ERROR -- error'
  infold=`pwd`
  fold=$1
  point=$2

  ofold=${fold}

  tval=`echo ${point} | tr ':' ' ' | awk '{print $1}'`
  zval=`echo ${point} | tr ':' ' ' | awk '{print $2}'`
  yval=`echo ${point} | tr ':' ' ' | awk '{print $3}'`
  xval=`echo ${point} | tr ':' ' ' | awk '{print $4}'`

  kvals=`echo ${kvalues} | tr ':' ' '`
  graphks=`echo ${graphkinds} | tr ':' ' '`

  if test ${file} -eq 1; then
    python ${PYhome}/nc_var.py -o netcdf_fold_concatenation -f wrfout_d01 -S ${infold}/${fold},Time -v ${allvars}

    if test $? -ne 0; then
      echo ${errormsg}
      echo "  "${main}"python has failed!!"
      echo "    python ${PYhome}/nc_var.py -o netcdf_fold_concatenation -f wrfout_d01\
      -S ${infold}/${fold},Time -v ${allvars}"
      exit
    fi
    mv netcdf_fold_concatenated.nc ${infold}/${fold}/wrfout_d01_check.nc
  fi

  if test ${cut} -eq 1; then
    infile=${infold}/${fold}/wrfout_d01_check.nc
    for kvn in ${kvals}; do
      echo "  kind of values: "${kvn}
      case $kvn in
        '2dh')
          pyvls='Time,'${tval}','${tval}',1@bottom_top,'${zval}','${zval}',1'
        ;;
        'hovz')
          pyvls='south_north,'${yval}','${yval}',1@west_east,'${xval}','${xval}',1'
        ;;
        'tevol')
          pyvls='bottom_top,'${zval}','${zval}',1@south_north,'${yval}','${yval}
          pyvls=${pyvls}',1@west_east,'${xval}','${xval}',1'
        ;;
      esac
      ofile=`python ${PYhome}/nc_var.py -o DataSetSection_multidims -f ${infile}     \
        -S ${pyvls} | grep succesfull | awk '{print $6}' | tr '"' ' '`
      if test $? -ne 0; then
        echo ${errmsg}
        echo "  "${main}": python failed!"
        echo "    python "${PYhome}/nc_var.py -o DataSetSection_multidims              \
          -f ${infile} -S ${pyvls}
        exit
      fi
      case $kvn in
        '2dh')
# Adding CF longitudes and latitudes
          python ${PYhome}/nc_var.py -f ${ofile} -o WRF_CFlonlat_creation            \
            -S longitude,latitude -v XLONG,XLAT
        ;;
        'tevol')
# Adding CF time
          python ${PYhome}/nc_var.py -f ${ofile} -o WRF_CFtime_creation              \
            -S 19491201000000,hours -v time
        ;;
        'hovz')
# Adding CF time
          python ${PYhome}/nc_var.py -f ${ofile} -o WRF_CFtime_creation              \
            -S 19491201000000,hours -v time
# Adding variable dimension from 'bottom_top'
          python ${PYhome}/nc_var.py -f ${ofile} -o dimVar_creation -S bottom_top
        ;;
      esac
    done
  fi

  if test $draw -eq 1; then
  echo "Drawing..."
  mkdir -p ${ofold}/figs
  infile=${fold}'/wrfout_d01_check'

  for kvn in ${kvals}; do
    case $kvn in
      '2dh')
        filetail='Time_B'${tval}'-E'${tval}'-I1_bottom_top_B'${zval}
        filetail=${filetail}'-E'${zval}'-I1'
        dh2varns=''
      ;;
      'hovz')
        filetail='south_north_B'${yval}'-E'${yval}'-I1_west_east_B'${xval}
        filetail=${filetail}'-E'${xval}'-I1'
        hovzvarns=''
      ;;
      'tevol')
        filetail='bottom_top_B'${zval}'-E'${zval}'-I1_south_north_B'${yval}'-E'
        filetail=${filetail}${yval}'-I1_west_east_B'${xval}'-E'${xval}'-I1'
        tevolvarns=''
      ;;
    esac
    drwfile=${infile}_${filetail}.nc
    var0s=`python ~/PY/nc_var.py -o ivars -f ${drwfile} | grep allvars | awk '{print $3}'`
    vars=`echo ${var0s} | tr ':' ' '`
    echo ${vars}
    for var in ${vars}; do
      if ! $(isin_list ${NOdraw} ${var}); then
      echo "  var: "${var}
      varvals=`python ${PYhome}/drawing.py -o variable_values -S ${var} | grep all_values`
      if test $? -ne 0; then
        echo "  "${main}": Variable '"${var}"' not defined !!!!"
        exit
      fi
      stdn=`echo ${varvals} | tr ',' ' ' | awk '{print $3}'`
      cbar=`echo ${varvals} | tr ',' ' ' | awk '{print $9}'`
      min=`echo ${varvals} | tr ',' ' ' | awk '{print $5}'`
      max=`echo ${varvals} | tr ',' ' ' | awk '{print $6}'`

      if test ${kvn} = 'hovz' && ! $(isin_list ${OneDvars} ${var}); then
# Shadded vertical Hovmoller diagrams
##
      hovzvarns=${hovzvarns}':'${stdn}
# Let's make it easy
      if $(isin_list ${ZEROcenteredvars} ${var}); then
        min='Spercentile@1'
        max='Spercentile@1'
        cbar='seismic'
      else
        min='Saroundminmax@0.8'
        max='Saroundminmax@0.8'
      fi 

      if test ${var} = 'LDQEVA'; then
        cbar='Reds'
      fi

      echo "    min: "${min}" max: "${max}

      values=${stdn}':Time|-1,bottom_top|-1:varDIMbottom_top:time:'${cbar}':'
      values=${values}${min}','${max}':vertical|temporal|evolution|:pdf:'
      values=${values}'transpose:time|hours!since!1949-12-01|exct,12,h|'
      values=${values}'$%d^{%H}$|date!($[DD]^{[HH]}$)|x:True'

      python ${PYhome}/drawing.py -o draw_2D_shad_time -f ${infile}_${filetail}.nc   \
        -S ${values} -v ${var},time
      if test $? -ne 0; then
        echo ${errormsg}
        echo "  "${main}": drawing.py fails!"
        echo python ${PYhome}/drawing.py -o draw_2D_shad_time \
          -f ${infile}_${filetail}.nc -S ${values} -v ${var},time
        exit
      else
        mv 2Dfields_shadow_time.pdf ${ofold}/figs/${stdn}_Hovmoller.pdf
        if test $vis -eq 1; then evince ${ofold}/figs/${stdn}_Hovmoller.pdf & fi
      fi
      
      elif test ${kvn} = 'tevol'; then
# Point evolution
      tevolvarns=${tevolvarns}':'${stdn}

      values=${stdn}':time|($[DD]^{[HH]}$):exct,12,h:$%d^{%H}$:time|evolution|at'
      values=${values}'|(x='${xval}',y='${yval}',z='${zval}'):1:pdf'
      python ${PYhome}/drawing.py -o draw_timeSeries -f ${infile}_${filetail}.nc     \
        -S ${values} -v ${var},time
      if test $? -ne 0; then
        echo ${errormsg}
        echo "  "${main}": drawing.py fails!"
        echo python ${PYhome}/drawing.py -o draw_timeSeries \
          -f ${infile}_${filetail}.nc -S ${values} -v ${var},time
        exit
      else
        mv TimeSeries_${stdn}.pdf ${ofold}/figs/${stdn}_TimeSeries.pdf
        if test $vis -eq 1; then evince ${ofold}/figs/${stdn}_TimeSeries.pdf & fi
      fi

      elif test ${kvn} = '2dh'; then
# Horizontal plot
      dh2varns=${dh2varns}':'${stdn}
      mapv='None'

      min=`echo ${varvals} | tr ',' ' ' | awk '{print $5}'`
      max=`echo ${varvals} | tr ',' ' ' | awk '{print $6}'`

      if $(isin_list ${ZEROcenteredvars} ${var}); then
        min='Spercentile@1'
        max='Spercentile@1'
        cbar='seismic'
      else
        min='Saroundminmax@0.97'
        max='Saroundminmax@0.97'
      fi 

      if test ${var} = 'LRH2M'; then
        max='100.'
      elif test ${var} = 'T2'; then
        min='270'
        max='300.'
      fi

      echo "    min: "${min}" max: "${max}

      values=${stdn}':south_north|-1,west_east|-1:longitude:latitude:'${cbar}':'
      values=${values}${min}','${max}':horizontal|section|at|(t='${tval}
      values=${values}',z='${zval}'):pdf:None:'${mapv}':True'

      python ${PYhome}/drawing.py -o draw_2D_shad -f ${infile}_${filetail}.nc        \
        -S ${values} -v ${var}
      if test $? -ne 0; then
        echo ${errormsg}
        echo "  "${main}": drawing.py fails!"
        echo python ${PYhome}/drawing.py -o draw_2D_shad \
          -f ${infile}_${filetail}.nc -S ${values} -v ${var}
        exit
      else
        mv 2Dfields_shadow.pdf ${ofold}/figs/${stdn}_hor.pdf
        if test $vis -eq 1; then evince ${ofold}/figs/${stdn}_hor.pdf & fi
      fi

#      exit
      fi
      fi
    done
  done
  fi

  if test ${tex} -eq 1; then
# Tex generation
    echo "  TeX file generation...."

    otex='concatenate.tex'

    icol=1
    while test ${icol} -le ${Ncoldraw}; do
      NcolS=${NcolS}'c'
      icol=`expr ${icol} + 1`
    done

    colfrac=`echo "1 "${Ncoldraw} | awk '{print $1/$2}'`
    maxlinegraph=`echo "2 "${Ncoldraw} | awk '{printf("%d",$2*sqrt($1))}'`
    echo "%concatenate LaTeX output" > ${otex}

    for graphk in ${graphks}; do
      cat << EOF >> ${otex}

\\begin{figure}
\\begin{center}
\\begin{tabular}{${NcolS}}
EOF
      echo "  graphk: "${graphk}

      case ${graphk} in
        'hor')
          vardrwns=${dh2varns}
          GraphName='Horizontal section at time='${tval}', zlev='${zval}
        ;;
        'TimeSeries')
          vardrwns=${tevolvarns}
          GraphName='Punctual temporal evolution at zlev='${zval}', ypt= '${yval}', xpt= '${xval}
        ;;
        'Hovmoller')
          vardrwns=${hovzvarns}
          GraphName='Hovmoeller vertical section evolution at ypt= '${yval}', xpt= '${xval}
        ;;
      esac
      echo "    vardrwns: "${vardrwns}
      vardns=`echo ${vardrwns} | tr ':' ' '`

      Ngraph=1
      Nline=1
      Ngraphk=1
      varsgraph=''
      for graph in ${vardns}; do
        varvals=`python ${PYhome}/drawing.py -o variable_values -S ${var} | grep all_values`
        stdn=`echo ${varvals} | tr ',' ' ' | awk '{print $3}'`

        graphn=${ofold}'/'figs'/'${graph}'_'${graphk}'.png'
        echo "\includegraphics[width="${colfrac}"\textwidth]{"${graphn}"}" >> ${otex}
        modNgraph=`expr ${Ngraph} % ${Ncoldraw}`
        if test ${modNgraph} -eq 0; then
          echo "\\\\" >> ${otex}
          Nline=`expr ${Nline} + 1`
          if test ${Nline} -gt ${maxlinegraph}; then
            Nline=1
            varsgraph=${varsgraph}', '${graph}
            cat << EOF >> ${otex}
\\end{tabular}
\\end{center}
\\caption{${GraphName} graphics from top left to bottom right of ${varsgraph}}
\\label{fig:${graphk}_${Ngraphk}}
\\end{figure}
\\clearpage

\\begin{figure}
\\begin{center}
\\begin{tabular}{${NcolS}}
EOF

           Ngraphk=`expr ${Ngraphk} + 1`
           varsgraph=''
          fi
        else
          echo "&" >> ${otex}
          varsgraph=${varsgraph}', '${graph}
        fi
        Ngraph=`expr ${Ngraph} + 1`
      done
      cat << EOF >> ${otex}
\\end{tabular}
\\end{center}
\\caption{${GraphName} graphics from top left to bottom right of ${varsgraph}}
\\label{fig:${graphk}_${Ngraphk}}
\\end{figure}
\\clearpage

EOF

    done
  fi

fi