RollingStat.h

/*
 * RollingStat.h
 *
 * This file is part of Mozzi.
 *
 * Copyright 2013-2024 Tim Barrass and the Mozzi Team
 *
 * Mozzi is licensed under the GNU Lesser General Public Licence (LGPL) Version 2.1 or later.
 *
 * WARNING: this code is incomplete, it doesn't work yet
 */
 
#ifndef ROLLINGSTAT_H
#define ROLLINGSTAT_H
 
#include "RollingAverage.h"
#include "mozzi_fixmath.h"
 
// adapted from RunningStat, http://www.johndcook.com/standard_deviation.html
/** @ingroup sensortools
WARNING: this class is work in progress, don't use it yet.
Calculates an approximation of the variance and standard deviation for a window of recent inputs.
@tparam T the type of numbers to use.  Choose unsigned int, int , uint8_t, int8_t, or float
@tparam WINDOW_LENGTH how many recent input values to include in the calculations.
*/
template <class T, int WINDOW_LENGTH>
class RollingStat
{
public:
 
    /** Constructor */
    RollingStat() : _previous_mean(0), _mean(0), _variance(0), WINDOW_LENGTH_AS_RSHIFT(trailingZerosConst(WINDOW_LENGTH))
    {}
 
 
    /** Update the mean and variance given a new input value.
    @param x the next input value
    @note timing for unsigned int 10us, int 22us
    */
    void update(T x) {
        _mean = rollingMean.next(x);
        _variance = (((long)x - _previous_mean)*((long)x - _mean))>>WINDOW_LENGTH_AS_RSHIFT; // should really be div by (WINDOW_LENGTH-1), but exact values are not important for this application
        _previous_mean = _mean;
    }
 
 
    /** Update the mean and variance given a new input value.
    @param x the next input value
    */
    void update(int8_t x) {
        _mean = rollingMean.next(x);
        _variance = (((int)x - _previous_mean)*((int)x - _mean))>>WINDOW_LENGTH_AS_RSHIFT; // should really be div by (WINDOW_LENGTH-1), but exact values are not important for this application
        _previous_mean = _mean;
    }
 
 
    /** Return the mean of the last WINDOW_LENGTH number of inputs.
    @return  mean
    */
    T getMean() const {
        return _mean;
    }
 
 
    /** Return the approximate variance of the last WINDOW_LENGTH number of inputs.
    @return  variance
    @note This should really be calculated using WINDOW_LENGTH-1, but sacrificing accuracy for speed
    we use the power of two value of WINDOW_LENGTH.
    */
    T getVariance() const {
        return _variance;
    }
 
    /** Return the approximate standard deviation of the last WINDOW_LENGTH number of inputs.
    @return  standard deviation.
    */
    T getStandardDeviation() const {
        return isqrt16(_variance);
    }
 
 
 
private:
    T _previous_mean, _mean, _variance;
    RollingAverage <T, WINDOW_LENGTH> rollingMean;
    const uint8_t WINDOW_LENGTH_AS_RSHIFT;
};
 
// no need to show the specialisations
/** @cond  */
 
/** float specialisation of RollingStat template */
template <int WINDOW_LENGTH>
class RollingStat <float, WINDOW_LENGTH>
{
public:
 
    /** Constructor */
    RollingStat() : _previous_mean(0), _mean(0), _variance(0), WINDOW_LENGTH_AS_RSHIFT(trailingZerosConst(WINDOW_LENGTH))
    {}
 
 
    /** Update the mean and variance given a new input value.
    @param x the next input value
    @note timing for float: 60 to 90us
    */
    void update(float x) {
        _mean = rollingMean.next(x);
        _variance = ((x - _previous_mean)*(x - _mean))/(WINDOW_LENGTH-1);
        _previous_mean = _mean;
    }
 
 
    /** Return the mean of the last WINDOW_LENGTH number of inputs.
    @return  mean
    */
    float getMean() const {
        return _mean;
    }
 
 
    /** Return the approximate variance of the last WINDOW_LENGTH number of inputs.
    @return  variance
    */
    float getVariance() const {
        return _variance;
    }
 
    /** Calculate and return the approximate standard deviation of the last WINDOW_LENGTH number of inputs.
    @return  standard deviation.
    @note this is probably too slow to use!
    */
    float getStandardDeviation() const {
        return sqrt(_variance);
    }
 
 
 
private:
    float _previous_mean, _mean, _variance;
    RollingAverage <float, WINDOW_LENGTH> rollingMean;
    const uint8_t WINDOW_LENGTH_AS_RSHIFT;
};
 
// no need to show the specialisations
/** @endcond  */
 
#endif        //  #ifndef ROLLINGSTAT_H