#!/bin/bash
export LC_ALL=C
ipipe=0
echo //////////////////////////////////////////////////
echo           Going to analyse a field
echo //////////////////////////////////////////////////
divacalclog='./output/divacalc.log'
echo ////////////////////////////////////////////////// > $divacalclog
echo           Going to analyse a field >> $divacalclog
echo ////////////////////////////////////////////////// >> $divacalclog
#
divarunning=$(ps -eafl | grep -i 'diva\.a' | wc -l)
if [ "$divarunning" != "0" ]
then
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! >> $divacalclog
echo ????????diva already running?????? >> $divacalclog
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! >> $divacalclog
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
echo ????????diva already running??????
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ps -eafl | grep -i 'diva\.a'
fi
iqc=0
if [ "$1" == "qc" ] 
then
echo =============== >> $divacalclog
echo Quality control >> $divacalclog
echo =============== >> $divacalclog
iqc=1
fi

ibii=0
if [ "$1" == "Bii" ] 
then

echo =============== >> $divacalclog
echo Bii calculation >> $divacalclog
echo =============== >> $divacalclog
ibii=1
ipipe=0
fi


Filepar=./input/param.par
{
read linecomment
read lc
read linecomment
read icoordchange
read linecomment
read ispec
read linecomment
read ireg
read linecomment
read xori
read linecomment
read yori
read linecomment
read dx
read linecomment
read dy
read linecomment
read nx
read linecomment
read ny
read linecomment
read valex
read linecomment
read snr
read linecomment
read varbak
} < $Filepar


if [ "$ibii" == "1" ]
then
ispec=-1000
fi


if [ -f ./input/valatxy.coord ] 
then
cp ./input/valatxy.coord ./divawork/fort.79
fi
xi=1
cd ./meshgenwork
File23=./fort.23

cp  ../output/meshvisu/fort.10 .
cp  ../output/meshvisu/fort.22 .
cp  ../output/meshvisu/fort.23 .

{
read nnt1
read nnint
read nelt
} < $File23
 
cp fort.22 ../divawork/fort.11
cd ../divawork
# try to estimate memory allocation....
# 
NREA=0
NINT=0
NDATA=$(cat ../input/data.dat | wc -l)
echo 'Number of data points:' ${NDATA} >> ../$divacalclog
echo ' '
let NREA=$NREA+4*$NDATA
let NREA=$NREA+${nx}*${ny}
let NINT=$NINT+8*$NDATA
let NINT=$NINT+$nelt
let NINT=$NINT+41
let NDDLT=3*$nnt1+$nnint
let NDDLE=16
let NNT=$nnt1+$nnint
let NINT=$NINT+2*$NNT
let NINT=$NINT+$NDDLT
let NINT=$NINT+16*$nelt
let NINT=$NINT+2*$NDDLE
let NTERM=$NDDLT*$NDDLT
# instead of 90, square root of $NDDLT
let NDDLTB=$NDDLT
#echo Square root of $NDDLTB
#sqrt() { local x=$1 s=$1 os=0;while ((s!=os));\
#do os=$s;s=$(((s+x/s)/2)); done;echo $s; }
##sqrt() { local x=$1 s=$1 os=0; \
##while [ "$((os-s))" -ge "1" -o "$((s-os))" -ge "1" ]; \
##do os=$s;s=$(((s+x/s)/2)); done;echo $s; }
##xx=$(sqrt $NDDLTB)
xx=$(echo $NDDLTB | awk '{print int(sqrt($1))}')
#echo root is $xx $xxb
let NTERM=$NDDLT*$xx
let NTERM=$NTERM+$NTERM/2
let NREA=$NREA+$nelt
let NREA=$NREA+2*$NTERM
let NREA=$NREA+2*$NDDLT
let NREA=$NREA+$NDDLE*$NDDLE
let NREA=$NREA+$NDDLE*$NDDLE
let NREA=$NREA+$NDDLE
let NREA=$NREA+38*$nelt
let NREA=$NREA+2*$nnt1
# For optimisation of element location, rought estimate:
let NINT=$NINT+$nelt*100
#
itcs=0
if [ -f ../input/constraint.dat ]
then
let itcs=$itcs+1

echo Advection requested >> ../$divacalclog

cat ../input/constraint.dat > ./fort.49
cp ../input/UVinfo.dat ./fort.93
cp ../input/Uvel.dat ./fort.91
cp ../input/Vvel.dat ./fort.92

#JMB2012 added support for source terms in advection constraint
# as many sources as whished, format: x y Q
if [ -f ../input/sources.dat ]
then
cp ../input/sources.dat ./fort.24
fi


let NREA=$NREA+2*$nnt1
{
read aaaaa
read bbbbb
read ccccc
read ddddd
read nxuuu
read nyuuu
} < ../input/UVinfo.dat
let NREA=$NREA+$nxuuu*$nyuuu
let NREA=$NREA+2*$nnt1
let NREA=$NREA+1
fi
if [ -f ../input/RL.dat ]
then
let itcs=$itcs+2
echo Variable L requested >> ../$divacalclog
cp ../input/RL.dat ./fort.94
cp ../input/RLinfo.dat ./fort.95
let NREA=$NREA+$nnt1
{
read aaaaa
read bbbbb
read ccccc
read ddddd
read nxlll
read nylll
} < ../input/RLinfo.dat
let NREA=$NREA+$nxlll*$nylll
let NREA=$NREA+$nnt1
let NREA=$NREA+1
fi

#creation of file fort.13
#
# 'griddef.a' replaced by awk command 
#------------------------------------
#echo $xori > ./fort.13
#echo $yori >> ./fort.13
#echo $dx >> ./fort.13
#echo $dy >> ./fort.13
#echo $nx >> ./fort.13
#echo $ny >> ./fort.13
#echo $valex >> ./fort.13
#../../bin/griddef.a >> diva.log

echo $(echo $xori $dx | awk '{print ($1-$2)}') >  ./fort.13
echo $(echo $yori $dy | awk '{print ($1-$2)}') >> ./fort.13
echo $dx >> ./fort.13
echo $dy >> ./fort.13
echo $nx >> ./fort.13
echo $ny >> ./fort.13
echo $valex >> ./fort.13


echo $(echo $lc | awk '{print (1/($1*$1*$1*$1))}') > fort.12
echo $(echo $lc $xi | awk '{print (2*$2/($1*$1))}') >> fort.12

if [ "$ireg" == "3" ]
then
aaaref=$(head -1 fort.12)
echo 0 > newalpha
# Test, not really zero, but small

##echo  1.E-6 $aaaref | ../../bin/multiply.a > newalpha
echo $(echo $aaaref | awk '{print ($1*1.E-6)}') >  newalpha

tail -1 fort.12 >> newalpha
mv newalpha fort.12
let ireg=1
fi
echo $varbak > ./fort.15
cd ..
if [ "$iqc" != "1" ]
then
if [ "$ispec" -lt "0" ] 
then
if [ "$varbak" != "0" ]
then
if [ "$1" != "-pfsum" ]
then
#dvkern
echo %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% >> $divacalclog
echo %  Going to call another diva execution    % >> $divacalclog
echo %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% >> $divacalclog


echo divacalc: ////////////////////////////////////////// >> ./output/diva.log
echo divacalc:         running dvcovariance >> ./output/diva.log
echo divacalc: ////////////////////////////////////////// >> ./output/diva.log
echo ' '  >> ./output/diva.log

. dvcovariance
    cat $divacalclog ./output/dvcovariance.log > bidon
    mv bidon $divacalclog

if [ "$ibii" == "1" ]
then
cd divawork
#pwd
../../bin/bintoasc.a | awk '{print $1,$2,$3/(1-$3)}' >../output/Bii.dat
cd ..
echo Finished bii  >> $divacalclog
exit
fi

if [ "$ipipe" == "1" ]
then
echo %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% >> $divacalclog
echo %  Have launched another diva execution    % >> $divacalclog
echo %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% >> $divacalclog
else
echo %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% >> $divacalclog
echo %  Have finished another diva execution    % >> $divacalclog
echo %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% >> $divacalclog
fi
fi
fi
fi
fi
if [ -f ./input/data.dat ] 
then 
cp ./input/data.dat ./divawork/fort.44
else
echo Need to provide data.dat in ./input !
##ls -l ./input
##pwd
fi
cd ./divawork
echo coord 1 > ./fort.10
echo $icoordchange >> ./fort.10
echo mathpr 1 >> ./fort.10
echo 2 >> ./fort.10
echo 0 >> ./fort.10
echo 2 >> ./fort.10
echo topolo 1 >> ./fort.10
echo $nnt1 >> ./fort.10
echo $nnint >> ./fort.10
echo $nelt >> ./fort.10
echo datapr 1 >> ./fort.10
echo $ireg >> ./fort.10
#JMB
if [ "$itcs" != "0" ] 
then 
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! >> ../$divacalclog
echo  Module constraint activated >> ../$divacalclog
echo $itcs
echo constr 1 >> ./fort.10
echo $itcs >> ./fort.10
fi
echo solver 1 >> ./fort.10
echo 0 >> ./fort.10
echo stores 1 >> ./fort.10
echo 3 >> ./fort.10
echo VARBAK: $varbak >> ../$divacalclog
if [ "$iqc" != "1" ]
then
if [ "$varbak" != "0" ]
then
if [ "$ispec" != "0" ]
then
	echo ' ' >> ../$divacalclog
echo Errors will be calculated >> ../$divacalclog
echo esterr 1 >> ./fort.10
echo $ispec >> ./fort.10
fi
fi
#echo $ >> ./fort.10
echo gcvfac 1 >> ./fort.10
echo 5 >> ./fort.10
echo stopex >> ./fort.10
else
#QC
echo dataqc 1 >> ./fort.10
echo stopex >>./fort.10
fi
#creation of file fort.20
nbccol=$(head -n 1 ./fort.44 | wc -w)
#echo $lc $snr $nbccol $xi | ../../bin/calcmu.a >> diva.log

# New version
multconst=$(echo $lc $snr $xi | ../../bin/calcmuc.a)
multvar=$(echo $multconst | awk '{print 1*$1}')

#echo multiply weight $multvar $multconst so what
awk -v multconst=$multvar '{ if (NF>2) {if ($4 > "") {$4=$4*multconst} else {$4=multconst}}; {print $1,$2,$3,$4}}' fort.44 > fort.20

##?? if [ $? -ne 0 ];then
##?? echo ' '
##?? echo --------------------------------------------
##?? echo A problem was encountered during execution !
##?? echo       divacalc:  calmu.a
##?? echo Check execution track
##?? echo --------------------------------------------

NDATAB=$(cat ./fort.20 | wc -l)
if [ "$NDATAB" -lt "$NDATA" ]
then
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! >> ../$divacalclog
echo ------------------------------------------- >> ../$divacalclog
echo ??????????????????????????????????????????? >> ../$divacalclog
echo ?????   Corrupted data file    ???????????? >> ../$divacalclog
echo From $NDATA lines in data.dat >> $divacalclog
echo Only $NDATAB have been feed into DIVA >> $divacalclog
echo ??????????????????????????????????????????? >> ../$divacalclog
echo ------------------------------------------- >> ../$divacalclog
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! >> ../$divacalclog
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
echo -------------------------------------------
echo ???????????????????????????????????????????
echo ?????   Corrupted data file    ????????????
echo From $NDATA lines in data.dat
echo Only $NDATAB have been feed into DIVA
echo ???????????????????????????????????????????
echo -------------------------------------------
echo !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
fi
if [ "$1" == "-pfsum" ]
then
echo Summing covariances >> ../$divacalclog
cd ..

echo divacalc: ////////////////////////////////////////// >> ./output/diva.log
echo divacalc:         running dvcovariancesum >> ./output/diva.log
echo divacalc: ////////////////////////////////////////// >> ./output/diva.log
echo ' '  >> ./output/diva.log
. ./dvcovariancesum
    cat $divacalclog ./output/dvcovariancesum.log > bidon
    mv bidon $divacalclog
exit
fi

echo divacalc: ////////////////////////////////////////// >> ../output/diva.log
echo divacalc:'        Diva calculation (diva.a)'  >> ../output/diva.log
echo divacalc: ////////////////////////////////////////// >> ../output/diva.log
echo ' '  >> ../output/diva.log

##echo 'REMOVE FORT.35 and FORT.43'
rm -f fort.43
rm -f fort.35
#execute diva
echo $NINT $NREA | ../../bin/diva.a >> ../output/diva.log

if [ $? -ne 0 ];then
echo ' '
echo --------------------------------------------
echo A problem was encountered during execution !
echo       divacalc:  diva.a
echo 'Check execution track in diva.log'
echo --------------------------------------------
echo ' ' >> ../$divacalclog
echo -------------------------------------------- >> ../$divacalclog
echo A problem was encountered during execution ! >> ../$divacalclog
echo       divacalc:  diva.a >> ../$divacalclog
echo 'Check execution track in diva.log' >> ../$divacalclog
echo -------------------------------------------- >> ../$divacalclog
echo ' ' >> ../output/diva.log
echo -------------------------------------------- >> ../output/diva.log
echo A problem was encountered during execution ! >> ../output/diva.log
echo       divacalc:   diva.a  >> ../output/diva.log
echo 'Check execution track in diva.log' >> ../output/diva.log
echo -------------------------------------------- >> ../output/diva.log
fi

dataonmesh=0
##ctroupin Oct 2012
# Check if fort.84 was generated AND if not empty...
if [ ! -s ./fort.84 ];
then
        echo ' ' >> ../$divacalclog
        echo '//////////////////////////////////////////////////' >> ../$divacalclog
        echo 'Problem in diva.a : analysis could not be finished' >> ../$divacalclog
        echo ' ' >> ../$divacalclog
        echo 'Check diva.log' >> ../$divacalclog
        echo '//////////////////////////////////////////////////' >> ../$divacalclog

        echo ' '
        echo '//////////////////////////////////////////////////'
        echo 'Problem in diva.a : analysis could not be finished'
        echo ' '
        echo 'Check diva.log'
        echo '//////////////////////////////////////////////////'
fi
##---




#if [ -f ./fort.35  ]; then
#dataonmesh=$(cat ./fort.35 | awk '{print $1-$2}')
#echo Data are outside the mesh, will produce empty files
#fi
#if [ "$dataonmesh" -lt "1" ]; then
#cd ..
#dvnil
#cd ./divawork
#else
#if [ -f ./fort.43 ] 
#then
#echo No data, will produce empty files
#echo $nx $ny $valex | ../../bin/ghervalex.a
#cp fort.84 fort.87
#fi
#fi

# New version 11 Jan 2012
#------------------------

if [ -f ./fort.35  ]; then
	dataonmesh=$(cat ./fort.35 | awk '{print $1-$2}')
fi
if [ "$dataonmesh" -lt "1" ]; then
	cd ..
	echo Data are outside the mesh, will produce empty files
	dvnil
	cd ./divawork
fi

rm -f fort.59
rm -f fort.61
if [ "$iqc" == "1" ] 
then
#
#output of results to user
cp fort.76 ../output/expectederroratdatapoint.anl
cp fort.71 ../output/fieldatdatapoint.anl
echo ' '
echo divacalc: Looking for outliers > ../$divacalclog
echo ' '


echo divacalc: ////////////////////////////////////////// >> ../output/diva.log
echo divacalc:       lookforoutliers   >> ../output/diva.log
echo divacalc: ////////////////////////////////////////// >> ../output/diva.log
echo ' '  >> ../output/diva.log

echo $valex | ../../bin/lookforoutliers.a >> ../output/diva.log
if [ $? -ne 0 ];then
echo ' '
echo --------------------------------------------
echo A problem was encountered during execution !
echo       divacalc:   lookforoutliers.a
echo 'Check execution track in diva.log'
echo --------------------------------------------
echo ' ' >> ../output/diva.log
echo -------------------------------------------- >> ../output/diva.log
echo A problem was encountered during execution ! >> ../output/diva.log
echo       divacalc:    lookforoutliers.a  >> ../output/diva.log
echo 'Check execution track in diva.log' >> ../output/diva.log
echo -------------------------------------------- >> ../output/diva.log
fi
echo ' '

echo ' ' >> ../$divacalclog
echo Copying outlier lists >> ../$divacalclog
echo ' '
if [ -f fort.66 ]
then
mv -v fort.66 ../output/outliers.dat
mv -v fort.67 ../output/outliers.normalized.dat
else
echo ' ' >> ../$divacalclog
echo -------------------------------------------- >> ../$divacalclog
echo A problem was encountered during execution ! >> ../$divacalclog
echo 'Check execution track in diva.log' >> ../$divacalclog
echo -------------------------------------------- >> ../$divacalclog
echo ' '
echo --------------------------------------------
echo A problem was encountered during execution !
echo 'Check execution track in diva.log'
echo --------------------------------------------
fi
else
#output as netcdf not yet done
echo ' ' >> ../$divacalclog
echo Output of results copied in ./output/ >> ../$divacalclog
echo ' ' >> ../$divacalclog
cp fort.84 ../output/fieldgher.anl
if [ -f fort.82 ]
then
cp fort.82 ../output/valatxyascii.anl
fi
cp fort.83 ../output/fieldascii.anl
if [ -f fort.87 ]
then
cp fort.87 ../output/errorfieldgher.anl
fi
if [ -f fort.73 ]
then
cp fort.73 ../output/erroratxyascii.anl
fi
if [ -f fort.86 ]
then
cp fort.86 ../output/errorfieldascii.anl
fi
cp fort.71 ../output/fieldatdatapoint.anl
if [ -f fort.72 ] 
then
cp fort.72 ../output/erroratdatapoint.anl
fi
if [ -f fort.77 ]
then
cp fort.77 ../output/gcvval.dat
fi

##echo ' '
##echo Creation of file GridInfo.dat
##echo ' '

echo $xori > ../output/ghertonetcdf/GridInfo.dat
echo $yori >> ../output/ghertonetcdf/GridInfo.dat
echo $dx >> ../output/ghertonetcdf/GridInfo.dat
echo $dy >> ../output/ghertonetcdf/GridInfo.dat
echo $nx >> ../output/ghertonetcdf/GridInfo.dat
echo $ny >> ../output/ghertonetcdf/GridInfo.dat
echo $valex >> ../output/ghertonetcdf/GridInfo.dat


#Iceedge case
if [ -f ../input/iceedge.dat ]
then
echo Ice edge to be used as mask
cd ..
dvicemask
cp ./output/fieldgher.anl ./divawork/fort.84
cd divawork
echo finished ice edge case
fi
cp fort.84 ../output/ghertonetcdf
if [ -f fort.87 ]
then 
cp  fort.87 ../output/ghertonetcdf
fi
iecal=0
if [ "$ispec" == "1" ]
then
iecal=1
fi
if [ "$ispec" == "3" ]
then
iecal=1
fi
if [ "$ispec" == "5" ]
then
iecal=1
fi
if [ "$ispec" == "7" ]
then
iecal=1
fi
if [ "$ispec" == "11" ]
then
iecal=1
fi
if [ "$ispec" == "13" ]
then
iecal=1
fi
if [ "$ispec" == "15" ]
then
iecal=1
fi
if [ "$ispec" == "17" ]
then
iecal=1
fi
if [ "$ispec" == "-1" ]
then
iecal=1
fi
if [ "$ispec" == "-3" ]
then
iecal=1
fi
if [ "$ispec" == "-5" ]
then
iecal=1
fi
if [ "$ispec" == "-7" ]
then
iecal=1
fi
if [ "$ispec" == "-11" ]
then
iecal=1
fi
if [ "$ispec" == "-13" ]
then
iecal=1
fi
if [ "$ispec" == "-15" ]
then
iecal=1
fi
if [ "$ispec" == "-17" ]
then
iecal=1
fi
# JMB ADAPTED FOR DIFFERENT CASES IF ISPEC
if [ "$varbak" != "0" ]
then
if [ "$iecal" == "1" ] 
then
echo Creating netcdf file for field and associated error >> ../$divacalclog
echo in directory ./output/ghertonetcdf/ >> ../$divacalclog
cd ../output/ghertonetcdf

echo divacalc: ////////////////////////////////////////// >> ../diva.log
echo divacalc:       netcdfoutput   >> ../diva.log
echo divacalc: ////////////////////////////////////////// >> ../diva.log
echo ' '  >> ../diva.log

../../../bin/netcdfoutput.a >> ../diva.log
if [ $? -ne 0 ];then
echo ' '
echo --------------------------------------------
echo A problem was encountered during execution !
echo       divacalc:   netcdfoutput.a
echo 'Check execution track in diva.log'
echo --------------------------------------------
echo ' ' >> ../$divacalclog
echo -------------------------------------------- >> ../$divacalclog
echo A problem was encountered during execution ! >> ../$divacalclog
echo       divacalc:   netcdfoutput.a >> ../$divacalclog
echo 'Check execution track in diva.log' >> ../$divacalclog
echo -------------------------------------------- >> ../$divacalclog
echo ' ' >> ../diva.log
echo -------------------------------------------- >> ../diva.log
echo A problem was encountered during execution ! >> ../diva.log
echo       divacalc:   netcdfoutput.a  >> ../diva.log
echo 'Check execution track in diva.log' >> ../diva.log
echo -------------------------------------------- >> ../diva.log
fi
fi
fi
if [ "$varbak" == "0"  -o  "$iecal" == "0" ] 
then
echo Creating netcdf file only for field  >> ../$divacalclog
echo since Varbak and ispec are $varbak $ispec >> ../$divacalclog

cd ../output/ghertonetcdf
../../../bin/netcdfoutputfield.a >> ../diva.log
if [ $? -ne 0 ];then
echo ' '
echo --------------------------------------------
echo A problem was encountered during execution !
echo       divacalc:   netcdfoutputfield.a
echo 'Check execution track in diva.log'
echo --------------------------------------------
echo ' ' >> ../$divacalclog
echo -------------------------------------------- >> ../../$divacalclog
echo A problem was encountered during execution ! >> ../../$divacalclog
echo       divacalc:   netcdfoutputfield.a >> ../../$divacalclog
echo 'Check execution track in diva.log' >> ../../$divacalclog
echo -------------------------------------------- >> ../../$divacalclog
echo ' ' >> ../diva.log
echo -------------------------------------------- >> ../diva.log
echo A problem was encountered during execution ! >> ../diva.log
echo       divacalc:    netcdfoutputfield.a  >> ../diva.log
echo 'Check execution track in diva.log' >> ../diva.log
echo -------------------------------------------- >> ../diva.log
fi
###rm -f ../error*.anl
fi
cd ../..

echo ' ' >> $divacalclog
echo -------------------- >> $divacalclog
echo Analysis is finished >> $divacalclog
echo -------------------- >> $divacalclog
echo ' '
echo --------------------
echo Analysis is finished
echo --------------------
fi
if [ "$ispec" -lt "0" ] 
then
if [ "$varbak" != "0" ]
then
rm -f ./divawork/sub/fort.*
fi
fi
