Subversion Repositories gelsvn

Rev

Rev 5 | Rev 113 | Go to most recent revision | Details | Compare with Previous | Last modification | View Log | RSS feed

Rev Author Line No. Line
2 bj 1
#include "ArithSqMat3x3Float.h"
2
#include "Mat3x3f.h"
3
 
4
namespace CGLA {
5
 
6
	template<class V, class M>
7
	M invert(const ArithSqMat3x3Float<V,M>& _a)
8
	{
9
		M a(_a[0], _a[1], _a[2]);
10
		M	b;
11
		b.identity();
12
 
45 jab 13
		unsigned int  i, j, i1;
2 bj 14
 
15
		for (j=0; j<3; j++) 
16
			{   
17
				i1 = j;                 // Row with largest pivot candidate
18
				for (i=j+1; i<3; i++)
19
					if (fabs(a[i][j]) > fabs(a[i1][j]))
20
						i1 = i;
21
 
22
				// Swap rows i1 and j in a and b to put pivot on diagonal
23
				V a_tmp = a[i1];
24
				a[i1] = a[j];
25
				a[j]  = a_tmp;
26
 
27
				V b_tmp = b[i1];
28
				b[i1] = b[j];
29
				b[j]  = b_tmp;
30
 
31
				// Scale row j to have a unit diagonal
5 jab 32
				if (a[j][j] == 0.0)
2 bj 33
					throw(Mat3x3fSingular("Tried to invert Singular matrix"));
34
 
35
				b[j] /= a[j][j];
36
				a[j] /= a[j][j];
37
 
38
				// Eliminate off-diagonal elems in col j of a, doing identical ops to b
39
				for (i=0; i<3; i++)
40
					if (i!=j) 
41
						{
42
							b[i] -= a[i][j] * b[j];
43
							a[i] -= a[i][j] * a[j];
44
						}
45
			}
46
		return b;
47
	}                                                                               
48
 
49
	template class ArithSqMat3x3Float<Vec3f,Mat3x3f>;
50
	template Mat3x3f invert(const ArithSqMat3x3Float<Vec3f,Mat3x3f>&);
51
}