#!/bin/bash

#. env.sh

local=$( pwd )

#set -vx

########################################################################
#
# Computing one metric from netcdf 2D z-t file
# --------------------------------------------
#
# Inheritating from extract_onemetric.sh
#
# Compute average btwn t1,t2,zmin,zmax for ncdf file
# Auteur: C Rio, F Hourdin
# Modified: F Couvreux
# Modified: 23/01/2019 N Villefranque from r110
# Modified: 23/09/2022 N Villefranque from r318
# Changed to get_one_metric_from_file.sh : 10/08/2024 F Hourdin & N Villefranque
#
########################################################################

##############################################################
# A faire
# Gestion des options.
# Notamment shift_time
# Pour l'instant le 3eme argument
##############################################################

shift_time=+0

while (( $# > 0 )) ; do
    case $1 in
       -h|--help) cat <<eod
Use : $0 [-h] metric_name dir_or_file [N_shift]
if dir_or_file is : a netcdf file, compute the metric for this file
                  : a directory, compute the metrics for all the files in it
dir_or_file is a path like MAIN/CASE/SUBCASE with main within WAVE1, 2... or LES or RAD or CTRL
N_shift : to shift by N_shift time steps the time axis when computing the metrics
eod
       shift ;;
       *) input_file_or_dir=$1 ; metric_full_name=$2 ; shift 2
          if [ $# -gt 0 ] ; then shift_time=$1 ; fi ; shift
    esac
done

if [ -d $input_file_or_dir ] ; then
    input_nc_files="$( ls $input_file_or_dir/*.nc )"
elif [ -f $input_file_or_dir ] ; then
    input_nc_files=$input_file_or_dir
else
    echo Error in $0 : $input_file_or_dir not a file nor a directory ; $0 -h ; exit 1
fi


sepm=_

# if radiative metric, metric_full_name is
# RAD_CASE_SUBCASETIME_metrics_SZA1_SZA2
# ex: RAD_RICO_REF005_uptoa_01_01

#dir=$2 # = LES/RAD or WAVEi


prefile_def=SCM

# metric_output_name : LES_cas_souscas_... ou RAD_RAD_cas_souscastemps_... ou WAVE1_cas_souscas ou WAVE1_RAD_cas_souscastemps_...
# metric_output_name=${dir}_$metric_full_name
casename=$case_basis/$name_subcase
# echo Cas $casename, metric $metric_name, time $tmin : $tmax, treating $dir
# echo $metric_output_name # same info as above + subcase (+ time if RAD)

case_basis=$( echo $input_file_or_dir | cut -d/ -f2 )
case $case_basis in
	ARMCU|RICO|IHOP) zmax_domain=3500 ;;
	SANDU) zmax_domain=2800 ;;
	GABLS1|GABLS4) zmax_domain=400 ;;
	*) zmax_domain=2000
esac

#########################################################################
export LC_ALL=C # set locale for computation with awk
#########################################################################

# only print the number without header or anything
NCKS="ncks -H -C -s %15.8f "

#########################################################################
# Decompisition of inputs
#########################################################################
#  IFS='_' read -a temp <<< $metric_full_name : TOO SOPHISTICATED FOR MACOS
temp=( $( echo $metric_full_name | sed -e 's/'${sepm}'/ /g' ) )
metric_name=${temp[0]}
tmin=${temp[1]}
tmax=${temp[2]} 


############# Vieux commentaires ########################################
#echo metric_full_name $metric_full_name
#echo $temp
# si metrique non radiative
# temp = ( CAS SOUSCAS METRIC T1 T2 )
# si metrique radiative
# temp = ( RAD CAS SOUSCASTIME METRIC SZA1 SZA2 )

# Decomposition soit en
# cas_souscas_metrique_t1_t2 => comparaison SCM / LES ou LES / LES
# RAD_cas_souscas_metrique_t2_t2 => comparaison ecRad / MC sur le champ de cas/souscas, un souscas par instant, t2 = i_sza
# Attention de ne pas appeler un cas "RAD"

#tmp=${temp[0]}
#temp=( $( echo $input_nc_file | sed -e 's/'${sepm}'/ /g' ) )
#if [ ${tmp} == RAD ] # 3 first characters
#then
#  REF=RAD
#  case_basis=${temp[1]}
#  name_subcase=${temp[2]:0:-3}
#  time=${temp[2]: -3}
#  prefile_def=RAD${time}
#  #INS=${REF}${name_subcase: -3}
#else
#  REF=LES
#  case_basis=${temp[0]}
#  name_subcase=${temp[1]}
#fi




#########################################################################
# nWAVE=${dir:4}

# a changer en fonction du repertoire de l'utilisateur
# DIRin=$dir/$casename
# if [ $REF == RAD -a $dir == RAD ] ; then
#   DIRin=$DIRin/$dir_rad
# fi
# ex: LES/ARMCU/REF ou RAD/ARMCU/REF ou WAVE1/ARMCU/REF

########################################################################
# list of files to be post processed
# LES or SCM
# trick to compute LES var when one LES only is available
# by shifting the time of the metrics by +1 or -1
# Default : no shift (shf=+0)
########################################################################
# cd $DIRin
# shfs="+0"
# #if [ -f LES0.nc ] ; then
#  prefile=LES
#  if [ ! -f LES1.nc ] ; then shfs="+0 -1 +1" ; fi
#elif [ -f LESLESSCM_${nWAVE}-001.nc ] ; then
#  prefile=LESLESSCM_${nWAVE}-
#elif [[ `pwd` == *"/RAD/"* ]] ; then # on est dans RAD/
#  if [ -f RAD${time}.nc ] ; then
#    prefile=RAD
#  else 
#    echo "error: NO REFERENCE FOR THE METRIC" $metric_full_name
#    exit 1
#  fi
#elif [ -f SCM.nc ] ; then
#  # Gestion assez horrible. A retravailler.
#  # Toutes lees autres simulations on des numéros à la fin.
#  # Mais pas le controle qui s'appelle SCM.nc
#  prefile="S"
#else
#  # if scm, prefile_def = SCM
#  # if rad, prefile_def = RAD$time (RAD005)
#  prefile=${prefile_def}-${nWAVE}-
#fi


#if [[ `pwd` == *"/RAD/"* ]] ; then
#  # many RADXXX in CASE/SUBCASE, use only one
#  ls ${prefile}${time}*nc > list
#else
#  ls ${prefile}*nc > list
#fi
#sed -e "s/$prefile//" -e "s/.nc//" list > list_nruns
#nruns=`wc -l list_nruns | awk '{print $1}'`
#nr=0
#cd $local

########################################################################
# Loop on metrics
########################################################################

#echo $metric_output_name
#file=$( echo $input_nc_file | sed -e 's/.nc//' )
for ncFile in $input_nc_files ; do
#for run in `cat $DIRin/list_nruns`; do
    #run3=${run}
    #file=$DIRin/$prefile${run3}
    #(( nr = $nr + 1 ))
    #################################################################
    #defini les bornes temporelles et verticales
    #################################################################
    t1=`echo $tmin | awk ' { print $1 - 1 '$shift_time' } '`
    t2=`echo $tmax | awk ' { print $1 - 1 '$shift_time' } '`
    t1R=`echo $tmin | awk ' { print $1 '$shift_time' } '`
    t2R=`echo $tmax | awk ' { print $1 '$shift_time' } '`

    #if [ $REF == RAD ]
    if [[ "${metric_name}" == "RAD-"* ]]
    then 
      metric_name=$( echo $metric_name | sed -e 's/RAD-//' )
      #if [ $dir == RAD ] 
      ncdump -h $ncFile | grep vertical_levels > /dev/null 2>&1
      if [ $? -eq 0 ]
      then # In Monte Carlo netCDF files
        vlev=vertical_levels
        vlev_toa=-1
        vlev_surf=0
        calc_std=1
      else # In ecRad netCDF files (reverse z order and head / tail are different )
        vlev=half_level
        vlev_toa=0
        vlev_surf=-1
        calc_std=0
      fi
      get_rad_met(){
        ilev=$1
        varname=$2
        $NCKS -d column,$t1,$t2 -d $vlev,$ilev -v $varname $ncFile|head -1|awk '{for (i=1; i<=NF; i++) sum+=$i; print sum/NF}'
      }
      if [ "$metric_name" = "lwu" -o "$metric_name" = "lwd" ] ; then
        if [ $calc_std -eq 1 ] ; then # MC
          std_lwtoa=$(get_rad_met $vlev_toa std_flux_up_lw)
          std_lwsurf=$(get_rad_met $vlev_surf std_flux_dn_lw)
        else # ecRad => select right height because vertical levels don't correspond !!! THIS IS BAD CODING !!!
          vlev_toa=30
          vlev_surf=-4
        fi
        lwsurf=$(get_rad_met $vlev_surf flux_dn_lw)
        lwtoa=$(get_rad_met $vlev_toa flux_up_lw)
      else # SW metrics
        dntoa=$(get_rad_met $vlev_toa flux_dn_sw)
        uptoa=$(get_rad_met $vlev_toa flux_up_sw)
        dnsurf=$(get_rad_met $vlev_surf flux_dn_sw)
        upsurf=$(get_rad_met $vlev_surf flux_up_sw)
        if [ $calc_std -eq 1 ] ; then
          std_dntoa=$(get_rad_met $vlev_toa std_flux_dn_sw)
          std_uptoa=$(get_rad_met $vlev_toa std_flux_up_sw)
          std_dnsurf=$(get_rad_met $vlev_surf std_flux_dn_sw)
          std_upsurf=$(get_rad_met $vlev_surf std_flux_up_sw)
        fi
      fi
    fi
    case ${metric_name:0:3} in
      "unt") # untrans(mitted)
        metric=`echo $dntoa $dnsurf | awk '{ print $1 - $2 }' `
        ;;
      "dns") metric=$dnsurf ;;
      "upt") metric=$uptoa  ;;
      "dnt") metric=$dntoa  ;;
      "abs") metric=`echo $dntoa $dnsurf $upsurf $uptoa | awk '{ print $1 - $2 + $3 - $4 }' ` ;;
      "lwu") metric=$lwtoa  ;;
      "lwd") metric=$lwsurf ;;

      #################################################################
      # Metrics computed in R with htune_netcdf2csvMetrics.R
      "neb"|"lwp"|"Ay-"|"net"|"rat"|"tra") 
      #################################################################
	      #metric=`Rscript --vanilla htune_netcdf2csvMetrics.R $ncFile $metric_name $t1R $t2R $zmax_domain | awk ' { print $2 } '`
	      metric=`Rscript htune_netcdf2csvMetrics.R $ncFile $metric_name $t1R $t2R $zmax_domain | awk ' { print $2 } '`
	      ;;

      #################################################################
      # Metrics computed in R with htune_netcdf2csvMetrics.R
      "noi")
      #################################################################
              tmp_=( `echo $metric_name | sed -e 's/-/ /g'` )
              var_=${tmp_[1]}
              z1_=${tmp_[2]}
              z2_=${tmp_[3]}
              # echo tmp_ ${tmp_[*]}
              # echo var $var_
	      # echo Rscript htune_noise_metrics.R $ncFile $var_ $t1R $t2R $z1_ $z2_
	      metric=`Rscript htune_noise_metrics.R $ncFile $var_ $t1R $t2R $z1_ $z2_ | awk ' { print $2 } '`
	      ;;

      #################################################################
      # case of z-t box averages ("zav")
      # Metrics computed with cdo
      "zav")
      #################################################################
	      z1=`echo $metric_name | awk -F- ' { print $2 } '`
	      z2=`echo $metric_name | awk -F- ' { print $3 } '`
	      met=`echo $metric_name | awk -F- ' { print $4} '`
	   
        rm -f XXXX_intermediate_extract*.nc
        timename=`ncdump -h $ncFile | grep -i time |head -1 | awk ' {print $1 } '`
        # defining a mask in altitude (named "zf" in nc file) between z1 and z2.
	ncap2 -s "mask= (zf >= ${z1} && zf <= ${z2})" $ncFile XXXX_intermediate_extract.nc
        #echo "execute ncwa to perform time and z average"
        ncwa -m mask -a ${timename} -d ${timename},${t1},${t2} -v ${met} XXXX_intermediate_extract.nc XXXX_intermediate_extract2.nc
        #echo "execute ncwa to get the only non zero value"
        ncwa -y avg XXXX_intermediate_extract2.nc XXXX_intermediate_extract3.nc
        # "put the value in type_metric_value"
        metric=`ncks -v ${met} -c -H XXXX_intermediate_extract3.nc | grep ${met} | grep '=' |  awk ' {print $3}'`

        # These lines should re removed in a near future once the real pb has been diagnosed ...
        if ls XXXX_*.tmp 2>/dev/null ; then
             echo Pourquoi reste t il des fichiers ncwa.tmp
             sleep 3 ; if ls XXXX_*.tmp 2>/dev/null ; then
                           echo Et 3 secondes apres ca craint. On arrete les frais
                           exit 1
                       fi
        fi
        ;;

      #################################################################
      *) echo Metrics $metric_name not available yet. ; echo 'Want to contribute ?' ; exit 1
      #################################################################

    esac
    echo $metric
done