matrix.cpp

adaboost code in matlab

Matrix Matrix::operator+ (const Matrix& m) const
{
	if(rows != m.row() || cols != m.col())
		reptError("Matrix sizes mismatch, addition failed!");

	Matrix temp(rows, cols);
	temp = *this;
	temp += m;

	return temp;
}

Matrix Matrix::operator- (const Matrix& m) const
{
	int r=rows;
	int c=m.col();
	if (c > 1)
	   c = c-1;
	else
	   c = c;
	Matrix temp(r, c);

	for (int i=0; i < r; i++)
	{
		for (int j=0; j < c; j++)
		{
			temp.element[i][j] = element[i][j] - m.element[i][j];
		}
	}

	return temp; 
 // if(rows != m.row() || cols != m.col())
 //		reptError("Matrix sizes mismatch, subtraction failed!");
    /*
	Matrix temp(rows, cols);
	temp = *this;
	temp -= m;

	return temp;
	*/
}

// compare the correct label matrix and scaler
Matrix Matrix::operator> (const T& c) const
{
	Matrix temp(rows, cols);
	temp = *this;
	temp >= c;

	return temp;
}

// combined scalar multiplication and assignment operator //
Matrix& Matrix::operator>= (const T& c)
{
	for (int j = 0; j < cols; j++){
		for (int i = 0; i < rows; i++){
			if (element[i][j] > c)
				element[i][j]=1;
			else
				element[i][j]=0;
		}
	}
	
	return *this;
}

// scale a matrix by a constant c //
Matrix Matrix::operator* (const T& c) const
{
	Matrix temp(rows, cols);
	temp = *this;
	temp *= c;

	return temp;
}

// two matrix multiplication // 
Matrix Matrix::operator* (const Matrix& m) const
{
//	if(cols != m.row())
//		reptError("Matrix dimensions error, multiplication failed!");

	int r=rows;
	int c=m.col();
	Matrix temp(r, c);

	for (int i=0; i < r; i++)
	{
		for (int j=0; j < c; j++)
		{
			for (int k=0; k < cols; k++)
			{
				temp.element[i][j] += element[i][k] * m.element[k][j];
				if(fabs(temp.element[i][j])<=AS_ZERO)
					temp.element[i][j]=0;			// round small # to zero
			}
		}
	}

	return temp;
}

void Matrix::scalarDivisonScalar(double scalar_1, double scalar_2, const Matrix& m) const
{
	int r=rows;;
	int c=m.col();
	for (int i=0; i < r; i++)
		for (int j=0; j < c; j++)
			 element[i][j] = scalar_1 / scalar_2;
}



void Matrix::matrixdDivisonScalar(double scalar, const Matrix& m) const
{
	int r=rows;
	int c=m.col();
//	Matrix temp(r, c);
	for (int i=0; i < r; i++){
		for (int j=0; j < c; j++){
			 element[i][j] = m.element[i][j] / scalar;
			//if(fabs(temp.element[i][j])<=AS_ZERO)
			//	temp.element[i][j]=0;			// round small # to zero
		}
	}
	
//	return temp;
}


// divison matrix by scaler 
Matrix Matrix::operator / (const Matrix& m) const
{
	//	if(cols != m.row())
	//		reptError("Matrix dimensions error, multiplication failed!");
	
	int r=rows;
	int c=m.col();
	Matrix temp(r, c);
	
	for (int i=0; i < r; i++)
	{
		for (int j=0; j < c; j++)
		{
			temp.element[i][j] = element[i][j] / m.element[i][j];
			if(fabs(temp.element[i][j])<=AS_ZERO)
				temp.element[i][j]=0;			// round small # to zero
		}
	}	
	return temp;
}


// set a matrix to an identity matrix //
void Matrix::toUnit(void)
{
	if (!isSquare())
		reptError("No identity matrix for a non-square matrix!");

	for (int i = 0; i < rows; i++)
		for (int j = 0; j < cols; j++)
			element[i][j]=(i==j) ? (T) 1 : (T) 0;
}

void Matrix::reptError (char* errorMsg) const
{
	cerr << errorMsg << endl;
	exit(1);
}

// compute the transpose matrix //
void Matrix::matrixTranspose (int r, int c, const Matrix& m) const
{
	for (int i = 0; i < r; i++)
		for (int j = 0; j < c; j++)
			element[j][i] = m.element[i][j];
}

// compute the transpose matrix //
Matrix Matrix::transpose (void) const
{
	int r=rows;
	int c=cols;
	Matrix temp(c, r);

	for (int i = 0; i < r; i++)
		for (int j = 0; j < c; j++)
	      temp.element[j][i] = element[i][j];

	return temp;
}
/*
void Matrix::TransposeMatrix (const Matrix& m) const
{
	int r=rows;
	int c=cols;
	Matrix temp(c, r);

	for (int i = 0; i < r; i++)
		for (int j = 0; j < c; j++)
	      temp.element[j][i] = element[i][j];

	return temp;
}
*/
// copy the matrix
void Matrix::copy (int tt, Matrix& m) const
{	
	int r=rows;
	int i = 0;
	int j = tt;
	double val;
	val = m.element[0][0];
	element[i][tt] = val;
}

// copy if the value is zero
void Matrix::ifValZero(double step,Matrix& m) const
{	
	double val;
	double one =1;
	val = m.element[0][0];
	if (val == 0.0){
		val = one/step;
	}
	if (val > 1.0 ){
        val = 0.999999;
	}
	element[0][0] = val;
}

int Matrix::ifValueIsMissing(int *missing, Matrix& m) const
{	
	int r=rows;
	int c=cols;
	double val=0.0;
	for (int i = 0; i < r; i++){
		for (int j = 0; j < c; j++){
	         val = m.element[i][j];
			 if (val == -999){
                (*missing)=1;
			     goto L10;
			 }
		}
	}
	L10:
	return (0);
}

// copy matrix to matrix
void Matrix::copyMatrixToMatrix (Matrix& m)
{	
	int r=rows;
	int c=cols;
	for (int i = 0; i < r; i++)
		for (int j = 0; j < c; j++)
	      element[i][j] = m.element[i][j];
}

// copy matrix to array
void Matrix::copyMatrixToArray (int nsample, int mdim, int array[], Matrix& m)
{	  
	int k=0;
	int j=0;
	int newmdim=mdim+1;
	double Val;
	for (k = 0; k < nsample; ++k) 
	{
		for (j = 0; j < mdim; ++j) 
		{
			Val = m.element[k][j];
			array[(j+1) + (k+1) * mdim - newmdim]= (int)Val;
		}
	 }
}

// copy array to array
void Matrix::copyArrayToArray (int nsample, int mdim, float array_1[], int array_2[])
{	
	for (int k = 1; k <= nsample; ++k) 
	{
		for (int j = 1; j <= mdim; ++j) 
		{
			//array_1[j + k * 36 - 37]= (float)array_2[j+k* 36 - 37];
			array_1[j + k * mdim - (mdim+1)]= (float)array_2[j+k* mdim - (mdim+1)];
		}
	 }
}
// copy array to matrix
void Matrix::copyArrayToMatrix (int nsample, int mdim, int array[], Matrix& m)
{	
    int k=0;
	int j=0;
	int newmdim=mdim+1;
	double Val;
	for (k = 0; k < nsample; ++k) 
	{
		for (j = 0; j <=0; ++j) 
		{
			Val = array[k];
			m.element[k][0]= Val;
		}
	 }
}

// copy array to matrix
void Matrix::copyVectorToMatrix (int k,Matrix& m)
{	
    int r=rows-1;
	for (int i = 0; i < r; i++)
	      element[i][k-1] = m.element[i][0];
}

void Matrix::copyToMatrix(int tt, Matrix& m) const
{
	double val;
	tt=tt-1;
	//for (int i=0; i < tt; i++){
	   val = m.element[0][tt];
	   element[tt][0] = val;
	//}
}


// minor() returns the minior of an element (r, c) //
void Matrix::minor(int r, int c, Matrix& m) const
{
	// precondition: a square matrix & the calling matrix with
	// dimension >= 3 //

	int k, l;
	for(int i = 0; i < rows-1; i++)
	{
		k=(i>=r) ? (i+1):i;
		for(int j = 0; j < cols-1; j++)
		{
			l=(j>=c) ? (j+1):j;
			m.element[i][j]=element[k][l];
		}
	}
}

// cofact() calculates the cofactor of an element at (r, c) //
T Matrix::cofact(int r, int c) const
{
	// precondition: a square matrix & the calling matrix with
	// dimension >= 3 //

	T cof;
	Matrix temp(rows-1, cols-1);

	minor(r, c, temp);
//	cout << endl;
//	temp.print();
	cof=(T) pow(-1, r+c)*temp.det();

	return cof;
}

// computes the adjoin matrix of a given matrix //
void Matrix::adjoin(Matrix& m) const
{
	int dim=rows;								// rows=cols

	if (!isSquare())
		reptError("No adjoin matrix for a non-square matrix!");
	if(m.row() != rows || m.col() != cols)
		reptError("Non square adjoin matrix defined!");

	for(int i = 0; i < dim; i++)
		for(int j = 0; j < dim; j++)
			m.element[i][j]=cofact(i, j);
//	cout << endl;
//	m.print();
	m=m.transpose();
//	cout << endl;
//	m.print();
}

// calculates the determinant of a matrix by first changing it
// to a upper triangular matrix //
T Matrix::det(void) const
{
	if (!isSquare())
		reptError("No determinant for a non-square matrix!");

	T determinant=1;
	Matrix temp(rows, cols);
	temp=*this;
//	cout << endl;
//	temp.print();

	temp.toUpper();
//	cout << endl;
//	temp.print();
	for (int i = 0; i < rows; i++)
		determinant *= temp.element[i][i];				// m(n, n) could be 0

	if(fabs(determinant)<=AS_ZERO)						// if too small, set to zero
		determinant=0;

	return determinant;
}

// covert to upper triangular matrix //
void Matrix::toUpper(void)
{
	// precondition: a square matrix //

	int i, j, k, dim=rows;
	for (j = 0; j < dim-1; j++)
	{
		// check if the diagonal element is zero //
		if(fabs(element[j][j])<=AS_ZERO)
		{
			element[j][j]=0;
			for(i = j+1; i < dim; i++)
			{
				if(element[i][j] != 0)
				{
					swapRow(i, j, dim);
					break;
				}
			}
			if(i>=dim)								// every element=0 //
			{
				cout << "Cann't convert to upper triangular matrix, ";
				cout << "Error in column No. j=" << j << endl;
				reptError("Every element in column j is zero!");
			}
		}

		// convert to upper triangular matrix //
		for (i = j+1; i < dim; i++)
		{
			if(element[i][j]==0)	continue;
			for(k = j+1; k < dim; k++)
				element[i][k]=element[i][k]-element[j][k]*element[i][j]/element[j][j];
			element[i][j]=0;
		}
	}
}

// swap 2 lines, each has c columns //
Matrix& Matrix::swapRow(int r1, int r2, int c)
{
	T temp;
	for (int j=0; j<c; j++)
	{
		temp=element[r1][j];
		element[r1][j]=element[r2][j];
		element[r2][j]=-temp;						// swapws, det()=-det();
	}

	return *this;									// have to return *this?
}

Matrix Matrix::inverse (void) const
{
	if (!isSquare())
		reptError("Not a square matrix, no inverse matrix!");
	if (det() == 0)
		reptError("Determinant=0, no inverse matrix!");

	int dim=rows;									// square matrix dimension
	Matrix temp(dim, dim);
	if(dim==1)
		temp=*this;
	else if(dim==2)
	{
		temp.element[0][0]=element[1][1];
		temp.element[0][1]=-element[1][1];
		temp.element[1][0]=-element[1][0];
		temp.element[1][1]=element[0][0];
		temp *= (1/det());
	}
	else
	{
		adjoin(temp);
		temp *= (1/det());
	}

	return temp;
}

void Matrix::readFile (char* file_name)
{
	// precondition: matrix dimension matches file format //

	ifstream is;

	is.open(file_name, ios::binary);
	if(is.fail())
	{
		cout << "File name = " << file_name << endl;
		reptError("Open the above file failed in readFile()!");
	}

	long i=0, j=0;
	while (is)
	{
		
		int ch = is.peek();
		if(ch==-1)									// end of file //
			break;
		if(ch==32 || (ch>=9 && ch<=13))				// skip whitespace //
		{
			is.get();
			continue;
		}
		if (ch=='#')
		{
			is.ignore(1000,'\n');					// skip comment line //
			continue;
		}

		for(j=0; j<cols; j++)
			is >> element[i][j];
		i++;
	}
	is.close();
//	if(i != rows)
//		reptError("File row# != rows!");
}

void Matrix::writeFile (char* file_name)
{
	ofstream fos;

	fos.open(file_name, ios::out, filebuf::sh_none);
	if(fos.fail())
	{
		cout << "File name = " << file_name << endl;
		reptError("Open file to write failed!");
	}

	fos << *this;
//	for(int i=0; i<rows; i++)
//	{
//		for(int j=0; j<cols; j++)
//			fos << element[i][j] << "\t";
//		fos << endl;
//	}
	fos.close();
}