// -------------------------------------------------------------------------
// @author Leonardo Florez-Valencia (florez-l@javeriana.edu.co)
// -------------------------------------------------------------------------

#ifndef __CPEXTENSIONS__ALGORITHMS__ITERATIVEGAUSSIANMODELESTIMATOR__H__
#define __CPEXTENSIONS__ALGORITHMS__ITERATIVEGAUSSIANMODELESTIMATOR__H__

#include <vector>
#include <itkConceptChecking.h>
#include <itkObject.h>
#include <vnl/vnl_matrix.h>

namespace cpExtensions
  {
    namespace Algorithms
    {
      /**
       */
      template< class S, unsigned int D >
      class IterativeGaussianModelEstimator
        : public itk::Object
      {
      public:
        typedef IterativeGaussianModelEstimator Self;
        typedef itk::Object                     Superclass;
        typedef itk::SmartPointer< Self >       Pointer;
        typedef itk::SmartPointer< const Self > ConstPointer;

        typedef S TScalar;
        itkStaticConstMacro( Dimension, unsigned int, D );

        // Begin concept checking
#ifdef ITK_USE_CONCEPT_CHECKING
        itkConceptMacro(
          ScalarTypeHasFloatResolution,
          ( itk::Concept::IsFloatingPoint< S > )
          );
#endif
        // End concept checking

        typedef vnl_matrix< S >        TMatrix;
        typedef std::vector< TMatrix > TMatrices;

      public:
        itkNewMacro( Self );
        itkTypeMacro( IterativeGaussianModelEstimator, itkObject );

      public:
        unsigned long GetNumberOfSamples( ) const
          { return( ( unsigned long )( this->m_N ) ); }
        const TMatrix& GetMu( ) const
          { return( this->m_M ); }
        const TMatrix& GetOmega( ) const
          { return( this->m_O ); }

        void SetNumberOfSamples( unsigned long n );
        void SetMu( const TMatrix& m );
        void SetOmega( const TMatrix& O );

        bool SaveModelToFile( const std::string& filename ) const;
        bool LoadModelFromFile( const std::string& filename );

        template< class V >
        S Probability( const V& sample ) const;

        S Probability( const S& s_x, const S& s_y, ... ) const;

        template< class V, class M >
        void GetModel( V& m, M& E ) const;

        void Clear( );

        template< class V >
        void AddSample( const V& sample );

        void AddSample( const S& s_x, const S& s_y, ... );

      protected:
        IterativeGaussianModelEstimator( );
        virtual ~IterativeGaussianModelEstimator( );

      protected:
        void _UpdateModel( ) const;

      private:
        // Purposely not implemented
        IterativeGaussianModelEstimator( const Self& other );
        void operator=( const Self& other );

      protected:
        S m_N;
        TMatrix m_M;
        TMatrix m_O;

        mutable bool m_Updated;
        mutable TMatrix m_Cov;
        mutable TMatrix m_Inv;
        mutable S m_Norm;
      };

    } // ecapseman

} // ecapseman

#include <cpExtensions/Algorithms/IterativeGaussianModelEstimator.hxx>

#endif // __CPEXTENSIONS__ALGORITHMS__ITERATIVEGAUSSIANMODELESTIMATOR__H__

// eof - $RCSfile$