#!/bin/bash 

echo
echo This bash script gets the SSTs and SICc from
echo hourly ERA5 data and prepares the input files
echo containing daily values that should be used to 
echo create limit.nc boundary file for LMDZ.
echo 
echo The SST field contains sstk over ocean and skt over land.
echo The SIC files contain either siconc field for the 
echo period 1979-2007 or ci field for the 2008-2017 period.
echo See below !!!
echo 
echo This script was tested on spirit2 machine.
echo 
echo To execute this script you need to give the year as input 
echo ./mk_prepa_daily-SST-SIC_ERA.sh $yr 
echo
echo For exemple for 2020 year :
echo ./mk_prepa_daily-SST-SIC_ERA.sh 2020
echo
echo Ionela Musat, 12/03/2024
echo

module load ferret

local=`pwd`

yr=$1

data=ERA5
res=GLOBAL_025
ret=hourly
dtype=AN_SF
suff=as1e5

echo din is the input directory containing ERA5 daily data 
echo for $yr year
###din=/bdd/ERA5/NETCDF/GLOBAL_025/hourly/AN_SF/${yr}
din=/bdd/${data}/NETCDF/${res}/${ret}/${dtype}/${yr}

MONTHS="01 02 03 04 05 06 07 08 09 10 11 12"

for mth in $MONTHS; do

if [ "${mth}" == "02" ]; then
### ndays=`mk_an_bissextile.sh ${yr}`

echo
echo Check if the year is a leap one or not
echo
echo A leap year is a year which is a multiple of 4,
echo except for years evenly divisible by 100, but not by 400.
echo
let " yrs4 = ${yr} / 4"
let " yrs100 = ${yr} / 100"
let " yrs400 = ${yr} / 400"

let "yrd4 = ${yrs4} *4 "
let "yrd100 = ${yrs100} *100 "
let "yrd400 = ${yrs400} *400 "

ndays=28
if [ ${yr} == ${yrd400} ]; then
ndays=29
elif [ ${yr} == ${yrd100} ]; then
ndays=28
elif [ ${yr} == ${yrd4} ]; then
ndays=29
fi

if [ "${ndays}" == "29" ]; then
 echo $yr is a leap year
else
 echo $yr is a noleap year
fi

fi

case $mth in 
             01|03|05|07|08|10|12) ndays=31 ;;
             04|06|09|11) ndays=30 ;;
esac

var=skt
cp ${din}/${var}.${yr}${mth}.${suff}.${res}.nc .
cdo daymean ${var}.${yr}${mth}.${suff}.${res}.nc ${var}_daily.${yr}${mth}.${res}.nc
chmod u+w ${var}.${yr}${mth}.${suff}.${res}.nc
\rm ${var}.${yr}${mth}.${suff}.${res}.nc 

var=sstk
cp ${din}/${yr}/${var}.${yr}${mth}.${suff}.${res}.nc .
cdo daymean ${var}.${yr}${mth}.${suff}.${res}.nc ${var}_daily.${yr}${mth}.${res}.nc
chmod u+w ${var}.${yr}${mth}.${suff}.${res}.nc
\rm ${var}.${yr}${mth}.${suff}.${res}.nc 

if [ -f tmp.ps ]; then
\rm tmp.ps
fi

cat<<eof>go_sstk_stk.jnl
cancel data/all
cancel view

use "skt_daily.$1.$3.nc"
use "sstk_daily.$1.$3.nc"

let sst = if sstk[d=2] then sstk[d=2] else skt[d=1]
shade sst[l=1]-skt[d=1,l=1]; go land

save/file=sst_$4.$1.nc sst[l=1]
repeat/l=2:$2 (save/append/file=sst_$4.$1.nc sst[l=`l`])
quit

eof

ferret -batch tmp.ps -script go_sstk_stk.jnl ${yr}${mth} ${ndays} ${res} ${data}

echo mth $mth ndays ${ndays}

done

varnew=sst
cdo cat ${varnew}_${data}.${yr}*.nc ${varnew}_daily_${data}.${yr}.nc

echo
echo Change dimensions and SST variable names and add the K unit.
echo This are mandatory for using the SSTs/SICs files as input files to the
echo create_eta0_limit.F90 main code.
echo 
var=sst
varM=SST
varnew=tosbcs
fin=${var}_daily_${data}.${yr}.nc
fout=${varnew}_daily_${data}.${yr}.nc
ncrename -v LONGITUDE,longitude -d LONGITUDE,longitude ${fin}
ncrename -v LATITUDE,latitude -d LATITUDE,latitude ${fin}
ncrename -v TIME,time -d TIME,time ${fin}
ncrename -v ${varM},${varnew} ${fin}
ncatted -a units,${varnew},a,c,"K" ${fin}
mv ${fin} ${fout}

echo 
echo Get SICs field : files named "ci" are containing the siconc or ci variable 
echo corresponding to sea-ice in 0-1 unit
echo
var=ci
for mth in $MONTHS; do

cp ${din}/${var}.${yr}${mth}.${suff}.${res}.nc .
cdo daymean ${var}.${yr}${mth}.${suff}.${res}.nc ${var}_daily.${yr}${mth}.${res}.nc
chmod u+w ${var}.${yr}${mth}.${suff}.${res}.nc
\rm ${var}.${yr}${mth}.${suff}.${res}.nc

done

cdo cat ${var}_daily.${yr}*.${res}.nc ${var}_daily_${data}.${yr}.nc

echo
echo Rename siconc or ci variable as sicbcs, multiply it by 100 and change unit
echo
if [ $yr -ge 2008 ] && [ $yr -le 2017 ]; then
varname=ci
else
varname=siconc
fi

fin=${var}_daily_${data}.${yr}.nc
varnew=sicbcs
mulc=100.
fout=${varnew}_daily_${data}.${yr}.nc
cdo mulc,$mulc ${fin} ${fout}
ncrename -v ${varname},${varnew} ${fout}
ncatted -a units,${varnew},o,c,"%" ${fout}

\rm ci*nc sst*nc skt*nc

exit
