$search
00001 #include "matrices.h" 00002 00003 void addMatrix2 (float A[2][2], float B[2][2], float C[2][2]){ 00004 /*C= A+B; 00005 * A,B and C need to be the same size (m by n); function does not check*/ 00006 int i, j; 00007 for(i=0; i<2; ++i){ 00008 for(j=0; j<2; ++j){ 00009 C[i][j] = A[i][j] + B[i][j]; 00010 }//inner for 00011 }//outer for 00012 }//add matrices 00013 00014 void addMatrix (matrix A, matrix B, matrix C, int m, int n){ 00015 /*C= A+B; 00016 * A,B and C need to be the same size (m by n); function does not check*/ 00017 int i, j; 00018 for(i=0; i<m; ++i){ 00019 for(j=0; j<n; ++j){ 00020 C[i][j] = A[i][j] + B[i][j]; 00021 }//inner for 00022 }//outer for 00023 }//add matrices 00024 00025 void addVector (float* A, float* B, float* C, int m){ 00026 /*C= A+B; 00027 * A,B and C need to be the same size (m); function does not check*/ 00028 int i; 00029 for(i=0; i<m; ++i){ 00030 C[i] = A[i] + B[i]; 00031 }//for 00032 }//add Vectors 00033 00034 void addVector2(float V[2], float V1[2], float V2[2]) 00035 { 00036 V2[0] = V[0] + V1[0]; 00037 V2[1] = V[1] + V1[1]; 00038 } 00039 00040 void subMatrix2 (float A[2][2], float B[2][2], float C[2][2]){ 00041 /*C= A-B; 00042 * A,B and C need to be the same size (m by n); function does not check*/ 00043 int i, j; 00044 for(i=0; i<2; ++i){ 00045 for(j=0; j<2; ++j){ 00046 C[i][j] = A[i][j] - B[i][j]; 00047 }//inner for 00048 }//outer for 00049 }//sub matrices 00050 00051 00052 void subMatrix (matrix A, matrix B, matrix C, int m, int n){ 00053 /*C= A-B; 00054 * A,B and C need to be the same size (m by n); function does not check*/ 00055 int i, j; 00056 for(i=0; i<m; ++i){ 00057 for(j=0; j<n; ++j){ 00058 C[i][j] = A[i][j] - B[i][j]; 00059 }//inner for 00060 }//outer for 00061 }//sub matrices 00062 00063 void subVector (float* A, float* B, float* C, int m){ 00064 /*C= A-B; 00065 * A,B and C need to be the same size (m); function does not check*/ 00066 int i; 00067 for(i=0; i<m; ++i){ 00068 C[i] = A[i] - B[i]; 00069 }//for 00070 }//add Vectors 00071 00072 void subVector2 (float A[2], float B[2], float C[2]){ 00073 /*C= A-B; 00074 * A,B and C need to be the same size (m); function does not check*/ 00075 int i; 00076 for(i=0; i<2; ++i){ 00077 C[i] = A[i] - B[i]; 00078 }//for 00079 }//add Vectors 00080 00081 00082 void negMatrix(matrix A, int m, int n){ 00083 /*negate matrix*/ 00084 int i,j; 00085 for(i=0; i<m; ++i){ 00086 for(j=0; j<n; ++j){ 00087 A[i][j] = -A[i][j]; 00088 }//inner for 00089 }//outer for 00090 }//negMAtrix 00091 00092 00093 void multMatrix(matrix A, matrix B, matrix C, int Am, int An, int Bn){ 00094 //C = A*B 00095 //A,B and C need to be of correct dimensions, function does not check 00096 int i,j,k; 00097 float sum; 00098 00099 for(i=0; i < Bn; i++){ 00100 for(j=0; j < Am; j++){ 00101 sum=0; 00102 for(k=0; k < An; k++) 00103 sum += A[j][k] * B[k][i]; 00104 C[j][i] = sum; 00105 }//inner for 00106 }//outer for 00107 } 00108 00109 void multMatrix2(float A[2][2], float B[2][2], float C[2][2]){ 00110 //C = A*B 00111 //A,B and C need to be of correct dimensions, function does not check 00112 int i,j,k; 00113 float sum; 00114 00115 for(i=0; i < 2; i++){ 00116 for(j=0; j < 2; j++){ 00117 sum=0; 00118 for(k=0; k < 2; k++) 00119 sum += A[j][k] * B[k][i]; 00120 C[j][i] = sum; 00121 }//inner for 00122 }//outer for 00123 } 00124 void multMatVec(matrix A, float* B, float* C, int Am, int Bn){ 00125 /*C = A*B 00126 *A,B and C need to be of correct dimensions, function does not check */ 00127 00128 00129 int j,k; 00130 float sum; 00131 for(j=0; j < Am; j++){ 00132 sum=0; 00133 for(k=0; k < Bn; k++) 00134 sum += A[j][k] * B[k]; 00135 C[j] = sum; 00136 }//for 00137 }//multMatVec 00138 00139 void multMatVec2(float A[2][2], float B[2], float C[2]){ 00140 /*C = A*B 00141 *A,B and C need to be of correct dimensions, function does not check */ 00142 00143 int j,k; 00144 float sum; 00145 for(j=0; j < 2; j++){ 00146 sum=0; 00147 for(k=0; k < 2; k++) 00148 sum += A[j][k] * B[k]; 00149 C[j] = sum; 00150 }//for 00151 }//multMatVec 00152 00153 00154 void transpose(matrix A, matrix B, int Am, int An){ 00155 /*B = A'*/ 00156 int i,j; 00157 for(i=0; i < Am; i++) 00158 for(j=0; j < An; j++) 00159 B[j][i] = A[i][j]; 00160 00161 } 00162 00163 void invert3(matrix A, matrix B){ 00164 /* B = A^-1 for 3x3 matrices*/ 00165 float det = (A[0][0]*A[1][1]*A[2][2]-A[0][0]*A[1][2]*A[2][1]-A[1][0]*A[0][1]*A[2][2]+A[1][0]*A[0][2]*A[2][1]+A[2][0]*A[0][1]*A[1][2]-A[2][0]*A[0][2]*A[1][1]); 00166 00167 B[0][0]=(A[1][1]*A[2][2]-A[1][2]*A[2][1])/det; 00168 B[0][1]=(-A[0][1]*A[2][2]+A[0][2]*A[2][1])/det; 00169 B[0][2]=(A[0][1]*A[1][2]-A[0][2]*A[1][1])/det; 00170 B[1][0]=(-A[1][0]*A[2][2]+A[1][2]*A[2][0])/det; 00171 B[1][1]=(A[0][0]*A[2][2]-A[0][2]*A[2][0])/det; 00172 B[1][2]=(-A[0][0]*A[1][2]+A[0][2]*A[1][0])/det; 00173 B[2][0]=(A[1][0]*A[2][1]-A[1][1]*A[2][0])/det; 00174 B[2][1]=(-A[0][0]*A[2][1]+A[0][1]*A[2][0])/det; 00175 B[2][2]=(A[0][0]*A[1][1]-A[0][1]*A[1][0])/det; 00176 }//invert3 00177 00178 void invert2(float A[2][2], float B[2][2]){ 00179 /* B = A^-1 for 2x2 matrices*/ 00180 float det = (A[0][0]*A[1][1]-A[0][1]*A[1][0]); 00181 00182 B[0][0]=(A[1][1])/det; 00183 B[0][1]=(-A[0][1])/det; 00184 B[1][0]=(-A[1][0])/det; 00185 B[1][1]=(A[0][0])/det; 00186 }//invert2 00187 00188 float determinant(matrix A,matrix tmp, int n){ 00189 /*returns determinant of n x n matrix*/ 00190 int i,j,k; 00191 float q; 00192 float dtr=1; 00193 /* copy A to tmp*/ 00194 for(i=0; i < n; i++) 00195 for(j=0; j < n; j++) 00196 tmp[i][j] = A[j][i]; 00197 00198 /*compute determinant*/ 00199 for(i=0;i < n;i++){ 00200 for(j=0;j < n;j++){ 00201 q=tmp[j][i]/tmp[i][i]; 00202 for(k=0;k < n;k++){ 00203 if(i==j) break; 00204 tmp[j][k]=(tmp[j][k]-tmp[i][k])*q; 00205 }//innermost for 00206 }//middle for 00207 }//outer for 00208 for(i=0;i < n;i++){ 00209 dtr=dtr*tmp[i][i]; 00210 }//for 00211 return(dtr); 00212 }//det 00213 00214 00215 void getMinor(matrix A, matrix B,int An, int m, int n){ 00216 /* use this function to generate a matrix B that is the mn minor of A; A must be square 00217 * m>=1, n>=1, An is the size of A*/ 00218 int i,j; 00219 int k=0; 00220 for(i=0; i < An; i++){ 00221 if(i!=(m-1)){ 00222 for(j=0; j < (An-1); j++){ 00223 if(j<(n-1)){ 00224 B[k][j] = A[i][j]; 00225 } 00226 else{ 00227 B[k][j] = A[i][j+1]; 00228 } 00229 }//inner for (column) 00230 k++; 00231 }//if 00232 }//outer for 00233 }//get minor 00234 00235 00236 void invt(matrix A, matrix B, matrix tmp1, matrix tmp2, int An){ 00237 /*B = A^-1 ; tmp1 & tmp2 must be at least of size (An-1)x(An-1) 00238 * function does not check if matrix is singular, tried, but even if det(A)==0, it did not escape, 00239 * just reutrned "nan" for some elements, and garbage for others 00240 * */ 00241 int i,j; 00242 float detA = determinant(A,tmp2,An); 00243 for(i=0; i < An; i++){ 00244 for(j=0; j < An; j++){ 00245 getMinor(A,tmp1,An,i+1,j+1); 00246 if(((i+j)%2)==0){ 00247 B[j][i]= determinant(tmp1,tmp2,(An-1))/detA; 00248 } 00249 else{ 00250 B[j][i]= -determinant(tmp1,tmp2,(An-1))/detA; 00251 } 00252 //problem if an element is zero (returns "nan" for zero element) 00253 }//inner for 00254 }//outer for 00255 }//invert 00256 00257 float BhQhBhT(matrix Qh, float* Bh, float* tmp){ 00258 //function returns the value for Bh*Qh*Bh', where Qh is an 3x3 matrix, and Bh is a vector of length 3 00259 00260 00261 00262 //Bh*Qh 00263 int i,j; 00264 float sum; 00265 for(i=0; i < 3; i++){ 00266 sum=0; 00267 for(j=0; j < 3; j++) 00268 sum += Bh[j]*Qh[j][i]; 00269 tmp[i] = sum; 00270 }//outer for 00271 00272 //tmp*BhT 00273 sum = 0; 00274 for(i=0; i < 3; i++) 00275 sum += tmp[i]*Bh[i]; 00276 00277 00278 00279 return sum; 00280 } 00281 00282 void VectdivSc(float* V, float scl){ 00283 //function returns the value of Vector/scalar where vector is 3x1 00284 00285 int i; 00286 for(i=0; i < 3; i++){ 00287 V[i] = V[i]/scl; 00288 }// for 00289 } 00290 void VectmultSc2(float V[2], float scl, float V1[2]) 00291 { 00292 V1[0] = V[0] * scl; 00293 V1[1] = V[1] * scl; 00294 }