From f7345ea71779db6bc35a9d5336813ad4d104d21e Mon Sep 17 00:00:00 2001
From: Antonin Portelli <antonin.portelli@me.com>
Date: Mon, 17 Feb 2014 18:48:44 +0000
Subject: [PATCH] improvements of StatArray and derived classes

---
 latan/Dataset.hpp   |   7 +--
 latan/Sample.cpp    |  29 ++++++++-
 latan/Sample.hpp    |  83 +++++++-----------------
 latan/StatArray.hpp | 149 +++++++++++++++++++++++++++++---------------
 4 files changed, 151 insertions(+), 117 deletions(-)

diff --git a/latan/Dataset.hpp b/latan/Dataset.hpp
index d712270..14bd062 100644
--- a/latan/Dataset.hpp
+++ b/latan/Dataset.hpp
@@ -39,6 +39,7 @@ private:
     typedef StatArray<T> Base;
 public:
     // constructor
+    using Base::Base;
     Dataset(void);
     Dataset(const std::string &listFileName, const std::string &dataName);
     template <typename Derived>
@@ -71,12 +72,6 @@ Dataset<T, FileType>::Dataset(const std::string &listFileName,
     load(listFileName, dataName);
 }
 
-template <typename T, typename FileType>
-template <typename Derived>
-Dataset<T, FileType>::Dataset(const Eigen::EigenBase<Derived> &dataset)
-: Base(dataset)
-{}
-
 // destructor //////////////////////////////////////////////////////////////////
 template <typename T, typename FileType>
 Dataset<T, FileType>::~Dataset(void)
diff --git a/latan/Sample.cpp b/latan/Sample.cpp
index b21dad2..58ac376 100644
--- a/latan/Sample.cpp
+++ b/latan/Sample.cpp
@@ -23,8 +23,33 @@
 using namespace Latan;
 using namespace std;
 
-template <>
-IoObject::IoType Sample<DMat>::getType(void) const
+/******************************************************************************
+ *                        DMatSample implementation                           *
+ ******************************************************************************/
+// constructors ////////////////////////////////////////////////////////////////
+DMatSample::DMatSample(void)
+: Sample<DMat>()
+{}
+
+DMatSample::DMatSample(const unsigned int nSample, const unsigned int nRow,
+                       const unsigned int nCol)
+: Sample<DMat>(nSample)
+{
+    resizeMat(nRow, nCol);
+}
+
+// destructor //////////////////////////////////////////////////////////////////
+DMatSample::~DMatSample(void)
+{}
+
+// resize all matrices /////////////////////////////////////////////////////////
+void DMatSample::resizeMat(const unsigned int nRow, const unsigned int nCol)
+{
+    this->unaryExpr([nRow, nCol](DMat &m){m.resize(nRow, nCol);});
+}
+
+// IO type /////////////////////////////////////////////////////////////////////
+IoObject::IoType DMatSample::getType(void) const
 {
     return IoType::dMatSample;
 }
diff --git a/latan/Sample.hpp b/latan/Sample.hpp
index 7453da8..739dfdf 100644
--- a/latan/Sample.hpp
+++ b/latan/Sample.hpp
@@ -27,16 +27,18 @@
 
 BEGIN_NAMESPACE
 
-const int central = -1;
+#define SAMPLE_OFFSET 1u
+
+const int central = -static_cast<int>(SAMPLE_OFFSET);
 
 /******************************************************************************
  *                              Sample class                                  *
  ******************************************************************************/
 template <typename T>
-class Sample: public StatArray<T>, public IoObject
+class Sample: public StatArray<T, SAMPLE_OFFSET>
 {
 private:
-    typedef StatArray<T> Base;
+    typedef StatArray<T, SAMPLE_OFFSET> Base;
 public:
     // constructors
     Sample(void);
@@ -45,74 +47,35 @@ public:
     Sample(const Eigen::EigenBase<Derived> &s);
     // destructor
     virtual ~Sample(void);
-    // operators
-    T& operator[](const int s);
-    // IO type
-    virtual IoType getType(void) const;
-private:
-    // index of the first element to take into account for statistics
-    virtual unsigned int getOffset(void) const;
 };
 
-template <>
-IoObject::IoType Sample<DMat>::getType(void) const;
-
-// specialization aliases
-typedef Sample<DMat> DMatSample;
+/******************************************************************************
+ *                          DMatSample class                                  *
+ ******************************************************************************/
+class DMatSample: public Sample<DMat>, public IoObject
+{
+public:
+    // constructors
+    DMatSample(void);
+    DMatSample(const unsigned int nSample, const unsigned int nRow,
+               const unsigned int nCol);
+    using Sample<DMat>::Sample;
+    // destructor
+    virtual ~DMatSample(void);
+    // resize all matrices
+    void resizeMat(const unsigned int nRow, const unsigned int nCol);
+    // IO type
+    virtual IoType getType(void) const;
+};
 
 /******************************************************************************
  *                    Sample class template implementation                    *
  ******************************************************************************/
-// constructor /////////////////////////////////////////////////////////////////
-template <typename T>
-Sample<T>::Sample(void)
-: Base(static_cast<typename Base::Index>(getOffset()))
-{}
-
-template <typename T>
-Sample<T>::Sample(const unsigned int nSample)
-: Base(static_cast<typename Base::Index>(nSample + getOffset()))
-{}
-
-template <typename T>
-template <typename Derived>
-Sample<T>::Sample(const  Eigen::EigenBase<Derived> &s)
-: Base(s)
-{}
-
 // destructor //////////////////////////////////////////////////////////////////
 template <typename T>
 Sample<T>::~Sample(void)
 {}
 
-// operators ///////////////////////////////////////////////////////////////////
-template <typename T>
-T& Sample<T>::operator[](const int s)
-{
-    if (s >= 0)
-    {
-        return Base::operator[](s + 1);
-    }
-    else
-    {
-        return Base::operator[](0);
-    }
-}
-
-// IO type /////////////////////////////////////////////////////////////////////
-template <typename T>
-IoObject::IoType Sample<T>::getType(void) const
-{
-    return IoType::noType;
-}
-
-// statistics //////////////////////////////////////////////////////////////////
-template <typename T>
-unsigned int Sample<T>::getOffset(void) const
-{
-    return 1u;
-}
-
 END_NAMESPACE
 
 #endif // Latan_Sample_hpp_
diff --git a/latan/StatArray.hpp b/latan/StatArray.hpp
index 6dfa08a..291abcf 100644
--- a/latan/StatArray.hpp
+++ b/latan/StatArray.hpp
@@ -22,13 +22,14 @@
 
 #include <latan/Global.hpp>
 #include <latan/Mat.hpp>
+#include <iostream>
 
 BEGIN_NAMESPACE
 
 /******************************************************************************
  *                     Array class with statistics                            *
  ******************************************************************************/
-template <typename T>
+template <typename T, unsigned int offset = 0>
 class StatArray: public Eigen::Array<T, Eigen::Dynamic, 1>
 {
 private:
@@ -41,101 +42,151 @@ public:
     StatArray(const Eigen::EigenBase<Derived> &s);
     // destructor
     virtual ~StatArray(void);
+    // access
+    unsigned int size(void) const;
+    // operators
+          T & operator[](const int s);
+    const T & operator[](const int s) const;
     // statistics
-    T mean(void)     const;
-    T variance(void) const;
-private:
-    // index of the first element to take into account for statistics
-    virtual unsigned int getOffset(void) const;
-    // operations for reduction in statistical computations
-    static inline T square(const T &a);
-    static inline T sum(const T &a, const T &b);
+    void bin(unsigned int binSize);
+    T    mean(const unsigned int pos, const unsigned int n) const;
+    T    mean(void) const;
+    T    variance(const unsigned int pos, const unsigned int n) const;
+    T    variance(void) const;
 };
 
-template <>
-inline DMat StatArray<DMat>::square(const DMat &a);
+// reduction operations
+namespace ReducOp
+{
+    template <typename T>
+    inline T square(const T &a);
+    template <typename T>
+    inline T sum(const T &a, const T &b);
+    template <>
+    inline DMat square(const DMat &a);
+}
 
 /******************************************************************************
  *                 StatArray class template implementation                    *
  ******************************************************************************/
 // constructors ////////////////////////////////////////////////////////////////
-template <typename T>
-StatArray<T>::StatArray(void)
-: Base(static_cast<typename Base::Index>(1))
+template <typename T, unsigned int offset>
+StatArray<T, offset>::StatArray(void)
+: Base(static_cast<typename Base::Index>(offset))
 {}
 
-template <typename T>
-StatArray<T>::StatArray(const unsigned int size)
-: Base(static_cast<typename Base::Index>(size))
-{}
-
-template <typename T>
-template <typename Derived>
-StatArray<T>::StatArray(const  Eigen::EigenBase<Derived> &s)
-: Base(s)
+template <typename T, unsigned int offset>
+StatArray<T, offset>::StatArray(const unsigned int size)
+: Base(static_cast<typename Base::Index>(size + offset))
 {}
 
 // destructor //////////////////////////////////////////////////////////////////
-template <typename T>
-StatArray<T>::~StatArray(void)
+template <typename T, unsigned int offset>
+StatArray<T, offset>::~StatArray(void)
 {}
 
+// access //////////////////////////////////////////////////////////////////////
+template <typename T, unsigned int offset>
+unsigned int StatArray<T, offset>::size(void) const
+{
+    return Base::size() - offset;
+}
+
+// operators ///////////////////////////////////////////////////////////////////
+template <typename T, unsigned int offset>
+T & StatArray<T, offset>::operator[](const int s)
+{
+    return Base::operator[](s + offset);
+}
+
+template <typename T, unsigned int offset>
+const T & StatArray<T, offset>::operator[](const int s) const
+{
+    return Base::operator[](s + offset);
+}
+
+
 // statistics //////////////////////////////////////////////////////////////////
-template <typename T>
-T StatArray<T>::mean(void) const
+template <typename T, unsigned int offset>
+void StatArray<T, offset>::bin(unsigned int binSize)
+{
+    unsigned int q = size()/binSize, r = size()%binSize;
+
+    for (unsigned int i = 0; i < q; ++i)
+    {
+        (*this)[i] = mean(i*binSize, binSize);
+    }
+    if (r != 0)
+    {
+        (*this)[q] = mean(q*binSize, r);
+        this->conservativeResize(offset + q + 1);
+    }
+    else
+    {
+        this->conservativeResize(offset + q);
+    }
+}
+
+template <typename T, unsigned int offset>
+T StatArray<T, offset>::mean(const unsigned int pos, const unsigned int n) const
 {
     T result;
-    unsigned int size = this->size() - getOffset();
     
-    if (size)
+    if (n)
     {
-        result = this->tail(size).redux(&StatArray<T>::sum);
+        result = this->segment(pos+offset, n).redux(&ReducOp::sum<T>);
     }
     
-    return result/static_cast<double>(size);
+    return result/static_cast<double>(n);
 }
 
-template <typename T>
-T StatArray<T>::variance(void) const
+template <typename T, unsigned int offset>
+T StatArray<T, offset>::mean(void) const
+{
+    return mean(0, size());
+}
+
+template <typename T, unsigned int offset>
+T StatArray<T, offset>::variance(const unsigned int pos, const unsigned int n) const
 {
     T s, sqs, result;
-    unsigned int size = this->size() - getOffset();
     
-    if (size)
+    if (n)
     {
-        s      = this->tail(size).redux(&StatArray<T>::sum);
-        sqs    = this->tail(size).unaryExpr(&StatArray<T>::square)
-                                 .redux(&StatArray<T>::sum);
-        result = sqs - square(s)/static_cast<double>(size);
+        s   = this->segment(pos+offset, n).redux(&ReducOp::sum<T>);
+        sqs = this->segment(pos+offset, n).unaryExpr(&ReducOp::square<T>)
+                                          .redux(&ReducOp::sum<T>);
+        result = sqs - ReducOp::square(s)/static_cast<double>(n);
     }
     
-    return result/static_cast<double>(size - 1);
+    return result/static_cast<double>(n - 1);
 }
 
+template <typename T, unsigned int offset>
+T StatArray<T, offset>::variance(void) const
+{
+    return variance(0, size());
+}
+
+// reduction operations ////////////////////////////////////////////////////////
 template <typename T>
-inline T StatArray<T>::sum(const T &a, const T &b)
+inline T ReducOp::sum(const T &a, const T &b)
 {
     return a + b;
 }
 
 template <typename T>
-inline T StatArray<T>::square(const T &a)
+inline T ReducOp::square(const T &a)
 {
     return a*a;
 }
 
 template <>
-inline DMat StatArray<DMat>::square(const DMat &a)
+inline DMat ReducOp::square(const DMat &a)
 {
     return a.cwiseProduct(a);
 }
 
-template <typename T>
-unsigned int StatArray<T>::getOffset(void) const
-{
-    return 0u;
-}
-
 END_NAMESPACE
 
 #endif // Latan_StatArray_hpp_