I'm trying to fit a parabola to a segment of data using the two end points and a point near the center.
I've decided to use Cramer's rule because I think that it's much
easier to only worry with determinants instead of having to find
the inverse of the matrix. I was wondering if there is a better or
faster way to go about this.
Also, since my determinants are only for 3x3 matrices would it be
a better idea to simply crank out the equation:
I realize that readability is in favor of the matrix_det function, but I'm concerned with theCode:double m00 = array[0][0]; double m01 = array[0][1]; double m02 = array[0][2]; double m10 = array[1][0]; double m11 = array[1][1]; double m12 = array[1][2]; double m20 = array[2][0]; double m21 = array[2][1]; double m22 = array[2][2]; double determinant = m00*m11*m22 + m01*m12*m20 + m02*m10*m21 - m02*m11*m20 - m01*m10*m22 - m00*m12*m21;
time factor. This code needs to be fairly efficient. This particular loop of code is
executed anywhere from 100 to 1000 times per minute.
My concerns stem from stepping through the code and my inexperience with pointer
objects and memory control. I understand how the array of pointers is the matrix
structure that is passed to the function, which is a similar version of Prelude's determinant
function.
When I step through the code it jumps through all kinds of strange code in malloc.c,
dbgheap.c, new.cpp, memset.asm . . . would it be any less time consuming to just create
an array and carefully create a mapping of a 2d array to a 1d array?
Ex: 2D[0][0]->1D[0] . . . 2D[2][2]->1D[9].
(Or is this concern pointless) I think that the readability will be even more convoluted,
but from a speed standpoint is this really any faster?
![]()
Here's the current version of the code:
Code://parameters[3] is set at 10 for now if ( consecutiveFlags[c] > parameters[3] ) { //*******Set up for solving the system of equations double **ARRAY, **ARRAY2; ARRAY = new double *[Q]; //Q defines the square matrix dimensions ARRAY2 = new double *[Q]; double YDATA[3]; for ( i = 0; i < Q; i++ ) { ARRAY[i] = new double[Q]; ARRAY2[i] = new double[Q]; }//end set up storage //set up for first row of matrix //grab pt @ beginning of segment ARRAY[0][0] = zlen[c]*zlen[c]; ARRAY[0][1] = zlen[c]; ARRAY[0][2] = 1.0; YDATA[0] = yadjusted[c]; //set up for second row of matrix //grab pt near middle of k = c + consecutiveFlags[c]/2-1; //failed segment ARRAY[1][0] = zlen[k]*zlen[k]; ARRAY[1][1] = zlen[k]; ARRAY[1][2] = 1.0; YDATA[1] = yadjusted[k]; //set up for third row of matrix //grab pt @ end of failed segment k = c + consecutiveFlags[c]-1; ARRAY[2][0] = zlen[k]*zlen[k]; ARRAY[2][1] = zlen[k]; ARRAY[2][2] = 1.0; YDATA[2] = yadjusted[k]; /*To find the coefficients of the parabola for the segment correction we must solve: [zlen1^2 zlen1 1][A] [y_front ] [zlen2^2 zlen2 1][B]= [y_middle] or ARRAY * Coefficients = YDATA [zlen3^2 zlen3 1][C] [y_end ] (Use Cramer's Rule) Step 1)Check for det of ARRAY, if zero break out of loop b/c division by zero Step 2)Set up ARRAY2 with the original ARRAY matrix data but substituting YDATA transpose for each column corresponding to the position of the desired Coefficient and then calculating the determinant. Step 3)Coefficients will be equal to to ratio of the determinant of ARRAY2 to the determinant of ARRAY. */ //Step 1: double detARRAY = matrix_det ( ARRAY, 3 ); if ( detARRAY == 0 ) break; //error else ;//nothing //continue //Step 2: double Coefficients[3]; for ( int m = 0; m < Q; m++ ) //m iterates the Coefficients { for ( j = 0; j < Q; j++ ) //j iterates the columns { for( i = 0; i < Q; i++ ) //i iterates the rows { if ( j==m ) ARRAY2[i][j] = YDATA[i]; //fill appropriate column with ydata else ARRAY2[i][j] = ARRAY[i][j];//copy original array in other columns } }//end for populate ARRAY2 double detARRAY2 = matrix_det ( ARRAY2, 3 ); Coefficients[m] = detARRAY2/detARRAY; }//end for populate Coefficients //created an array object to model a matrix //must delete it to avoid a memory leak for ( i = 0; i < Q; i++ ) { delete [] ARRAY[i]; delete [] ARRAY2[i]; }//delete storage containers delete [] ARRAY; delete [] ARRAY2; }//end if find the parabola for the segment double matrix_det ( double **in_matrix, int n ) { int i, j, k; double **matrix; double det = 1; matrix = new double *[n]; for ( i = 0; i < n; i++ ) //setting up pointers for rows i, [0,n-1] matrix[i] = new double[n]; for ( i = 0; i < n; i++ ) { for ( j = 0; j < n; j++ ) matrix[i][j] = in_matrix[i][j]; //setting up pointers for columns j, [0,n-1] }//end pointer setup for ( k = 0; k < n; k++ ) { if ( matrix[k][k] == 0 ) //check for zeros along the diagonal { bool ok = false; for ( j = k; j < n; j++ ) //ensures that if a zero exists along the diagonal { //at least one value in that column has to be non-zero if ( matrix[j][k] != 0 ) ok = true; }//end for check if det will be zero if ( !ok ) //if column of zeros exist, then the matrix has det=0 return 0; for ( i = k; i < n; i++ ) { double temp; temp = matrix[i][k]; matrix[i][k] = matrix[i][j]; matrix[i][j] = temp; //std::swap ( matrix[i][j], matrix[i][k] ); }//end for swap column k for column j, the last column det = -det; }//end if a zero occurs along diagonal det *= matrix[k][k]; if ( k + 1 < n ) { for ( i = k + 1; i < n; i++ ) { for ( j = k + 1; j < n; j++ ) { matrix[i][j] = matrix[i][j] - matrix[i][k] * matrix[k][j] / matrix[k][k]; }//end for }//end for }//end if }//end for k's for ( i = 0; i < n; i++ ) delete [] matrix[i]; delete [] matrix; return det;



LinkBack URL
About LinkBacks



