pca.hpp

/* $Id$ */

#ifndef JMATH_PCA_T_HPP
#define JMATH_PCA_T_HPP

#include "jafarConfig.h"

#ifdef HAVE_BOOST_SANDBOX
#ifdef HAVE_LAPACK

#include "jmath/jblas.hpp"
#include "kernel/jafarDebug.hpp"

namespace jafar {
  namespace jmath {

    template<typename NUMTYPE>
    class PCA_T {
      public:

      enum UFlag {
        increase, 
        preserve, 
        ore 
      };


      PCA_T(bool basisOnly_=false) : basis_only(basisOnly_) {};

      PCA_T(const ublas::matrix<NUMTYPE>& X_, int dim_=0, bool basisOnly_ = false) {
        basis_only = basisOnly_;
        batchPCA(X_,dim_);
      };

      PCA_T(PCA_T const & pca_) {
        mean = pca_.mean;
        eigenvalues = pca_.eigenvalues;
        eigenvectors = pca_.eigenvectors;
        coefficients = pca_.coefficients;
      }

      PCA_T& operator= (PCA_T const & pca_) {
        mean = pca_.mean;
        eigenvalues = pca_.eigenvalues;
        eigenvectors = pca_.eigenvectors;
        coefficients = pca_.coefficients;
        return *this;
      };

      //Accessors
      ublas::vector<NUMTYPE>& getMean() {
        JFR_PRECOND(mean.size() != 0, "PCA::getMean: no results available");
        return mean;
      };

      ublas::matrix<NUMTYPE>& getEigenVectors() {
        JFR_PRECOND(eigenvectors.size2() != 0, "PCA::getEigenVectors: no results available");
        return eigenvectors;
      };

      ublas::vector<NUMTYPE>& getEigenValues() {
        JFR_PRECOND(eigenvalues.size() != 0, "PCA::getEigenValues: no results available");
        return eigenvalues;
      }

      ublas::matrix<NUMTYPE>& getCoefficients() {
        JFR_PRECOND(coefficients.size2() != 0, "PCA::getEigenValues: no results available");
        return coefficients;
      };


      void batchPCA(const ublas::matrix<NUMTYPE>& X_, int dim_=0);


      void updatePCA(const ublas::vector<NUMTYPE>& I_, UFlag f_=preserve, NUMTYPE thd_=0.25);

      ublas::vector<NUMTYPE> project(const ublas::vector<NUMTYPE>& I_) const;

      ublas::vector<NUMTYPE> reconstruct(const ublas::vector<NUMTYPE>& I_) const;

      private:

#ifdef USE_JMATH_SERIALIZATION
      friend class boost::serialization::access;
      template<class Archive>
        void serialize(Archive & ar, const unsigned int version) {
        ar & BOOST_SERIALIZATION_NVP(basis_only);
        ar & BOOST_SERIALIZATION_NVP(eigenvectors);
        ar & BOOST_SERIALIZATION_NVP(coefficients);
        ar & BOOST_SERIALIZATION_NVP(mean);
        ar & BOOST_SERIALIZATION_NVP(eigenvalues);
      };
#endif
      bool basis_only;
      ublas::matrix<NUMTYPE> eigenvectors,coefficients;
      ublas::vector<NUMTYPE> mean, eigenvalues;

    }; // class PCA_T
  } // namespace jmath
} // namespace jafar

#include <jmath/pca.hxx>

namespace jafar{
  namespace jmath {
    typedef PCA_T<double> PCA;
  }
}
#endif // HAVE_LAPACK
#endif // HAVE_BOOST_SANDBOX
#endif // JMATH_PCA_T_HPP