#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;
|
}
|