#ifndef __FILE_SIMPLE_H_SEEN__
#define __FILE_SIMPLE_H_SEEN__ 

#include "libscl.h"

struct stats {
  REAL mean;
  REAL sdev;
  REAL var;
  REAL skew;
  REAL kurt;
  INTEGER nobs;
};
  
template <class P> stats simple(const P begin, const P end)
{
  stats ret;
  ret.mean = 0.0;
  ret.nobs = 0;
  
  P x = begin;
  while (x != end) {
    ret.mean += REAL(*x++);
    ++ret.nobs;
  }

  if (ret.nobs<2) 
    scl::error("Error, simple, not enough data, need two or more obs");

  ret.mean /= REAL(ret.nobs);

  REAL adev = 0.0;
  REAL r,p;
  ret.var = ret.skew = ret.kurt = 0.0;

  x = begin;
  while (x != end) {
    adev += (r = REAL(*x++) - ret.mean) > 0 ? r : -r;
    ret.var  += (p = r*r);
    ret.skew += (p *= r);
    ret.kurt += (p *= r);
  }

  adev /= REAL(ret.nobs);
  ret.var  /= REAL(ret.nobs-1);
  ret.sdev = sqrt(ret.var);
  if (ret.var) {
    ret.skew /= (REAL(ret.nobs) * ret.var * ret.sdev);
    ret.kurt =  ret.kurt /(REAL(ret.nobs) * ret.var * ret.var)-3.0;
  } 
  else {
    ret.skew = ret.kurt = REAL_MAX;
  }

  return ret;
}

#endif