ResInsight/ApplicationLibCode/ResultStatisticsCache/RigStatisticsMath.h

/////////////////////////////////////////////////////////////////////////////////
//
//  Copyright (C) 2011-2012 Statoil ASA, Ceetron AS
//
//  ResInsight is free software: you can redistribute it and/or modify
//  it under the terms of the GNU General Public License as published by
//  the Free Software Foundation, either version 3 of the License, or
//  (at your option) any later version.
//
//  ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY
//  WARRANTY; without even the implied warranty of MERCHANTABILITY or
//  FITNESS FOR A PARTICULAR PURPOSE.
//
//  See the GNU General Public License at <http://www.gnu.org/licenses/gpl.html>
//  for more details.
//
/////////////////////////////////////////////////////////////////////////////////
#pragma once

#include "RigStatisticsTools.h"

#include <cmath>
#include <cstddef>
#include <expected>
#include <set>
#include <string>
#include <vector>

class RigStatisticsMath
{
public:
    enum class PercentileStyle
    {
        REGULAR,
        SWITCHED
    };

    static void
        calculateBasicStatistics( const std::vector<double>& values, double* min, double* max, double* sum, double* range, double* mean, double* dev );

    static void calculateStatisticsCurves( const std::vector<double>& values,
                                           double*                    p10,
                                           double*                    p50,
                                           double*                    p90,
                                           double*                    mean,
                                           PercentileStyle            percentileStyle );
    static std::expected<std::vector<double>, std::string>
        calculatePercentiles( const std::vector<double>& values, const std::vector<double>& quantiles, PercentileStyle percentileStyle );

    static std::expected<std::vector<double>, std::string> calculateNearestRankPercentiles( const std::vector<double>& inputValues,
                                                                                            const std::vector<double>& percentiles,
                                                                                            PercentileStyle            percentileStyle );

    static std::expected<std::vector<double>, std::string> calculateInterpolatedPercentiles( const std::vector<double>& inputValues,
                                                                                             const std::vector<double>& percentiles,
                                                                                             PercentileStyle            percentileStyle );

    static double calculateMean( const std::vector<double>& values );
};

//==================================================================================================
/// Class to calculate a histogram, and histogram based p-value estimates
//==================================================================================================

class RigHistogramCalculator
{
public:
    enum class BinningMode
    {
        LINEAR,
        LOGARITHMIC
    };

    enum class OutOfRangeHandling
    {
        EXCLUDE,
        INCLUDE_IN_BOUNDARY_BINS
    };

    RigHistogramCalculator( double               min,
                            double               max,
                            size_t               nBins,
                            std::vector<size_t>* histogram,
                            BinningMode          binningMode        = BinningMode::LINEAR,
                            OutOfRangeHandling   outOfRangeHandling = OutOfRangeHandling::EXCLUDE );

    void addData( const std::vector<double>& data );
    void addData( const std::vector<float>& data );

    void addValue( double value );

    /// Calculates the estimated percentile from the histogram.
    /// the percentile is the domain value at which pVal of the observations are below it.
    /// The estimate is only meaningful when the histogram covers the full range of the observed values.

    double calculatePercentil( double pVal, RigStatisticsMath::PercentileStyle percentileStyle );

private:
    size_t               m_maxIndex;
    double               m_range;
    double               m_min;
    double               m_max;
    BinningMode          m_binningMode;
    OutOfRangeHandling   m_outOfRangeHandling;
    size_t               m_observationCount;
    std::vector<size_t>* m_histogram;
};

class MinMaxAccumulator
{
public:
    MinMaxAccumulator( double initMin = HUGE_VAL, double initMax = -HUGE_VAL )
        : max( initMax )
        , min( initMin )
    {
    }

    void addData( const std::vector<double>& values )
    {
        for ( double val : values )
        {
            addValue( val );
        }
    }

    void addData( const std::vector<float>& values )
    {
        for ( float val : values )
        {
            addValue( val );
        }
    }

    void addValue( double value )
    {
        if ( RigStatisticsTools::isValidNumber<double>( value ) )
        {
            if ( value < min )
            {
                min = value;
            }

            if ( value > max )
            {
                max = value;
            }
        }
    }

    double max;
    double min;
};

class PercentilAccumulator
{
public:
    PercentilAccumulator() {}

    void addData( const std::vector<double>& values )
    {
        for ( double val : values )
        {
            addValue( val );
        }
    }

    void addData( const std::vector<float>& values )
    {
        for ( float val : values )
        {
            addValue( val );
        }
    }

    void addValue( double value ) { values.push_back( value ); }

    void computep10p90( double& p10, double& p90 )
    {
        double mean = HUGE_VAL;
        double p50  = HUGE_VAL;

        RigStatisticsMath::calculateStatisticsCurves( values, &p10, &p50, &p90, &mean, RigStatisticsMath::PercentileStyle::SWITCHED );
    }

    std::vector<double> values;
};

class PosNegAccumulator
{
public:
    PosNegAccumulator( double initPos = HUGE_VAL, double initNeg = -HUGE_VAL )
        : pos( initPos )
        , neg( initNeg )
    {
    }

    void addData( const std::vector<double>& values )
    {
        for ( double val : values )
        {
            addValue( val );
        }
    }

    void addData( const std::vector<float>& values )
    {
        for ( float val : values )
        {
            addValue( val );
        }
    }

    void addValue( double value )
    {
        if ( RigStatisticsTools::isValidNumber<double>( value ) )
        {
            if ( value < pos && value > 0 )
            {
                pos = value;
            }

            if ( value > neg && value < 0 )
            {
                neg = value;
            }
        }
    }

    double pos;
    double neg;
};

class SumCountAccumulator
{
public:
    SumCountAccumulator( double initSum = 0.0, size_t initCount = 0 )
        : valueSum( initSum )
        , sampleCount( initCount )
    {
    }

    void addData( const std::vector<double>& values )
    {
        for ( double val : values )
        {
            addValue( val );
        }
    }

    void addData( const std::vector<float>& values )
    {
        for ( float val : values )
        {
            addValue( val );
        }
    }

    void addValue( double value )
    {
        if ( RigStatisticsTools::isValidNumber<double>( value ) )
        {
            valueSum += value;
            ++sampleCount;
        }
    }

    double valueSum;
    size_t sampleCount;
};

class UniqueValueAccumulator
{
public:
    UniqueValueAccumulator() {}

    void addData( const std::vector<double>& values )
    {
        for ( double val : values )
        {
            addValue( val );
        }
    }

    void addData( const std::vector<float>& values )
    {
        for ( float val : values )
        {
            addValue( val );
        }
    }

    void addValue( double value )
    {
        if ( RigStatisticsTools::isValidNumber<double>( value ) )
        {
            uniqueValues.insert( static_cast<int>( value ) );
        }
    }

    std::set<int> uniqueValues;
};