home *** CD-ROM | disk | FTP | other *** search
/ gondwana.ecr.mu.oz.au/pub/ / Graphics.tar / Graphics / voglw.zip / minv4.c < prev    next >
Text File  |  1997-02-13  |  1KB  |  61 lines
  1. /*
  2.  * minv4
  3.  *
  4.  *    find the inverse of the 4 by 4 matrix b using gausian elimination
  5.  * and return it in a.
  6.  * 
  7.  *    (We don't actually use this yet in VOGL - maybe one day).
  8.  */
  9. minv4(a, b)
  10.     Matrix    a, b;
  11. {
  12.     float    val, val2;
  13.     int    i, j, k, ind;
  14.     Matrix    tmp;
  15.  
  16.     identmatrix(a);
  17.  
  18.     copymatrix(tmp, b);
  19.  
  20.     for (i = 0; i != 4; i++) {
  21.  
  22.         val = tmp[i][i];        /* find pivot */
  23.         ind = i;
  24.         for (j = i + 1; j != 4; j++) {
  25.             if (fabs(tmp[j][i]) > fabs(val)) {
  26.                 ind = j;
  27.                 val = tmp[j][i];
  28.             }
  29.         }
  30.  
  31.         if (ind != i) {            /* swap columns */
  32.             for (j = 0; j != 4; j++) {
  33.                 val2 = a[i][j];
  34.                 a[i][j] = a[ind][j];
  35.                 a[ind][j] = val2;
  36.                 val2 = tmp[i][j];
  37.                 tmp[i][j] = tmp[ind][j];
  38.                 tmp[ind][j] = val2;
  39.             }
  40.         }
  41.  
  42.         if (val == 0.0)
  43.             fatal("art: singular matrix in minv4.\n");
  44.  
  45.         for (j = 0; j != 4; j++) {
  46.             tmp[i][j] /= val;
  47.             a[i][j] /= val;
  48.         }
  49.  
  50.         for (j = 0; j != 4; j++) {    /* eliminate column */
  51.             if (j == i)
  52.                 continue;
  53.             val = tmp[j][i];
  54.             for (k = 0; k != 4; k++) {
  55.                 tmp[j][k] -= tmp[i][k] * val;
  56.                 a[j][k] -= a[i][k] * val;
  57.             }
  58.         }
  59.     }
  60. }
  61.