#include<math.h>
#include<string.h>
#include<stdlib.h>
#include<stdio.h>
#include"LeastSquare.h"



// fitting Çϱâ À§ÇØ ÇÊ¿äÇÑ ÇÔ¼öµé..................

//********************************************************************************************************************************************************************************
//					¿ªÇà·ÄÀ» ±¸Çϱâ À§ÇÑ ÇÔ¼öµé...
//********************************************************************************************************************************************************************************

// transepose matrix¸¦ ±¸ÇÏ°í ±× Çà·ÄÀ» TrMatrix ¿¡ ³Ö´Â ÇÔ¼ö.....
//row : Çà column: ¿­...
void InverseFns::TranseposeMatrix(float** Matrix , float** TrMatrix , int row , int column) 
{
	int i ;
	int j ;
	
	for( i = 0 ; i < row ; i++ )
		for( j = 0 ; j < column ; j++ )
		TrMatrix[j][i] = Matrix[i][j] ;	
}


// ÀÓÀÇÀÇ µÎ Çà·ÄÀ»   °öÇÏ°í ¼¼¹ø° Çà·Ä¿¡ ³Ö´Â ÇÔ¼ö.......
// ¾Õ Çà ¿­ µÚ ¿­....
void InverseFns::MultiplyMatrixAB(float** MatrixA , float** MatrixB , float** MatrixAB , int row , int column , int bcolumn) 
{
	int i ; 
	int j ;
	int m ;
	float sum ;


	for( i = 0 ; i < row ; i++ )
	{
		for( j = 0 ; j < bcolumn ; j++ )
		{
			sum = 0.0 ;
		
			for( m = 0 ; m < column  ; m++ )
			{
				sum += MatrixA[i][m]*MatrixB[m][j] ;
			}

			MatrixAB[i][j] = sum;
		}		
	}	
}




// inverse matrix¸¦ ±¸ÇÏ´Â ÇÔ¼ö.... nSize : Â÷¼ö¸¦ ¸»ÇÑ´Ù....
int InverseFns::CalculatingInverseMatrix(float** Matrix , int nSize , float** InverseMatrix)
{
	int i , j , k , a;
	float c , temp;		

	// ´ÜÀ§ Çà·Ä·Î ¸¸µç´Ù....
	for( i = 0 ; i < nSize ; i++ )
		memset(InverseMatrix[i] , 0x00 , sizeof(float)*nSize);


	for( i = 0 ; i < nSize ; i++ )
		InverseMatrix[i][i] = 1;


	// ÀÌÁ¦ ½ÇÁ¦ ¿ªÇà·ÄÀ» ±¸ÇÑ´Ù.....
	for( j = 0 ; j < nSize ; j++ )
	{
		for( i = j ; i < nSize ; i++)
			if( Matrix[i][j] != 0.0 ) // óÀ½À¸·Î 0.0 ÀÌ ¾ÈµÇ´Â °÷À» ã´Â´Ù...
				break ;
		
		if( i == nSize )	//0.0 ÀÌ ¾ÈµÇ´Â °÷ÀÌ ¾øÀ¸¸é ¿ªÇà·ÄÀÌ ¾ø´Ù°í ÇÑ´Ù..
			return 0 ;		// ¿ªÇà·ÄÀÌ Á¸ÀçÇÏÁö ¾Ê´Â´Ù..... 
		

		c = Matrix[i][j] ;

		// óÀ½À¸·Î 0.0 ÀÌ ¾Æ´Ñ°÷À» 1·Î ¸¸µé±â À§ÇØ ³ª´©¾îÁØ´Ù... °öÇؼ­ °°Àº ¼ö·Î ¸¸µé¾î »©±â ½±°Ô ÇÏ·Á°í...
		if( c != 1.0 )
		{
			for( a = 0 ; a < nSize ; a++ )
			{
				Matrix[i][a] = Matrix[i][a]/c;
				InverseMatrix[i][a] = InverseMatrix[i][a]/c;
			}
		}

		// 0 ÀÌ ¾Æ´Ï¶ó°í ãÀº °÷À» j ¹ø° Çà°ú À§Ä¡¸¦ ¹Ù²ÙÀÚ...
		for( a = 0 ; a < nSize ; a++ )
		{ 
			temp = Matrix[j][a];
			Matrix[j][a] = Matrix[i][a];
			Matrix[i][a] = temp;

			temp = InverseMatrix[j][a];
			InverseMatrix[j][a] = InverseMatrix[i][a];
			InverseMatrix[i][a] = temp;
		}

		for( k = 0 ; k < nSize ; k++)
		{
			if(Matrix[k][j] == 0.0 || k == j)
				continue ;

			c = Matrix[k][j] ;

			for( a = 0 ; a < nSize ; a++ )
			{
				Matrix[k][a] = Matrix[k][a] - c*Matrix[j][a];
				InverseMatrix[k][a] = InverseMatrix[k][a] - c*InverseMatrix[j][a];
			}		
		}
	}// end of outer for loop ... j  ¿­ ¹øÈ£.... 

	//¿ªÇà·ÄÀ» ¼º°øÀûÀ¸·Î ±¸Çϸé...
	return 1 ;
}


//coeff = (Tr(zMatrix)*zMatrix)ÀÇ ¿ªÇà·Ä*Tr(zMatrix)*Y(½ÇÁ¦ µ¥ÀÌÅÍÀÇ y °ª.. ¿©±â¼± z °ªÀÌ µÈ´Ù... )
//½ÇÁ¦·Î coeffcient ¸¦ ±¸ÇÏ´Â ÇÔ¼ö....
void LeastSquareFns::MultiplyMatrixAzThreeD(float** MatrixA , float* zThreeD , float aCoeff[] , int mDegree , int dataSize) 
{
	int i ; 
	int m ;
	float sum ;

	
	for( i = 0 ; i < mDegree ; i++ )
	{
		sum = 0.0 ;

		for( m = 0 ; m < dataSize ; m++ )
		{
			sum += MatrixA[i][m]*zThreeD[m] ;
		}

		aCoeff[i] = sum;
	}
}


//***********************************************************************************************************************
// - end -
//************************************************************************************************************************


//********************************************************************************************************************************************************************************
// Least square regression À¸·Î fitting ÇÏ´Â µ¥ ÇÊ¿äÇÑ ÇÔ¼öµé....
//********************************************************************************************************************************************************************************

// ³Ö¾îÁØ °ªÀ» nTimes ¸¸Å­ °öÇÏ´Â ÇÔ¼ö.....
float LeastSquareFns::GetNTimesValue(float value , int nTimes)
{
	
	if( nTimes == 0 )
		return 1.0;
	else
		return value*GetNTimesValue(value , nTimes - 1);
}

//Least square·Î  polynomial fitting ÇÏ´Â ÇÔ¼ö...  fitting ÀÌ ¾ÈµÇ¸é 1À» µÇ¸é 0À» ¸®ÅÏÇÑ´Ù...
//0 ºÎÅÍ ½ÃÀÛÇؼ­ Â÷¼ö +1 ¸¸Å­ ¸Þ¸ð¸®¸¦ Àâ¾Ò´Ù....
// a0 + a1*x + a2*x*x + .....
int LeastSquareFns::GetRegressionCoeffInPoly(float* xPoint , float* yPoint , float aCoeff[] , int mDegree , int nSize)
{
	int i , j , result ;
	float** zMatrix ;
	float** trzMatrix;
	float** trzzMatrix ;
	float** inversetrzzMatrix ;
	float** resultMatrix ;


	zMatrix = (float** ) malloc(nSize*sizeof(float*)) ;
	trzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	trzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	inversetrzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	resultMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;

	  

	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
	}

	for( i = 0 ; i < (mDegree+1) ; i++ )
	{
		trzMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
		trzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;  
		inversetrzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
		resultMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
	}
	
	for( i = 0 ; i < nSize ; i++ )
	{
		for( j = 0 ; j < mDegree + 1 ; j++ )
		{
			zMatrix[i][j] = GetNTimesValue( xPoint[i] , j ) ;		
		}// end of inner for loop 
	}

	// ÀüÄ¡ Çà·ÄÀ» ±¸ÇÑ´Ù.
	invFnsClass.TranseposeMatrix( zMatrix , trzMatrix , nSize , mDegree +1 ) ;
	// ÀüÄ¡Çà·ÄÀ» °öÇÑ´Ù.
	invFnsClass.MultiplyMatrixAB( trzMatrix , zMatrix , trzzMatrix , mDegree+1 , nSize , mDegree + 1) ;
	// ¿ªÇà·ÄÀ» ±¸ÇÑ´Ù.
	result = invFnsClass.CalculatingInverseMatrix( trzzMatrix , mDegree+1 , inversetrzzMatrix ) ;

	if( result == 1 )
	{
		invFnsClass.MultiplyMatrixAB(inversetrzzMatrix , trzMatrix , resultMatrix , mDegree+1 , mDegree+1 ,  nSize) ;
		MultiplyMatrixAzThreeD(resultMatrix , yPoint , aCoeff  , mDegree+1 , nSize) ; 

		for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

		for( i=0 ; i < (mDegree+1) ; i++ )
		{
			free( trzMatrix[i] ) ;
			free( trzzMatrix[i]) ;
			free( inversetrzzMatrix[i] ) ;
			free( resultMatrix[i] ) ;
		}

		free( zMatrix ) ;
		free( trzMatrix ) ;
		free( trzzMatrix) ;
		free( inversetrzzMatrix ) ;
		free( resultMatrix ) ;
		return 1 ;

	}

	for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

	for( i=0 ; i < (mDegree+1) ; i++ )
	{
		free( trzMatrix[i] ) ;
		free( trzzMatrix[i]) ;
		free( inversetrzzMatrix[i] ) ;
		free( resultMatrix[i] ) ;
	}

	free( zMatrix ) ;
	free( trzMatrix ) ;
	free( trzzMatrix) ;
	free( inversetrzzMatrix ) ;
	free( resultMatrix ) ;
	
	return 0;
}


// polynomial fitting ÈÄ ±× ÇÔ¼öÀÇ y°ªÀ» ¾Ë¾Æ³»´Â ÇÔ¼ö...
float LeastSquareFns::GetRegressionFnValue(float aCoeff[] ,float xValue , int mDegree) 
{
	float nTimesxValue = 1.0 ;
	int i ;
	if( mDegree == 0 )
		return aCoeff[0] ;
	else
	{
		for( i = 1 ; i <= mDegree ; i++ )
			nTimesxValue *= xValue ;
		return ( aCoeff[mDegree]* nTimesxValue + GetRegressionFnValue( aCoeff , xValue , mDegree - 1 ) );
	}
}



// Least square ·Î  Æò¸éÀ¸·Î fitting ÇÏ´Â ÇÔ¼ö...  fitting ÀÌ ¾ÈµÇ¸é 1À» µÇ¸é 0À» ¸®ÅÏÇÑ´Ù...
// 0 ºÎÅÍ ½ÃÀÛÇϹǷΠÂ÷¼ö +1 ¸¸Å­ ¸Þ¸ð¸®¸¦ Àâ¾Ò´Ù....
// z = ax + by + c  ÇüÅÂÀÇ Æò¸éÀÌ´Ù...
int LeastSquareFns::GetRegressionCoeffInPlane(float* xPoint , float* yPoint , float* zPoint, float aCoeff[] , int nSize)
{
	int i, result ;
	float** zMatrix ;
	float** trzMatrix;
	float** trzzMatrix ;
	float** inversetrzzMatrix ;
	float** resultMatrix ;
	int mDegree = 2;			// »ó¼ö¸¦ »« ³ª¸ÓÁö Â÷¼ö...

	zMatrix = (float** ) malloc(nSize*sizeof(float*)) ;
	trzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	trzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	inversetrzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	resultMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;

	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
	}

	for( i = 0 ; i < (mDegree+1) ; i++ )
	{
		trzMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
		trzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;  
		inversetrzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
		resultMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
	}
	
	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i][0] = xPoint[i];
		zMatrix[i][1] = yPoint[i];
		zMatrix[i][2] = 1;
	}

	invFnsClass.TranseposeMatrix( zMatrix , trzMatrix , nSize , mDegree +1 ) ;
	invFnsClass.MultiplyMatrixAB( trzMatrix , zMatrix , trzzMatrix , mDegree+1 , nSize , mDegree + 1) ;
	result = invFnsClass.CalculatingInverseMatrix( trzzMatrix , mDegree+1 , inversetrzzMatrix ) ;

	if( result == 1 )
	{
		invFnsClass.MultiplyMatrixAB( inversetrzzMatrix , trzMatrix , resultMatrix , mDegree+1 , mDegree+1 ,  nSize) ;
		MultiplyMatrixAzThreeD( resultMatrix , zPoint , aCoeff, mDegree+1 , nSize ) ; 

		for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

		for( i=0 ; i < (mDegree+1) ; i++ )
		{
			free( trzMatrix[i] ) ;
			free( trzzMatrix[i]) ;
			free( inversetrzzMatrix[i] ) ;
			free( resultMatrix[i] ) ;
		}

		free( zMatrix ) ;
		free( trzMatrix ) ;
		free( trzzMatrix) ;
		free( inversetrzzMatrix ) ;
		free( resultMatrix ) ;
		return 1;

	}

	for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

	for( i=0 ; i < (mDegree+1) ; i++ )
	{
		free( trzMatrix[i] ) ;
		free( trzzMatrix[i]) ;
		free( inversetrzzMatrix[i] ) ;
		free( resultMatrix[i] ) ;
	}

	free( zMatrix ) ;
	free( trzMatrix ) ;
	free( trzzMatrix) ;
	free( inversetrzzMatrix ) ;
	free( resultMatrix ) ;
	
	return 0;
}

//Least square·Î  Æò¸éÀ¸·Î fitting ÇÏ´Â ÇÔ¼ö...  fitting ÀÌ ¾ÈµÇ¸é 1À» µÇ¸é 0À» ¸®ÅÏÇÑ´Ù...
//0 ºÎÅÍ ½ÃÀÛÇϹǷΠÂ÷¼ö +1 ¸¸Å­ ¸Þ¸ð¸®¸¦ Àâ¾Ò´Ù....
// z = ax + by + c  ÇüÅÂÀÇ Æò¸éÀÌ´Ù...
// z = a0 + a1*2*x*y + a2(x*x + y*y) + a2(y*y - x*x)
int LeastSquareFns::GetRegressionCoeffInZernikePoly( float* xPoint , float* yPoint , float* zPoint, float aCoeff[] , int nSize )
{
	int i, result ;
	float** zMatrix ;
	float** trzMatrix;
	float** trzzMatrix ;
	float** inversetrzzMatrix ;
	float** resultMatrix ;
	int mDegree = 3;			// »ó¼ö¸¦ »« ³ª¸ÓÁö Â÷¼ö...

	zMatrix = (float** ) malloc(nSize*sizeof(float*)) ;
	trzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	trzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	inversetrzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	resultMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;

	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
	}

	for( i = 0 ; i < (mDegree+1) ; i++ )
	{
		trzMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
		trzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;  
		inversetrzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
		resultMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
	}
	
	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i][0] = 1;
		zMatrix[i][1] = 2*xPoint[i]*yPoint[i];
		zMatrix[i][2] = 2*(xPoint[i]*xPoint[i] + yPoint[i]*yPoint[i]);
		zMatrix[i][3] = yPoint[i]*yPoint[i] - xPoint[i]*xPoint[i];
//		zMatrix[i][4] = xPoint[i];
//		zMatrix[i][5] = yPoint[i];	
	}

	invFnsClass.TranseposeMatrix( zMatrix , trzMatrix , nSize , mDegree +1 ) ;
	invFnsClass.MultiplyMatrixAB( trzMatrix , zMatrix , trzzMatrix , mDegree+1 , nSize , mDegree + 1) ;
	result = invFnsClass.CalculatingInverseMatrix( trzzMatrix , mDegree+1 , inversetrzzMatrix ) ;

	if( result == 1 )
	{
		invFnsClass.MultiplyMatrixAB( inversetrzzMatrix , trzMatrix , resultMatrix , mDegree+1 , mDegree+1 ,  nSize) ;
		MultiplyMatrixAzThreeD( resultMatrix , zPoint , aCoeff  , mDegree+1 , nSize ) ; 

		for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

		for( i=0 ; i < (mDegree+1) ; i++ )
		{
			free( trzMatrix[i] ) ;
			free( trzzMatrix[i]) ;
			free( inversetrzzMatrix[i] ) ;
			free( resultMatrix[i] ) ;
		}

		free( zMatrix ) ;
		free( trzMatrix ) ;
		free( trzzMatrix) ;
		free( inversetrzzMatrix ) ;
		free( resultMatrix ) ;
		return 0 ;

	}

	for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

	for( i=0 ; i < (mDegree+1) ; i++ )
	{
		free( trzMatrix[i] ) ;
		free( trzzMatrix[i]) ;
		free( inversetrzzMatrix[i] ) ;
		free( resultMatrix[i] ) ;
	}

	free( zMatrix ) ;
	free( trzMatrix ) ;
	free( trzzMatrix) ;
	free( inversetrzzMatrix ) ;
	free( resultMatrix ) ;
	
	return 1 ;
}


//Least square·Î  Æò¸éÀ¸·Î fitting ÇÏ´Â ÇÔ¼ö...  fitting ÀÌ ¾ÈµÇ¸é 1À» µÇ¸é 0À» ¸®ÅÏÇÑ´Ù...
//0 ºÎÅÍ ½ÃÀÛÇϹǷΠÂ÷¼ö +1 ¸¸Å­ ¸Þ¸ð¸®¸¦ Àâ¾Ò´Ù....
// z = ax + by + c  ÇüÅÂÀÇ Æò¸éÀÌ´Ù...
int LeastSquareFns::GetRegressionCoeffInZernikePoly( float* xPoint , float* yPoint , int nSize)
{
	int i, j, result ;
	float** zMatrix ;
	float** trzMatrix;
	float** trzzMatrix ;
	float** inversetrzzMatrix ;
	float** resultMatrix ;
	int mDegree = 1;	// »ó¼ö¸¦ »« ³ª¸ÓÁö Â÷¼ö...  ¹è¿­ÀÇ ¸¶Áö¸· À妽º¸¦ ³ÖÀ¸¸é µÈ´Ù...


	FILE* fp;

	fp = fopen("LeastSquareMatrix.txt", "w");



	zMatrix = (float** ) malloc(nSize*sizeof(float*)) ;
	trzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	trzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	inversetrzzMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;
	resultMatrix = (float** ) malloc((mDegree+1)*sizeof(float*)) ;

	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
	}

	for( i = 0 ; i < (mDegree+1) ; i++ )
	{
		trzMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
		trzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;  
		inversetrzzMatrix[i] = (float* ) malloc((mDegree+1)*sizeof(float)) ;
		resultMatrix[i] = (float* ) malloc(nSize*sizeof(float)) ;
	}
	

/*
	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i][0] = 1;
		zMatrix[i][1] = xPoint[i];
		zMatrix[i][2] = yPoint[i];
		zMatrix[i][3] = 2*xPoint[i]*yPoint[i];
		zMatrix[i][4] = (2*xPoint[i]*xPoint[i] + 2*yPoint[i]*yPoint[i] - 1);
		zMatrix[i][5] = yPoint[i]*yPoint[i] - xPoint[i]*xPoint[i];	
	}
*/
	
	
/*	
	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i][0] = 1;
		zMatrix[i][1] = 2*xPoint[i]*yPoint[i];
		zMatrix[i][2] = 2*(xPoint[i]*xPoint[i] + yPoint[i]*yPoint[i]);
		zMatrix[i][3] = yPoint[i]*yPoint[i] - xPoint[i]*xPoint[i];	
	}
*/

/*
	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i][0] = 2*xPoint[i]*yPoint[i];
		zMatrix[i][1] = (2*xPoint[i]*xPoint[i] + 2*yPoint[i]*yPoint[i]);
		zMatrix[i][2] = yPoint[i]*yPoint[i] - xPoint[i]*xPoint[i];	
	}
*/

	
/*
	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i][0] = 1;
		zMatrix[i][1] = xPoint[i];
		zMatrix[i][2] = yPoint[i];
	}
*/


	for( i = 0 ; i < nSize ; i++ )
	{
		zMatrix[i][0] = 1;
//		zMatrix[i][1] = xPoint[i];
//		zMatrix[i][2] = yPoint[i];
		zMatrix[i][1] = xPoint[i]*xPoint[i] + yPoint[i]*yPoint[i];
	}

	invFnsClass.TranseposeMatrix( zMatrix , trzMatrix , nSize , mDegree +1 ) ;
	invFnsClass.MultiplyMatrixAB( trzMatrix , zMatrix , trzzMatrix , mDegree+1 , nSize , mDegree + 1) ;
	result = invFnsClass.CalculatingInverseMatrix( trzzMatrix , mDegree+1 , inversetrzzMatrix ) ;

	if( result == 1 )
	{
		invFnsClass.MultiplyMatrixAB( inversetrzzMatrix , trzMatrix , resultMatrix , mDegree+1 , mDegree+1 ,  nSize) ;

		
		for(i = 0; i < mDegree + 1; i++)
		{
			for(j = 0; j < nSize; j++)
			{
				fprintf(fp, "%f, ", resultMatrix[i][j]);
			}

			fprintf(fp, "\n");
		}




		fprintf(fp, "integer version\n\n\n\n");
		for(i = 0; i < mDegree + 1; i++)
		{
			for(j = 0; j < nSize; j++)
			{
				fprintf(fp, "%5d, ", (int)(resultMatrix[i][j]*100000.0));
			}

			fprintf(fp, "\n");
		}
		fclose(fp);




		
//		MultiplyMatrixAzThreeD( resultMatrix , zPoint , aCoeff  , mDegree+1 , nSize ) ; 

		for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

		for( i=0 ; i < (mDegree+1) ; i++ )
		{
			free( trzMatrix[i] ) ;
			free( trzzMatrix[i]) ;
			free( inversetrzzMatrix[i] ) ;
			free( resultMatrix[i] ) ;
		}

		free( zMatrix ) ;
		free( trzMatrix ) ;
		free( trzzMatrix) ;
		free( inversetrzzMatrix ) ;
		free( resultMatrix ) ;
		return 0 ;

	}

	for( i = 0 ; i < nSize ; i++ )
			free(zMatrix[i] ) ;

	for( i=0 ; i < (mDegree+1) ; i++ )
	{
		free( trzMatrix[i] ) ;
		free( trzzMatrix[i]) ;
		free( inversetrzzMatrix[i] ) ;
		free( resultMatrix[i] ) ;
	}

	free( zMatrix ) ;
	free( trzMatrix ) ;
	free( trzzMatrix) ;
	free( inversetrzzMatrix ) ;
	free( resultMatrix ) ;
	
	return 1 ;
}

int GetRegressionCoff_CircleFit(float* xps, float *yps, int nData)
{
	// xx, yy, ax, by, c
	return 0;
}

//Least square·Î  polynomial fitting ÇÏ´Â ÇÔ¼ö...  fitting ÀÌ ¾ÈµÇ¸é 1À» µÇ¸é 0À» ¸®ÅÏÇÑ´Ù...
//0 ºÎÅÍ ½ÃÀÛÇؼ­ Â÷¼ö +1 ¸¸Å­ ¸Þ¸ð¸®¸¦ Àâ¾Ò´Ù....
// a0 + a1*x + a2*x*x + .....
int LeastSquareFns::GetRegressionCoeffInCircleFitting(float* xPoint, float* yPoint, int dataNum, float *cx, float *cy, float *radius)
{
	// (x*x + y*y + ax + by + c)¿¡¼­ ÀÌ°ªÀÇ Á¦°öµéÀÇ ÇÕÀ» Á¦ÀÏ ÀÛ°ÔÇÏ´Â a, b, c ¸¦ ±¸ÇÏÀÚ.
	float xxx, yyy, xxy, yyx, xx, yx, xy, yy, x, y;
	int dataCnt, anyCnt;
	float **AMat;
	float **A_Mat;
	float *BMat;
	int result;
	float a, b, c;

	AMat = new float*[3];
	A_Mat = new float*[3];
	BMat = new float[3];

	for(anyCnt = 0; anyCnt < 3; anyCnt++)
	{
		AMat[anyCnt] = new float[3];
		A_Mat[anyCnt] = new float[3];
	}

	xxx = yyy = xxy = yyx = xx = yx = xy = yy = x = y = 0.0;
	
	// data °³¼ö ¸¸Å­ À§ ½ÄÀ» ¹ÌºÐÇÑ °ªµéÀ» ´õÇÏÀÚ..
	for(dataCnt = 0; dataCnt < dataNum; dataCnt++)
	{
		xxx += (xPoint[dataCnt]*xPoint[dataCnt]*xPoint[dataCnt]);
		yyy += (yPoint[dataCnt]*yPoint[dataCnt]*yPoint[dataCnt]);
		xxy += (xPoint[dataCnt]*xPoint[dataCnt]*yPoint[dataCnt]);
		yyx += (yPoint[dataCnt]*yPoint[dataCnt]*xPoint[dataCnt]);
		xx += (xPoint[dataCnt]*xPoint[dataCnt]);
		yx += (yPoint[dataCnt]*xPoint[dataCnt]);
		xy += (xPoint[dataCnt]*yPoint[dataCnt]);
		yy += (yPoint[dataCnt]*yPoint[dataCnt]);
		x += xPoint[dataCnt];
		y += yPoint[dataCnt];
	}

	AMat[0][0] = xx;
	AMat[0][1] = yx;
	AMat[0][2] = x;
	AMat[1][0] = xy;
	AMat[1][1] = yy;
	AMat[1][2] = y;
	AMat[2][0] = x;
	AMat[2][1] = y;
	AMat[2][2] = dataNum;

	BMat[0] = -xxx - yyx;
	BMat[1] = -xxy - yyy;
	BMat[2] = -xx - yy;

	result = invFnsClass.CalculatingInverseMatrix(AMat , 3, A_Mat);

	//ÀÌ °æ¿ì fittingÀÌ ¾ÈµÈ´Ù..
	if(result == 0)
	{
		for(anyCnt = 0; anyCnt < 3; anyCnt++)
		{
			delete [] AMat[anyCnt];
			delete [] A_Mat[anyCnt];
		}

		delete [] AMat;
		delete [] A_Mat;
		delete [] BMat;
		
		return 0;
	}

	a = b = c = 0.0;

	for(anyCnt = 0; anyCnt < 3; anyCnt++)
	{
		a += A_Mat[0][anyCnt]*BMat[anyCnt];
		b += A_Mat[1][anyCnt]*BMat[anyCnt];
		c += A_Mat[2][anyCnt]*BMat[anyCnt];
	}

	*cx = -a/2.0;
	*cy = -b/2.0;
	*radius = sqrt((*cx)*(*cx) + (*cy)*(*cy) - c);

	for(anyCnt = 0; anyCnt < 3; anyCnt++)
	{
		delete [] AMat[anyCnt];
		delete [] A_Mat[anyCnt];
	}

	delete [] AMat;
	delete [] A_Mat;
	delete [] BMat;

	return 1;
}