#ifndef __FILE_SIMPLE_H_SEEN__
#define __FILE_SIMPLE_H_SEEN__ 

/*-----------------------------------------------------------------------------

Copyright (C) 2012

A. Ronald Gallant
Post Office Box 659
Chapel Hill NC 27514
USA

This program 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 2 of the License, or
(at your option) any later version.

This program 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 for more details.

You should have received a copy of the GNU General Public License along
with this program; if not, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.

-----------------------------------------------------------------------------*/

#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