good_points

#!/usr/bin/env bash

if [ $# -lt 1 ]
then
clear
echo ""
echo " +----------------------------------------------------+"
echo " |   This program generates a file with 'good' data   |"
echo " |       based on values from specified column        |"
echo " |                                                    |"
echo " |  * Version 2017-02-23                              |"
echo " |  * Licensed under the MIT license:                 |"
echo " |    http://opensource.org/licenses/MIT              |"
echo " |  * Copyright (c) 2017 Przemysław Bruś              |"
echo " +----------------------------------------------------+"
echo " |   Program calculates an avarage x and a standard   |"
echo " |deviation d. Then it rejects lines containing values|"
echo " | differ by +/-A*d from x. It repeats this procedure |"
echo " |    N times. An output file is named as in.good     |"
echo " +----------------------------------------------------+"
echo " | 1. in.ext - input file                             |"
echo " | 2. C - column with data                            |"
echo " | 3. A - coefficient of standard deviation           |"
echo " | 4. N - number of iterations                        |"
echo " | 5. M - min. amount of points to calculate d value  |"
echo " |                                                    |"
echo " | Lines starting with # and empty lines are ignored  |"
echo " | Default values: C = 1, A = 2.5, N = 4, M = 30      |"
echo " +----------------------------------------------------+"
echo " |                      Usage:                        |"
echo " |./good_points in.ext C A N M                        |"
echo " |./good_points in.ext 1 3 5 15                       |"
echo " +----------------------------------------------------+"
echo ""
exit 0
fi

data_file="$1"
fext=`echo ${data_file} | cut -d '.' -f2`
fname=`basename ${data_file} .${fext}`
text="tmp"
awk '{if ($0 !~ /^#/ && NF) print $0}' ${data_file} > ${fname}.${text}

# unique name of temporary file prevents to
# overwrite files in the current directory
i=1
tmp=${fname}_${i}.${text}
while [ -e ${tmp} ]
do
    i=$((i+1))
    tmp=${fname}_${i}.${text}
done

if [ $# -le 5 ]
then
    C=1
    A=2.5
    N=4
    M=30
else
    clear
    echo ""
    echo "  Too many arguments!"
    echo ""
    exit
fi

if [ $# -eq 5 ]
then
    C=$2
    A=$3
    N=$4
    M=$5
elif [ $# -eq 4 ]
then
    C=$2
    A=$3
    N=$4
elif [ $# -eq 3 ]
then
    C=$2
    A=$3
elif [ $# -eq 2 ]
then
    C=$2
fi

# c4 depends on M value
array=(
    0.79788 0.88623 0.92132
    0.93999 0.95153 0.95937
    0.96503 0.96931 0.97266
    0.97535 0.97756 0.97941
    0.98097 0.98232 0.98348
    0.98451 0.98541 0.98621
    0.98693 0.98758 0.98817
    0.9887  0.98919 0.98964
    0.99005 0.99043 0.99079
    0.99111 0.99142 0.9917
)

coeff_c4()
{
    idx=`wc -l ${fname}.${text} | awk '{printf("%d", $1-2)}'`
    c4=1

    if [ $idx -ge 0 ] && [ $idx -lt ${#array[@]} ]
    then
        c4=${array[$idx]}
    fi
}

coeff_c4_it()
{
    idx=`wc -l ${fname}.${text}.it-${j} | awk '{printf("%d", $1-2)}'`
    c4=1

    if [ $idx -ge 0 ] && [ $idx -lt ${#array[@]} ]
    then
         c4=${array[$idx]}
    fi
}

state=1


    for (( i=1; i<=$N; i++ ))
    do
        j=$((i-1))
        if [ $i -eq 1 ]
        then
            coeff_c4
            x=`awk -v col=$C '{sum += $col} END {print sum/NR}' ${fname}.${text}`
            d=`awk -v col=$C -v ave=$x -v c4=$c4 '{sum += ($col - ave)*($col - ave)} END {print sqrt(sum/(NR-1))/c4}' ${fname}.${text}`
            awk -v col=$C -v ave=$x -v a=$A -v sig=$d '{if (($col > ave-a*sig) && ($col < ave+a*sig)) print $0}' ${fname}.${text} > ${tmp}
            oldline=`wc -l ${fname}.${text} | awk '{print $1}'`
        else
            coeff_c4_it
            x=`awk -v col=$C '{sum += $col} END {print sum/NR}' ${fname}.${text}.it-${j}`
            d=`awk -v col=$C -v ave=$x -v c4=$c4 '{sum += ($col - ave)*($col - ave)} END {print sqrt(sum/(NR-1))/c4}' ${fname}.${text}.it-${j}`
            awk -v col=$C -v ave=$x -v a=$A -v sig=$d '{if (($col > ave-a*sig) && ($col < ave+a*sig)) print $0}' ${fname}.${text}.it-${j} > ${tmp}
            oldline=`wc -l ${fname}.${text}.it-${j} | awk '{print $1}'`
        fi

        howmany=`wc -l ${tmp} | awk '{print $1}'`

        if [ $howmany -lt $M ]
        then
            echo "  Too little data. Goodbye"
            rm ${tmp}
            i=$j
            break
        elif [ $state -eq 1 ]
        then
            state=0
            echo ""
            echo "   iteration   average   stand.deviation   rejected"
            echo ""
        fi

        newline=`wc -l ${tmp} | awk '{print $1}'`
        mv ${tmp} ${fname}.${text}.it-${i}
        reject=$((oldline - newline))
        awk -v i=$i -v x=$x -v d=$d -v re=$reject 'BEGIN {printf("%8d %13lf %13lf %12d\n", i, x, d, re)}'
    done


if [ $howmany -ge $M ]
then
    j=$((j + 1))
fi

if [ -e ${fname}.${text}.it-${j} ]
then
    cp ${fname}.${text}.it-${j} ${fname}.good
else
    cp ${fname}.${text} ${fname}.good
fi
rm -f *.it-* ${fname}.${text}