#!/bin/bash
#
#-----------------------------------------------------------------------------
# Requires : netcdf, ferret, nco
#-----------------------------------------------------------------------------
# Check&modify section "User choices", then run the script with ./era2gcm_uv10m.sh
#-----------------------------------------------------------------------------
# When running LMDZ with SPLA (SiMPle Aerosol model), 10m-winds on grille_T are needed to compute sea-salt emissions
# ATTENTION : for a nudged simulation, use the same reanalysis as for nudging winds (same "rea" variable here and in era2gcm.sh)
#-----------------------------------------------------------------------------
# ( A VERIFIER PLUS TARD : le regrid pourrait se faire avec CDO, mais JEscribano le deconseillait; pourquoi ? resultat identique a la precision 1e-13 :
    #Pre-requis pour ERA5, fichiers 1h au lieu de 6h:
    #cdo selhour,0/18/6 ERAfile_1h.nc ERAfile.nc 
    # general
    #cdo selvar,grille_s ./INIT/grilles_gcm.nc  ./INIT/grille_s.nc
    #cdo remapbil,INIT/grille_s.nc ERAfile.nc regridded_fileNS.nc
    ##latitudes are N->S as in ERA files, not S->N, as written by Ferret, alors on fait invertlat (ou on changera le code)
    #cdo invertlat regridded_fileNS.nc regridded_file.nc 
#

#------------------------------------------------
# User choices
#------------------------------------------------
# Periode :
mth_i=200601
mth_f=200601
#
# Choix des fichiers de guidage ("rea"nalyses) : ERA5, ERAI, OPERA
# Contact pour demander acces aux fichiers sur Jean-Zay : 
#    Sophie Bouffiès-Cloché, IPSL : Sophie.Bouffies-Cloche@ipsl.fr
rea=ERA5
#------------------------------------------------

GRID_DIR=./INIT
#-----------------------------------------------------------------------------
#Utilite du block suivant a re-examiner :
#-----------------------------------------------------------------------------
# Rajouter "rea" aux les parametres d'entree a fournir au script ?
# Ou comme pour era2gcm.sh, c'est plus facile d'editer&modifier le script, puis lancer avec ./era2gcm_uv10m.sh
while (($# > 0))
   do
   case $1 in
     "-h") cat <<........fin
           ./era2gcm_uv10m.sh [-rea ERAreanalysis] [-mthini initial_month] [-mthend final_month] [-grid_dir directory_containing_grille_gcm.nc]
........fin
     exit ;;
     "-grille_dir") GRID_DIR=$2 ; shift ; shift ;;
     "-mthini") mth_i=$2 ; shift ; shift ;;
     "-mthend") mth_f=$2 ; shift ; shift ;;
     "-rea") rea=$2 ; shift ; shift ;;
     *) $0 -h ; exit
   esac
done


tmin=1
resol=grilles_gcm.nc

GRID_FI=${GRID_DIR}/${resol}
if [ ! -f "$GRID_FI" ] ; then
   echo Le fichier $GRID_DIR/$resol est necessaire
   echo creer le fichier $GRID_FI avec grilles_gcm_netcdf.e
   exit
fi
mth=$mth_i


#####################################################################
while (( $mth <= $mth_f )) ; do
#####################################################################
   echo mth $mth
   mois=`echo $mth | cut -c 5-`
   an=`echo $mth | cut -c -4`
   ndays=( 31 28 31 30 31 30 31 31 30 31 30 31 )
   months=( jan feb mar apr may jun jul aug sep oct nov dec )
   if [ $(( $an % 4 )) = 0 ] ; then ndays[1]=29 ; fi
   imois=$(( ${mois#0} - 1 ))
   month=${months[$imois]}
   nday=${ndays[$imois]}
   tmax=$(( $nday * 4 ))
   echo PARAMETRES CALENDAIRES $imois $month $nday $tmax


   iip1=`ncdump -h $GRID_FI | grep lonu | head -1 | awk ' { print $3 } '`
   jjm=`ncdump -h $GRID_FI | grep latv | head -1 | awk ' { print $3 } '`
   (( jjp1 =  $jjm + 1 ))
   (( iim  = $iip1 - 1 ))
    \rm u10m.nc v10m.nc
#####################################################################
# Choix de la periode temporelle
#####################################################################

   t0="l=$tmin"
   t1tn="l=${tmin}:${tmax}"

#####################################################################
# Lien avec les fichiers netcdf contenant les d0 ECMWF 
#####################################################################
   if [ "$rea" = "ERA5" ] ; then
      #ANA_DIR='/gpfsstore/rech/psl/rpsl376/ERA5/NETCDF/GLOBAL_025/4xdaily'
      #suf="2e5.GLOBAL_025"
      ANA_DIR='/gpfsstore/rech/psl/rpsl376/ERA5/NETCDF/GLOBAL_025/hourly'
      suf="1e5.GLOBAL_025"
      #"Neutral wind at 10 m" :
      varu10m=u10n 
      varv10m=v10n
      fileu1h="$ANA_DIR/AN_SF/$an/u10n.$an$mois.as$suf.nc"
      filev1h="$ANA_DIR/AN_SF/$an/v10n.$an$mois.as$suf.nc"
      fileu10m="$SCRATCH/u10n.$an$mois.as$suf.nc"
      filev10m="$SCRATCH/v10n.$an$mois.as$suf.nc"
      cdo selhour,0/18/6 $fileu1h $fileu10m
      cdo selhour,0/18/6 $filev1h $filev10m
   elif [ "$rea" = "ERAI" ] ; then
      ANA_DIR='/gpfsstore/rech/psl/rpsl376/ergon/ERAI/NETCDF/GLOBAL_1125/4xdaily'
      suf="hei.GLOBAL_1125"
      varu10m=u10
      varv10m=v10
      fileu10m="$ANA_DIR/AN_SF/$an/u10.$an$mois.as$suf.nc"
      filev10m="$ANA_DIR/AN_SF/$an/v10.$an$mois.as$suf.nc"
   else  #opera
      #ANA_DIR='http://prodn.idris.fr/thredds/dodsC/ipsl_private/rpsl376/OPERA/NETCDF/4xdaily'
      ANA_DIR='/gpfsstore/rech/psl/rpsl376/ergon/OPERA/NETCDF/4xdaily'
      suf=".GLOBAL_1125"
      varu10m=u10
      varv10m=v10
      fileu10m="$ANA_DIR/AN_SF/$an/u10.$an$mois.ash$suf.nc"
      filev10m="$ANA_DIR/AN_SF/$an/v10.$an$mois.ash$suf.nc"
   fi

ls $fileu10m
ls $filev10m

outd=${rea}_uv10m/$an/$mois
mkdir -p $outd

# effacer lien vers repertoire des vents_10m ERA interpoles (utilise dans script_SIMU) s'il existe deja
# recreer lien vers le repertoire "rea" choisi
rm -f ERA10m
ln -s ${rea}_uv10m ERA10m

# utile dans "save/file..." si souci de memoire (sinon, avec ou sans $iji, meme resultat)
#ij="i=1:$iip1,j=1:$jjp1"

###################################################################3
# script ferret pour interpolation
###################################################################3

cat <<eod>| tmp10m.jnl
set memory/size=512

use "$GRID_FI"

define axis/t="1-${month}-${an}:00:00":"${nday}-${month}-${an}:18:00":6/units=hours thour

! On traite une variable et un pas de temps a la fois pour raison de memoire 
! NOTE : @asn est le seul possible pour regrid time, avec thour=axe commune pour les reanalyses,
!     dont les axes de temps ont des origines differentes (1957 pour ERAI, 1800 pour ERA5)

! --- u10m ---
use "$fileu10m"

! Pour ecrire en simple precision, comme la variable dans le fichier ERA : 
! utiliser directement $var, et non pas une variable u10m definie dans ferret (let u10m=$varu10m), qui serait en double precision
! A la fin, on change le nom de variable $var* en celui attendu par le code (voir ncrename ci-dessous)
save/clobber/file="$outd/u10m.nc" $varu10m[d=2,gxy=grille_s[d=1],$t0,gt=thour@asn]
repeat/$t1tn save/file="$outd/u10m.nc"/append $varu10m[d=2,gxy=grille_s[d=1],gt=thour@asn]

! --- v10m --- 
cancel d 2
use "$filev10m"

save/clobber/file="$outd/v10m.nc" $varv10m[d=2,gxy=grille_s[d=1],$t0,gt=thour@asn]
repeat/$t1tn save/file="$outd/v10m.nc"/append $varv10m[d=2,gxy=grille_s[d=1],gt=thour@asn]
eod

ferret -nojnl <<eod
go tmp10m.jnl
quit
eod

ncrename -v $varu10m,U10M $outd/u10m.nc
ncrename -v $varv10m,V10M $outd/v10m.nc

echo AN MTH $an $mois
(( mth = $mth + 1 ))
if [ $mois = 12 ] ; then
   (( an = $an + 1 ))
   mth=${an}01
fi

done
