[Home] [News] [API] [Example usage] [FAQ] [Roadmap] [Download] [Contact]

/***************************************************************************
                          cmatrix.h  -  description
                             -------------------
    begin                : Thu Feb 10 2000
    copyright            : (C) 2000 by Rudolph Pienaar
    email                : pienaar@bme.ri.ccf.org
 ***************************************************************************/
//
//				NOTE
// Portions of this library are based on source code orginally written by
// Bruce Eckel. These portions are copyright either Osborne/McGraw-Hill or
// copyright Bruce Eckel. See http://www.mindview.net
//
// Additional contributors include:
// Antonio	(alrl1@alu.um.es)
//
//
/***************************************************************************
 *                                                                         *
 *   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.                                   *
 *                                                                         *
 ***************************************************************************/
//
// NAME
//
//	cmatrix.h	$Id: cmatrix.h,v 1.10 2002/02/11 01:28:39 pienaar Exp $
//
// DESCRIPTION
//
//	`cmatrix.h' contains the class declarations for the associated
//	source cmatrix.cpp.
//
// HISTORY
//
// < 6-10-1999
// o Unused, gathering dust on the filesystem. After some testing several
//   years ago, the code was found to run about 2 -> 3 times as slow as
//   conceptually comparable and similarly designed `c' code.
//
// 6-10-1999
// o Resurrected for usage in reinforcement learning problems, as well as
//   with an eye toward possible future java implementation.
//
// 6-14-1999
// o public deepcopy routine added.
//
// 6-16-1999
// o After much thought and deliberation, I decided to re-code the copy
//   constructor and assignment overload with a "deepcopy" approach. This will
//   arguably be slower than the current "point to data" and obliviates
//   the need for the refcnt element.
//
// o Hmmm... I'm definitely oscillating here. It seems that the original copy
//   constructor will still work. Since matrices are built dynamically, heap
//   space remains "alive" between scope changes. It is thus safe to
//   copy-by-point.
//
// 6-18-1999
// o Matrix structural manipulations routines added and debugged.
//
// 6-24-1999
// o scale operation added
//
// 7-7-1999
// o Added randomly_fill
//
// 7-8-1999
// o Added rms
//
// 7-13-1999
// o Expanded print routine
//
// 7-23-1999
// o Added findAndReplace
//
// 8-16-1999
// o Fixed rms routine
// o Added `sum' routine
//
// 9-01-1999
// o Added quantize routine
//
// 9-13-1999
// o Expanded print routine
//
// 2-18-2000
// o Added replicate routine
//
// 2-24-2000
// o Added base10 conversion routines
// o Added max/min_find() routines
// o Some miscellaneous cleaning of code
//
// 2-29-2000
// o Added `innerProd'
//
// 04 Aug 2000
// o Added 'quantize_linearly'
//
// 06 Aug 2000
// o Added 'find_quantized'
//
// 08 Aug 2000
// o Added continuous->quantized method
//
// 23 Aug 2000
// o Added stepPyramid_distribtion method
//
// 23 September 2000
// o Added diagonalise
//
// 25 September 2000
// o Added sqrdDistance
//
// 07 November 2000
// o Expanded normalisation to allow for different types
// o Added innerSum
//
// 10 November 2000
// o Added dot()
//
// 14 November 2000
// o Changed return types of is_() methods from int to bool
// o Added sort()
//
// 02 December 2000
// o Began fixing "local" return values in *_remove() methods
//
// 17 January 2001
// o Cleaned code for -Wall compile
//
// 20 March 2001
// o Added sprint()
//
// 02 April 2001
// o Added perimeter_set()
//
// 04 April 2001
// o Fixed up find_quantized()
//
// July 2001
// o Added statistical routines
// o Changed references of `long int' to `int'
//
// 12 September 2001
// o Added quantise_linearlyShifted()
//
// 4 October 2001
// o Added isOutOfBounds()
//
// 10 February 2002
// o Added _CMDATA as typedef placeholder for matrix data type
//


#ifndef __CMATRIX_H__
#define __CMATRIX_H__

#ifndef _CMDATA
#define _CMDATA double
#endif

#include 

typedef enum {
    e_matrix, e_rowVector, e_columnVector, e_square
} e_MATRIXTYPE;

typedef enum {
    e_row, e_column
} e_DOMINANCE;

typedef enum {
    e_MAGNITUDE, e_MAXIMUM, e_SUM, e_MEANSTD
} e_NORMALISATIONTYPE;

typedef enum {
    e_ascending, e_descending
} e_SORTTYPE;

class CMatrix {

    // Data members - declared as protected so that derived classes can
    // also access information

 protected:
    struct matstruct {
	_CMDATA 	**data;
	int 	rows, cols;
	int 	refcnt;
    } *ps_mat;

    // Public member functions

 public:

    //
    // Constructor, destructor, and miscellaneous initialization
    // routines
    //
    CMatrix(    int        	mrows=1,
                int        	mcols=1,
                _CMDATA         initvalue=0.0);
    CMatrix(    int        	mrows,
                int        	mcols,
                CMatrix         V,
                e_DOMINANCE     dominance);
    CMatrix(    CMatrix         M,
                e_DOMINANCE     dominance,
                e_MATRIXTYPE    type);
    CMatrix(    int        	mrows,
                int        	mcols,
                const _CMDATA*	data);
    CMatrix(    char*           flag,                   // create an identity
                int        	dimension);             //      matrix
    CMatrix(    char*           matfile);               // read from a CMATRIX
                                                        //      file
    CMatrix(    CMatrix& x);			        // copy initializer
    ~CMatrix();                                         // destructor

    CMatrix	reinitialise(	const _CMDATA* apf_data,
    				int a_elements = 0);	// overwrite internal
    							//	data with
    							//	new array
    int		positionVectorise(CMatrix& aM_dim);	// given a vector of
    							//	dimensions and
    							//	sizes, implements
    							//	a N->1 positional
    							//	mapping
    void        randomize(      _CMDATA	lower,          // randomize matrix
                                _CMDATA	upper,          //      cells with
		                int     asInt=0);       //      values in given
		                                        //      range

    void        randomly_fill(  _CMDATA	af_fillValue,   // fill a matrix with
                                _CMDATA	af_density);    //      values at given
                                                        //      density
                                                        //

    int         hardlimit(      _CMDATA   low,		// a hardlimiting
                                _CMDATA   high,		//      filter
		                _CMDATA   tolow,
		                _CMDATA   tohigh);
    void        quantize(       _CMDATA   f_grain); 	// quantization
    bool        quantizeIndex(  CMatrix & M_lookup,     // quantize a vector
                                CMatrix *&pM_index);	//     into a lookup
                                                        //      table and
                                                        //      return the
                                                        //      lookup index
    void        quantize_linearlyShifted(
    					_CMDATA af_start, // quantize a vector
                                        _CMDATA af_stop   //      over the given
                                                        //      range
                                );
    void        quantize_linearly(      _CMDATA af_start, // quantize a vector
                                        _CMDATA af_stop   //      over the given
                                                        //      range
                                );
    void	delta_add(      _CMDATA   af_delta,	// adds a delta value
                                _CMDATA   af_lower,	//      to matrix cells
                                _CMDATA   af_upper);	//      in the defined
                                                        //      range


    // Access and set internal values
    int         rows_get()
                                const {return ps_mat->rows;};
    int         cols_get()
                                const {return ps_mat->cols;};
    _CMDATA** 	data_get()
                                const {return ps_mat->data;};
    _CMDATA*	dataAsArray_get(	_CMDATA*	apf_data,
    					int		a_size = 0);
    int	        get_rows()
                                const {return ps_mat->rows;};
    int	        get_cols()
                                const {return ps_mat->cols;};
    _CMDATA**	get_data()
                                const {return ps_mat->data;};
    _CMDATA &	val(int row, int col);
                                                        // element selection:
                                                        //      read or
                                                        //      write

    // Output routines
    void	print(  char*	apch_msg 	= "ans",// print a matrix in
			int 	a_precision 	= 6,	//	a variety
			int	a_width		= 12,	//	of ways
			int	a_userFormat	= 0,
			int 	a_marginOffset 	= 0,
			char 	ach_left 	= (char) 0,
			int 	a_tabOffset 	= 0,
			bool 	ab_fancy 	= 1);
    string	sprint (
			char*	apch_msg 	= "ans",
			int 	a_precision 	= 6,
			int	a_width		= 12,
			int	a_userFormat	= 0,
			int 	a_marginOffset 	= 0,
			char 	ach_left 	= (char) 0,
			int 	a_tabOffset 	= 0,
			bool 	ab_fancy 	= 1);
    void	fprint(
			string	astr_filename,
			char*	apch_msg	= "",
			int 	a_precision 	= 6,
			int	a_width		= 12,
			int	a_userFormat	= 0
    );
    void	write_file(     char*   filename,       // write a matrix
                                char*   msg     ="");   //      to file

    // Matrix information routines
    bool		is_outOfBounds(			// checks if this
    		CMatrix&	aM_boundaryUpper,	// 	exceeds upper
		CMatrix&	aM_boundaryLower,	//	or lower
		bool&		ab_violate		// 	boundaries.
			);
    bool		is_vector();                    // checks if object is
                                                        //      a vector
    bool		is_column_vector();             // checks if object is
                                                        //      column vector
    bool		is_row_vector();	        // checks if object is
                                                        //      row vector
    bool		is_square();                    // checks if matrix is
                                                        //      square
    CMatrix*            size();                         // returns the matrix
                                                        //      size
    int                 size1D();                       // returns the 1D size
                                                        //      of the matrix,
                                                        //      i.e. the rows
                                                        //      multiplied by
                                                        //      column size
    e_MATRIXTYPE	matrix_type();	                // returns the matrix
                                                        //      type
    bool		compatible(     int     dim1,   // checks if *this is
                                        int     dim2);	//      compatible
    bool		compatible(const CMatrix& M);	//      with passed
                                                        //      dimensions

    // Matrix structural manipulation routines
    CMatrix	copy(   const   CMatrix &source);       // explicit deepcopy
    CMatrix	matrix_remove(  CMatrix*&       apM,
                                int             a_trow, // removes a submatrix
                                int             a_tcol,	//      from a base
			        int             a_rows,
			        int             a_cols);
    CMatrix	matrix_replace( int             a_trow,  // replaces part of
                                int             a_tcol,  //      base matrix
			        CMatrix&        aM_replacement);
    CMatrix	row_remove(     CMatrix*&       apV,    // removes a row
                                int             a_rownum);
    CMatrix	row_insert(     CMatrix&        aV,      // inserts a row
                                int             a_rownum);
    CMatrix	row_replace(    CMatrix&        aV,      // replaces a row
                                int             a_rownum);
    CMatrix	col_remove(     CMatrix*&       apV,    // removes a column
                                int             colnum);
    CMatrix	col_insert(     CMatrix&        aV,     // inserts a column
                                int             a_colnum);
    CMatrix	col_replace(    CMatrix&        aV,     // replaces a column
                                int             tcol);
    CMatrix	diagonal_replace(CMatrix & v,           // replaces a diagonal
				 int    dominant = 1);	//      defaults to
					                //      dominant
					                //      diagonal (i.e
					                //      top left to
					                //      bottom right)
    CMatrix     diagonalise(    _CMDATA   af_offdiag = 0.0);// creates an NxN matrix
                                                        //      with diagonal `v'
                                                        //      and off_diagonal
                                                        //      value `f_offdiag'
    CMatrix	replicate(      int             times,  // replicate a matrix
                                e_DOMINANCE     dir);
    CMatrix     sort(   e_SORTTYPE      e_sorttype      // sorts a matrix
                                = e_ascending,
                        e_DOMINANCE     e_dominance
                                = e_row);

    // Search/replace/fill routines
    CMatrix	findAll(        _CMDATA what,     int& occurences);
                                                        // searches for all
                                                        //      instances of
                                                        //      `what'
    int		findAndReplace( _CMDATA target,   _CMDATA source);
    bool        find_quantized( _CMDATA af_what,
                                CMatrix*        pM_index);
                                                        // find the position
                                                        //      corresponding
                                                        //      to the
                                                        //      quantum of
                                                        //      `af_what'
    void	perimeter_set(	_CMDATA 	af_value,
    				int	offset = 0);	// sets the perimeter
    							//	of the matrix
    							//	to af_value

    // Simple statistics
    _CMDATA	mean();					// mean of *all*
                                                        //      elements
    _CMDATA	mean(	e_DOMINANCE	ae_type,	// mean of a particular
    			int		a_index);	//	row or column
    CMatrix	mean(	e_DOMINANCE	ae_type);	// mean value on a
    							//	per row/col
    							//	basis
    _CMDATA	std(	bool	ab_sample = false);	// std dev of *all*
                                                        //      elements
    _CMDATA	std(	e_DOMINANCE	ae_type,	// std dev of a particular
    			int		a_index,	//	row or column
    			bool	ab_sample = false);
    CMatrix	std(	e_DOMINANCE	ae_type,	// std dev value on a
    			bool	ab_sample = false);	//	per row/col
    							//	basis
    _CMDATA	sum(	bool            ab_abs = false);// finds the sum of all
                                                        //      elements
    _CMDATA	sum(	e_DOMINANCE	ae_type,	// sum of a particular
    			int		a_index,	//	row or column
    			bool		ab_abs = false);//
    CMatrix	sum(	e_DOMINANCE	ae_type,	// sum value on a
    			bool		ab_abs = false);//	per row/col
    							//	basis
    _CMDATA	min();                                  // find the minimum
    _CMDATA	min(	e_DOMINANCE	ae_type,	// min of a particular
    			int		a_index);	//	row or column
    CMatrix	min(	e_DOMINANCE	ae_type);	// min value on a
    							//	per row/col
    							//	basis
    _CMDATA	max_filter(     CMatrix & M_mask);      // find the max
                                                        //      (filtering out
                                                        //      M_mask)
    _CMDATA	min_filter(     CMatrix & M_mask);      // find the min
                                                        //      (filtering out
                                                        //      M_mask)
    _CMDATA	min_epsilon(    _CMDATA   af_center,	// find the epsilon
                                _CMDATA   af_range);	//      min
    _CMDATA	min_find(       CMatrix& where);        // find the location of
                                                        //      the minimum
    _CMDATA	minabs();				// find the absolute
                                                        //      minimum
    _CMDATA	max();					// find the maximum
    _CMDATA	max(	e_DOMINANCE	ae_type,	// max of a particular
    			int		a_index);	//	row or column
    CMatrix	max(	e_DOMINANCE	ae_type);	// max value on a
    							//	per row/col
    							//	basis
    _CMDATA	max_epsilon(	_CMDATA   af_center,	// find the epsilon
                                _CMDATA   af_range);	//      max
    _CMDATA	max_find(       CMatrix& where);        // find the location of
                                                        //      maximum
    _CMDATA	maxabs();                               // find the absolute
                                                        //      maximum

    // Mathematical considerations
    CMatrix	functionApply(	double			// Apply a function
    				(*func)(double));	//	to each element
    _CMDATA	dot(            const CMatrix& rval);	// dot product
    CMatrix     add(CMatrix& M_target);                 // add M_target to *this
    _CMDATA	sqrdDistance(CMatrix & M_to);           // determines the squared
                                                        //      distance between
                                                        //      *this and M_to
    _CMDATA	distance(CMatrix & M_to);               // determines the distance
                                                        //      between *this and
                                                        //      M_to
    _CMDATA	determinant();                          // find the determinant
    CMatrix	inverse();                              // find the invese
    CMatrix     normalise(      e_NORMALISATIONTYPE     // normalise contents
                                ae_type = e_MAGNITUDE); //      of matrix
    CMatrix     normalise(      e_NORMALISATIONTYPE     // normalise contents
                                ae_type,                //      of matrix in either
                                e_DOMINANCE             //      a row or column
                                ae_dominance);          //      dominant manner
    void        flip(           _CMDATA   af_prob,        // flips (logical inverse)
                                _CMDATA   af_false = 0.0, //      matrix contents
                                _CMDATA   af_true = 1.0); //
    _CMDATA	mag();                                  // magnitude of matrix
    _CMDATA	mag(	e_DOMINANCE	ae_type,	// mag of a particular
    			int		a_index);	//	row or column
    CMatrix	mag(	e_DOMINANCE	ae_type);	// mag value on a
    							//	per row/col
    							//	basis
    CMatrix	abs();	                                // absolute value
    _CMDATA	variance();				// statistical variance
    CMatrix	lu_decompose(   CMatrix&        indx,   // L-U decomposition
                                int&       	d);
    void	lu_back_subst(  CMatrix&        indx,   // uses L-U decomp
                                CMatrix&        b);     //      for matrix
                                                        //      inverse
    CMatrix     scale(  const _CMDATA		rval);  // scales a matrix
    CMatrix     scale(  CMatrix& M);                    // does element for element
                                                        //      multiplication
    _CMDATA	innerProd();                            // finds the inner
                                                        //      product
                                                        //      (product of all
                                                        //      cells)
    _CMDATA	innerSum();                             // finds the inner
                                                        //      sum
                                                        //      (sum of all
                                                        //      cells)
    _CMDATA	rms(    CMatrix rval);                  // rms between *this
                                                        //      and rval
    int		b10_convertTo(  int     radix);         // convert a row vector
                                                        //      to a base 10
                                                        //      number
    CMatrix	b10_convertFrom(int     num,            // convert from number
                                int     radix,          //      in base 10 to
                                int     forcelength=0); //      radix, with optional
                                                        //      length
    bool        equal(  CMatrix &       aM_B);          // am I (numerically)
                                                        //      equal to B?
    _CMDATA	distanceTo(     CMatrix &       aM_B);  // distance from *this
                                                        //      to B

    // Operator overloads
    CMatrix	operator=(      const CMatrix& rval);   // matrix assignment
    CMatrix	operator+(      const CMatrix& rval);   // matrix addition
    CMatrix	operator+(      const _CMDATA rval);	// element by element
                                                        //      addition
    CMatrix	operator-(      const CMatrix& rval);   // matrix subtraction
    CMatrix	operator-(      const _CMDATA rval);	// element by element
                                                        //      subtraction
    CMatrix	operator-();				// unary minus
    CMatrix	operator*(      const CMatrix& rval);   // matrix mult
    CMatrix	operator*(      const _CMDATA rval);	// scalar mult
    CMatrix	operator/(      const _CMDATA rval);	// scalar division
    CMatrix	operator*=(     const CMatrix& rval);	// element by element
                                                        //      mult
    CMatrix	operator/=(     const CMatrix& rval);	// element by element
                                                        //      division
    CMatrix	operator>>=(    _CMDATA rval);		// element by element
                                                        //      comparison
                                                        //      with rval
    CMatrix	operator~();				// writes the elements
                                                        //      of rval into
                                                        //      a vector
    CMatrix	operator!();				// transpose a matrix
    _CMDATA	operator++();				// inner sum

    // Explicit mathematical function calls
    //	These were developed for test purposes
    //	and may later become obsolete
    void        multStore(
                        const   CMatrix&        M_A,    // multiply A and B, and
                        const   CMatrix&        M_B     //      store result in
                        );                              //      *this

    //
    // Private member functions
    //
 private:

    void	_error(         char*   	apch_proc,
                                char*   	apch_msg,
                                int     	code            = 1);
    void	_column_switch( int        	col1,
                                int        	col2);
    void	_column_copy(   CMatrix&        m,
                                int        	from_col,
                                int        	to_col);
    CMatrix	_scale();
    void	_deepcopy(      CMatrix&        from,
                                CMatrix&        to);	// deep copy matrix
    _CMDATA &	_mval(  int        row,
                        int        col)
                        const {return(ps_mat->data[row][col]);}
                                                        // used by matrix
                                                        //      functions that
                                                        //      know they
                                                        //      aren't
                                                        //      exceeding the
                                                        //      boundaries
};

//
// Non-class functions
//

CMatrix		col_insert(     CMatrix         mm,
                                CMatrix         vv,
                                int        	colnum);
CMatrix		row_insert(     CMatrix         mm,
                                CMatrix         vv,
                                int             rownum);
CMatrix		col_remove(     CMatrix         mm,
                                int 		colnum);
CMatrix		row_remove(     CMatrix         mm,
                                int             rownum);
CMatrix		insert_diagonal(CMatrix         mm,
                                CMatrix         vv);
int		dominant_eigen( CMatrix         mm,
                                _CMDATA 	tol,
                                _CMDATA* 	eigenvalue,
                                CMatrix*        eigenvector);
CMatrix		form_matrix(    CMatrix         mm,
                                CMatrix         vv);
int		b10_convertTo(	CMatrix M_num,	int radix);
CMatrix		b10_convertFrom(int num, 	int radix);

CMatrix         leastSquaresError_find(                 // return the least square
                                CMatrix&        aM_A,   //      error between A
                                CMatrix&        aM_b);  //      and b
CMatrix         linearSystem_solve(                     // return x, the solution
                                CMatrix&        aM_A,   //      to the linear sys
                                CMatrix&        aM_b);  //      Ax = b


#endif //__CMATRIX_H__