ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/common/invmat4.c
Revision: 2.3
Committed: Tue Feb 25 02:47:21 2003 UTC (21 years, 2 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R4, rad5R2, rad4R2P2, rad5R0, rad5R1, rad3R7P2, rad3R7P1, rad4R2, rad4R1, rad4R0, rad3R5, rad3R6, rad3R6P1, rad3R8, rad3R9, rad4R2P1, rad5R3, HEAD
Changes since 2.2: +1 -56 lines
Log Message:
Replaced inline copyright notice with #include "copyright.h"

File Contents

# User Rev Content
1 greg 2.1 #ifndef lint
2 greg 2.2 static const char RCSid[] = "$Id$";
3 greg 2.1 #endif
4     /*
5     * invmat4 - computes the inverse of mat into inverse. Returns 1
6     * if there exists an inverse, 0 otherwise. It uses Gaussian Elimination
7     * method.
8     */
9    
10 greg 2.3 #include "copyright.h"
11 greg 2.2
12 greg 2.1 #include "mat4.h"
13    
14    
15     #define ABS(x) ((x)>=0 ? (x) : -(x))
16    
17     #define SWAP(a,b,t) (t=a,a=b,b=t)
18    
19    
20 greg 2.2 int
21     invmat4(inverse, mat)
22 greg 2.1 MAT4 inverse, mat;
23     {
24     MAT4 m4tmp;
25     register int i,j,k;
26     register double temp;
27    
28     copymat4(m4tmp, mat);
29     /* set inverse to identity */
30 greg 2.2 setident4(inverse);
31 greg 2.1
32     for(i = 0; i < 4; i++) {
33     /* Look for row with largest pivot and swap rows */
34     temp = FTINY; j = -1;
35     for(k = i; k < 4; k++)
36     if(ABS(m4tmp[k][i]) > temp) {
37     temp = ABS(m4tmp[k][i]);
38     j = k;
39     }
40     if(j == -1) /* No replacing row -> no inverse */
41     return(0);
42     if (j != i)
43     for(k = 0; k < 4; k++) {
44     SWAP(m4tmp[i][k],m4tmp[j][k],temp);
45     SWAP(inverse[i][k],inverse[j][k],temp);
46     }
47    
48     temp = m4tmp[i][i];
49     for(k = 0; k < 4; k++) {
50     m4tmp[i][k] /= temp;
51     inverse[i][k] /= temp;
52     }
53     for(j = 0; j < 4; j++) {
54     if(j != i) {
55     temp = m4tmp[j][i];
56     for(k = 0; k < 4; k++) {
57     m4tmp[j][k] -= m4tmp[i][k]*temp;
58     inverse[j][k] -= inverse[i][k]*temp;
59     }
60     }
61     }
62     }
63     return(1);
64     }