#!/bin/bash

# This script creates netcdf files used by sparring to drive LPJG offline from raw IFS output.
#
# Created By Klaus Wyser, some modifications by Etienne Tourigny

set -xuve

# ==============================================================================
# edit these variables

expname=a0l3; npts=88838
#expname=AHIL; npts=35718

y1=1990
y2=1990
y0=1990

srcdir=/gpfs/scratch/bsc32/bsc32051/lpjg/a0l3/in
outdir=/gpfs/scratch/bsc32/bsc32051/lpjg/a0l3/out

#tmpdir=$SNIC_TMP
tmpdir=/gpfs/scratch/bsc32/bsc32051/lpjg/tmp
#nprocs=$SLURM_NPROCS
nprocs=1
# ==============================================================================

mkdir -p $outdir
mkdir -p $tmpdir
cd $tmpdir

# create grib_filter script

# add 169 (ssrd) and 175 (strd) for surface downward radiation
# LPJG currently uses surface net radiation
sstring=
for c in 39 40 41 42 139 167 170 183 236; do
    sstring+=" || marsParam is \"$c.128\""
done
sstring=${sstring:4:${#sstring}}

fstring=
for c in 142 143 144 176 177; do
    fstring+=" || marsParam is \"$c.128\""
done
fstring=${fstring:4:${#fstring}}

cat > gf << EOT
if ( $sstring ) { write "x_state" ; }
if ( $fstring ) { write "x_flux" ; }
EOT


# first timestep
ym=$((y1-1))12
if (( y1 == $y0 ))
then
    grib_filter gf $srcdir/output/ifs/001/ICMGG${expname}+000000
else
    grib_filter gf $srcdir/output/ifs/$(printf %03d $((y1-y0)))/ICMGG${expname}+$ym
fi
mv x_state x_state_$ym
cdo shifttime,"-6hour" x_flux x_flux_$ym
rm -f x_flux

njobs=0
for y in $(seq $y1 $y2); do
    for m in {1..12}; do
        (  (( njobs+=1 ))
           ym=$((100*y+m))

           mkdir -p $tmpdir/data_$ym
           cd $tmpdir/data_$ym
           grib_filter ../gf \
              $srcdir/output/ifs/$(printf %03d $((y-$((y0-1)))))/ICMGG${expname}+$ym
           mv x_state ../x_state_$ym
           cdo shifttime,"-6hour" x_flux ../x_flux_$ym
           cd -
           rm -r $tmpdir/data_$ym ) &
        if (( njobs%nprocs == 0 )); then wait; fi
    done
done
wait

njobs=0
for y in $(seq $y1 $y2); do
    ( (( njobs+=1 ))
      cdo cat x_state_$((y-1))12 x_state_$y* y_state_$y
      cdo splitcode -daymean -selyear,$y y_state_$y z_${y}_

      cdo cat x_flux_$y* y_flux_$y
      cdo splitcode -daymean -selyear,$y y_flux_$y* z_${y}_ ) &
    if (( njobs%nprocs ==0 )); then wait; fi
done
wait

for f in $(/bin/ls z_*grb ); do
    fnew=$(basename $f .grb | cut -d _ -f 3 )
    year=$(basename $f .grb | cut -d _ -f 2 )

    cdo -f nc4c -z zip setgrid,g${npts}x1 $f $outdir/var${fnew}_${year}_dayavg.nc
done

exit 0