sglMat4.hpp Source File

00001 /*****************************************************************************
00002  * SGL: A Scene Graph Library
00003  *
00004  * Copyright (C) 1997-2001  Scott McMillan   All Rights Reserved.
00005  *
00006  * This library is free software; you can redistribute it and/or
00007  * modify it under the terms of the GNU Library General Public
00008  * License as published by the Free Software Foundation; either
00009  * version 2 of the License, or (at your option) any later version.
00010  *
00011  * This library is distributed in the hope that it will be useful,
00012  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00013  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014  * Library General Public License for more details.
00015  *
00016  * You should have received a copy of the GNU Library General Public
00017  * License along with this library; if not, write to the Free
00018  * Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
00019  *****************************************************************************
00020  *     File: sglMat4.hpp
00021  *  Created: 8 March 1999
00022  *  Summary: 4x4 templated matrix class
00023  *****************************************************************************/
00024 
00025 #ifndef __SGL_MAT4_HPP
00026 #define __SGL_MAT4_HPP
00027 
00028 #include <sgl.h>
00029 #include <string.h> // memcpy
00030 #include <sglVector.hpp>
00031 #include <sglVec4.hpp>
00032 #include <sglVec3.hpp>
00033 
00034 // ---------------------------------------------------------------------------
00035 template <class T>
00036 class sglMat4
00037 {
00038 public:
00039     // constructors
00040    sglMat4() {};
00041 
00042    sglMat4(T m00, T m01, T m02, T m03,
00043            T m10, T m11, T m12, T m13,
00044            T m20, T m21, T m22, T m23,
00045            T m30, T m31, T m32, T m33)
00046       {
00047          set(m00, m01, m02, m03, m10, m11, m12, m13,
00048              m20, m21, m22, m23, m30, m31, m32, m33);
00049       }
00050 
00051    template <class S>
00052    sglMat4(const sglMat4<S>& mat)
00053       {
00054          set(mat[0][0], mat[0][1], mat[0][2], mat[0][3],
00055              mat[1][0], mat[1][1], mat[1][2], mat[1][3],
00056              mat[2][0], mat[2][1], mat[2][2], mat[2][3],
00057              mat[3][0], mat[3][1], mat[3][2], mat[3][3]);
00058       }
00059 
00060    template <class S>
00061    sglMat4(S mat[4][4])
00062       {
00063          set(mat[0][0], mat[0][1], mat[0][2], mat[0][3],
00064              mat[1][0], mat[1][1], mat[1][2], mat[1][3],
00065              mat[2][0], mat[2][1], mat[2][2], mat[2][3],
00066              mat[3][0], mat[3][1], mat[3][2], mat[3][3]);
00067       }
00068 
00069    // destructor
00070    ~sglMat4() {};  // not virtual, do not subclass
00071 
00072    inline T *getPtr() const { return (T *)m_matrix; }
00073 
00074    // indexing operators
00075    inline sglVec4<T>& operator[](unsigned int i)
00076       {
00077          return *(sglVec4<T> *) m_matrix[i];
00078       }
00079    inline const sglVec4<T>& operator[](unsigned int i) const
00080       {
00081          return *(sglVec4<T> *) m_matrix[i];
00082       }
00083 
00084    // set/get
00085    template <class S>
00086    inline void set(S m00, S m01, S m02, S m03,
00087                    S m10, S m11, S m12, S m13,
00088                    S m20, S m21, S m22, S m23,
00089                    S m30, S m31, S m32, S m33)
00090       {
00091          m_matrix[0][0] = (T)m00; m_matrix[0][1] = (T)m01;
00092          m_matrix[0][2] = (T)m02; m_matrix[0][3] = (T)m03;
00093          m_matrix[1][0] = (T)m10; m_matrix[1][1] = (T)m11;
00094          m_matrix[1][2] = (T)m12; m_matrix[1][3] = (T)m13;
00095          m_matrix[2][0] = (T)m20; m_matrix[2][1] = (T)m21;
00096          m_matrix[2][2] = (T)m22; m_matrix[2][3] = (T)m23;
00097          m_matrix[3][0] = (T)m30; m_matrix[3][1] = (T)m31;
00098          m_matrix[3][2] = (T)m32; m_matrix[3][3] = (T)m33;
00099       }
00100 
00101    template <class S>
00102    inline void setRow(const sglVec4<S>& vec, unsigned int i)
00103       {
00104          m_matrix[i][0] = (T)vec[0]; m_matrix[i][1] = (T)vec[1];
00105          m_matrix[i][2] = (T)vec[2]; m_matrix[i][3] = (T)vec[3];
00106       }
00107    template <class S>
00108    inline void setCol(const sglVec4<S>& vec, unsigned int i)
00109       {
00110          m_matrix[0][i] = (T)vec[0]; m_matrix[1][i] = (T)vec[1];
00111          m_matrix[2][i] = (T)vec[2]; m_matrix[3][i] = (T)vec[3];
00112       }
00113 
00114    template <class S>
00115    inline void getRow(sglVec4<S>& vec, unsigned int i)
00116       {
00117          vec[0] = (S)m_matrix[i][0]; vec[1] = (S)m_matrix[i][1];
00118          vec[2] = (S)m_matrix[i][2]; vec[3] = (S)m_matrix[i][3];
00119       }
00120 
00121    template <class S>
00122    inline void getCol(sglVec4<S>& vec, unsigned int i)
00123       {
00124          vec[0] = (S)m_matrix[0][i]; vec[1] = (S)m_matrix[1][i];
00125          vec[2] = (S)m_matrix[2][i]; vec[3] = (S)m_matrix[3][i];
00126       }
00127 
00128    // assignment operators
00129    template <class S>
00130    inline sglMat4& operator=(const sglMat4<S>& mat)
00131       {
00132          set(mat[0][0], mat[0][1], mat[0][2], mat[0][3],
00133              mat[1][0], mat[1][1], mat[1][2], mat[1][3],
00134              mat[2][0], mat[2][1], mat[2][2], mat[2][3],
00135              mat[3][0], mat[3][1], mat[3][2], mat[3][3]);
00136          return *this;
00137       }
00138 
00139    // boolean operators
00140    template <class S>
00141    inline bool almostEqual(const sglMat4<S>& mat, double eps) const
00142       {
00143          return ((*this)[0].almostEqual(mat[0], eps) &&
00144                  (*this)[1].almostEqual(mat[1], eps) &&
00145                  (*this)[2].almostEqual(mat[2], eps) &&
00146                  (*this)[3].almostEqual(mat[3], eps));
00147       }
00148 
00149    template <class S>
00150    inline bool operator==(const sglMat4<S>& mat) const
00151       {
00152          return ((*this)[0] == mat[0] && (*this)[1] == mat[1] &&
00153                  (*this)[2] == mat[2] && (*this)[3] == mat[3]);
00154       }
00155 
00156    template <class S>
00157    inline bool operator!=(const sglMat4<S>& mat) const
00158       {
00159          return ((*this)[0] != mat[0] || (*this)[1] != mat[1] ||
00160                  (*this)[2] != mat[2] || (*this)[3] != mat[3]);
00161       }
00162 
00163    // math operators
00164    template <class S>
00165    inline sglMat4<T>& operator+=(const sglMat4<S>& mat)
00166       {
00167          (*this)[0] += mat[0]; (*this)[1] += mat[1];
00168          (*this)[2] += mat[2]; (*this)[3] += mat[3];
00169          return *this;
00170       }
00171    template <class S>
00172    inline sglMat4<T>& operator-=(const sglMat4<S>& mat)
00173       {
00174          (*this)[0] -= mat[0]; (*this)[1] -= mat[1];
00175          (*this)[2] -= mat[2]; (*this)[3] -= mat[3];
00176          return *this;
00177       }
00178    template <class S>
00179    inline sglMat4<T>& operator*=(const sglMat4<S>& mat)
00180       {
00181          mul(*this, mat);   return *this;
00182       }
00183 
00184    inline sglMat4<T>& operator*=(const T s)
00185       {
00186          scale(*this, s, s, s, s);   return *this;
00187       }
00188 
00189    template <class S>
00190    inline const sglMat4<T> operator+(const sglMat4<S> &rhs) const
00191       {
00192          return sglMat4(m_matrix[0][0] + rhs[0][0],
00193                         m_matrix[0][1] + rhs[0][1],
00194                         m_matrix[0][2] + rhs[0][2],
00195                         m_matrix[0][3] + rhs[0][3],
00196 
00197                         m_matrix[1][0] + rhs[1][0],
00198                         m_matrix[1][1] + rhs[1][1],
00199                         m_matrix[1][2] + rhs[1][2],
00200                         m_matrix[1][3] + rhs[1][3],
00201 
00202                         m_matrix[2][0] + rhs[2][0],
00203                         m_matrix[2][1] + rhs[2][1],
00204                         m_matrix[2][2] + rhs[2][2],
00205                         m_matrix[2][3] + rhs[2][3],
00206 
00207                         m_matrix[3][0] + rhs[3][0],
00208                         m_matrix[3][1] + rhs[3][1],
00209                         m_matrix[3][2] + rhs[3][2],
00210                         m_matrix[3][3] + rhs[3][3]);
00211       }
00212    template <class S>
00213    inline const sglMat4<T> operator-(const sglMat4<S> &rhs) const
00214       {
00215          return sglMat4(m_matrix[0][0] - rhs[0][0],
00216                         m_matrix[0][1] - rhs[0][1],
00217                         m_matrix[0][2] - rhs[0][2],
00218                         m_matrix[0][3] - rhs[0][3],
00219 
00220                         m_matrix[1][0] - rhs[1][0],
00221                         m_matrix[1][1] - rhs[1][1],
00222                         m_matrix[1][2] - rhs[1][2],
00223                         m_matrix[1][3] - rhs[1][3],
00224 
00225                         m_matrix[2][0] - rhs[2][0],
00226                         m_matrix[2][1] - rhs[2][1],
00227                         m_matrix[2][2] - rhs[2][2],
00228                         m_matrix[2][3] - rhs[2][3],
00229 
00230                         m_matrix[3][0] - rhs[3][0],
00231                         m_matrix[3][1] - rhs[3][1],
00232                         m_matrix[3][2] - rhs[3][2],
00233                         m_matrix[3][3] - rhs[3][3]);
00234       }
00235    template <class S>
00236    inline const sglMat4<T> operator*(const sglMat4<S> &rhs) const
00237       {
00238          sglMat4<T> tmp;
00239          for (unsigned int i=0; i<4; i++)
00240          {
00241             for (unsigned int j=0; j<4; j++)
00242             {
00243                tmp[i][j] = (m_matrix[i][0]*rhs[0][j] +
00244                             m_matrix[i][1]*rhs[1][j] +
00245                             m_matrix[i][2]*rhs[2][j] +
00246                             m_matrix[i][3]*rhs[3][j]);
00247             }
00248          }
00249          return tmp;
00250       }
00251 
00252    inline const sglMat4<T> operator*(T scalar) const
00253       {
00254          return sglMat4(m_matrix[0][0]*scalar,
00255                         m_matrix[0][1]*scalar,
00256                         m_matrix[0][2]*scalar,
00257                         m_matrix[0][3]*scalar,
00258 
00259                         m_matrix[1][0]*scalar,
00260                         m_matrix[1][1]*scalar,
00261                         m_matrix[1][2]*scalar,
00262                         m_matrix[1][3]*scalar,
00263 
00264                         m_matrix[2][0]*scalar,
00265                         m_matrix[2][1]*scalar,
00266                         m_matrix[2][2]*scalar,
00267                         m_matrix[2][3]*scalar,
00268 
00269                         m_matrix[3][0]*scalar,
00270                         m_matrix[3][1]*scalar,
00271                         m_matrix[3][2]*scalar,
00272                         m_matrix[3][3]*scalar);
00273       }
00274 
00275    // explicit math operations (for backward compatibility)
00276    template <class S1, class S2>
00277    inline void add(const sglMat4<S1>& mat1, const sglMat4<S2>& mat2)
00278       {
00279          (*this)[0].add(mat1[0],mat2[0]);
00280          (*this)[1].add(mat1[1],mat2[1]);
00281          (*this)[2].add(mat1[2],mat2[2]);
00282          (*this)[3].add(mat1[3],mat2[3]);
00283       }
00284    template <class S1, class S2>
00285    inline void sub(const sglMat4<S1>& mat1, const sglMat4<S2>& mat2)
00286       {
00287          (*this)[0].sub(mat1[0],mat2[0]);
00288          (*this)[1].sub(mat1[1],mat2[1]);
00289          (*this)[2].sub(mat1[2],mat2[2]);
00290          (*this)[3].sub(mat1[3],mat2[3]);
00291       }
00292    template <class S>
00293    inline void scale(const sglMat4<S>& mat, T sx, T sy, T sz, T sw=1)
00294       {
00295          (*this)[0].scale(mat[0], sx);
00296          (*this)[1].scale(mat[1], sy);
00297          (*this)[2].scale(mat[2], sz);
00298          // scale does not apply to the translation line
00299          (*this)[3].scale(mat[3], sw);
00300       }
00301    template <class S>
00302    inline void translate(const sglMat4<S>& mat, T tx, T ty, T tz)
00303       {
00304          *this = mat;
00305          for (unsigned int i=0; i<3; i++)
00306          {
00307             (*this)[3][i] +=
00308                mat[0][i]*tx +
00309                mat[1][i]*ty +
00310                mat[2][i]*tz;
00311          }
00312       }
00313 
00314 // matrix product
00315    template <class S1, class S2>
00316    void mul(const sglMat4<S1>& mat1, const sglMat4<S2>& mat2)
00317       {
00318          T t[4][4];
00319          for (unsigned int i=0; i<4; i++)
00320             for (unsigned int j=0; j<4; j++)
00321                t[i][j] = (mat1[i][0]*mat2[0][j] +
00322                           mat1[i][1]*mat2[1][j] +
00323                           mat1[i][2]*mat2[2][j] +
00324                           mat1[i][3]*mat2[3][j]);
00325 
00326          memcpy(m_matrix, t, sizeof(m_matrix));
00327       }
00328 
00329    // methods to build a matrix
00330    inline void buildIdentity()
00331       {
00332          set(1, 0, 0, 0,   0, 1, 0, 0,   0, 0, 1, 0,   0, 0, 0, 1);
00333       }
00334 
00335    inline void buildScale(T sx, T sy, T sz, T sw = 1)
00336       {
00337          buildIdentity();
00338          m_matrix[0][0] = sx; m_matrix[1][1] = sy; m_matrix[2][2] = sz;
00339          m_matrix[3][3] = sw;
00340       }
00341 
00342    inline void buildTranslation(T tx, T ty, T tz)
00343       {
00344          buildIdentity();
00345          m_matrix[3][0] = tx; m_matrix[3][1] = ty; m_matrix[3][2] = tz;
00346       }
00347 
00348    void buildRotation(float x, float y, float z, double angle)  // radians
00349       {
00350          double mag = x*x + y*y + z*z;
00351 
00352          if ((angle == 0.0) || (mag == 0.0))
00353          {
00354             buildIdentity();
00355             return;
00356          }
00357 
00358          mag = (T)sqrt( mag );
00359          x /= mag;  y /= mag;  z /= mag;
00360 
00361          double s = sin(angle), c = cos(angle);
00362          double a = 1 - c;
00363 
00364          double xx = x*x, yy = y*y, zz = z*z;
00365          double xy = x*y, yz = y*z, zx = z*x;
00366          double sx = x*s, sy = y*s, sz = z*s;
00367 
00368          m_matrix[0][0] = (T)(a*xx + c);
00369          m_matrix[0][1] = (T)(a*xy + sz);
00370          m_matrix[0][2] = (T)(a*zx - sy);
00371          m_matrix[0][3] = 0;
00372 
00373          m_matrix[1][0] = (T)(a*xy - sz);
00374          m_matrix[1][1] = (T)(a*yy + c);
00375          m_matrix[1][2] = (T)(a*yz + sx);
00376          m_matrix[1][3] = 0;
00377 
00378          m_matrix[2][0] = (T)(a*zx + sy);
00379          m_matrix[2][1] = (T)(a*yz - sx);
00380          m_matrix[2][2] = (T)(a*zz + c);
00381          m_matrix[2][3] = 0;
00382 
00383          m_matrix[3][0] = m_matrix[3][1] = m_matrix[3][2] = 0;
00384          m_matrix[3][3] = 1;
00385       }
00386 
00392    void buildRotation(double h, double p, double r)  // radians
00393       {
00394          T sr = sin( p );
00395          T cr = cos( p );
00396 
00397          T sh = sin( r );
00398          T ch = cos( r );
00399 
00400          T sp = sin( h );
00401          T cp = cos( h );
00402 
00403          m_matrix[0][0] = -sp*sr*sh + cp*ch;
00404          m_matrix[0][1] =  cp*sr*sh + sp*ch;
00405          m_matrix[0][2] = -cr*sh;
00406          m_matrix[0][3] = 0;
00407 
00408          m_matrix[1][0] = -sp*cr;
00409          m_matrix[1][1] =  cp*cr;
00410          m_matrix[1][2] =  sr;
00411          m_matrix[1][3] =  0;
00412 
00413          m_matrix[2][0] =  sp*sr*ch + cp*sh;
00414          m_matrix[2][1] = -cp*sr*ch + sp*sh;
00415          m_matrix[2][2] =  cr*ch;
00416          m_matrix[2][3] =  0;
00417 
00418          m_matrix[3][0] = m_matrix[3][1] = m_matrix[3][2] = 0;
00419          m_matrix[3][3] = 1;
00420       }
00421 
00422    // other operations
00423    inline void transpose(const sglMat4& mat)
00424       {
00425          // diagonal first
00426          m_matrix[0][0] = mat[0][0];
00427          m_matrix[1][1] = mat[1][1];
00428          m_matrix[2][2] = mat[2][2];
00429          m_matrix[3][3] = mat[3][3];
00430          // then upper triangle
00431          T tmp;
00432          tmp = mat[0][1]; m_matrix[0][1] = mat[1][0]; m_matrix[1][0] = tmp;
00433          tmp = mat[0][2]; m_matrix[0][2] = mat[2][0]; m_matrix[2][0] = tmp;
00434          tmp = mat[0][3]; m_matrix[0][3] = mat[3][0]; m_matrix[3][0] = tmp;
00435          tmp = mat[1][2]; m_matrix[1][2] = mat[2][1]; m_matrix[2][1] = tmp;
00436          tmp = mat[1][3]; m_matrix[1][3] = mat[3][1]; m_matrix[3][1] = tmp;
00437          tmp = mat[2][3]; m_matrix[2][3] = mat[3][2]; m_matrix[3][2] = tmp;
00438       }
00439 
00440 // cofactors, determinants, inverse
00441    inline T cofactor(unsigned int r, unsigned int c) const
00442       {
00443          static const unsigned int index[4][3] = { {1,2,3}, {0,2,3},
00444                                                    {0,1,3}, {0,1,2} };
00445          unsigned int i0=index[r][0], i1=index[r][1], i2=index[r][2];
00446          unsigned int j0=index[c][0], j1=index[c][1], j2=index[c][2];
00447 #define __SUBDET2(r,c,i0,i1,j0,j1) (((r+c)&1) ? -1.f: 1.f) * \
00448               (m_matrix[i0][j0]*m_matrix[i1][j1] - \
00449                m_matrix[i0][j1]*m_matrix[i1][j0])
00450       return (m_matrix[i0][j0] * __SUBDET2(i0,j0,i1,i2,j1,j2)+
00451               m_matrix[i0][j1] * __SUBDET2(i0,j1,i1,i2,j0,j2)+
00452               m_matrix[i0][j2] * __SUBDET2(i0,j2,i1,i2,j0,j1));
00453 #undef __SUBDET2
00454    }
00455 
00456 // reducing methods which result is a T
00457    inline T det() const
00458       {
00459          return (m_matrix[0][0] * cofactor(0,0) +
00460                  m_matrix[0][1] * cofactor(0,1) +
00461                  m_matrix[0][2] * cofactor(0,2) +
00462                  m_matrix[0][3] * cofactor(0,3));
00463       }
00464 
00465    template <class S>
00466    void adjoint(const sglMat4<S>& mat)
00467       {
00468          T t[4][4];
00469          for (unsigned int i=0; i<4; i++)
00470             for (unsigned int j=0; j<4; j++)
00471                t[i][j] = mat.cofactor(i,j);
00472          set(t[0][0], t[0][1], t[0][2], t[0][3],
00473              t[1][0], t[1][1], t[1][2], t[1][3],
00474              t[2][0], t[2][1], t[2][2], t[2][3],
00475              t[3][0], t[3][1], t[3][2], t[3][3]);
00476       }
00477 
00478    template <class S>
00479    bool inverse(const sglMat4<S>& mat)
00480       {
00481          unsigned int indxc[4], indxr[4], ipiv[4];
00482          unsigned int i,j,k,l,ll;
00483          unsigned int icol = 0;
00484          unsigned int irow = 0;
00485          T temp, pivinv, dum, big;
00486 
00487          // copy in place this may be unnecessary
00488          *this = mat;
00489 
00490          for (j=0; j<4; j++) ipiv[j]=0;
00491 
00492          for(i=0;i<4;i++)
00493          {
00494             big=(T)0.0;
00495             for (j=0; j<4; j++)
00496                if (ipiv[j] != 1)
00497                   for (k=0; k<4; k++)
00498                   {
00499                      if (ipiv[k] == 0)
00500                      {
00501                         if (SGL_ABS(m_matrix[j][k]) >= big)
00502                         {
00503                            big = SGL_ABS(m_matrix[j][k]);
00504                            irow=j;
00505                            icol=k;
00506                         }
00507                      }
00508                      else if (ipiv[k] > 1)
00509                         return false;
00510                   }
00511             ++(ipiv[icol]);
00512             if (irow != icol)
00513                for (l=0; l<4; l++) SGL_SWAP(m_matrix[irow][l],
00514                                             m_matrix[icol][l],
00515                                             temp);
00516 
00517             indxr[i]=irow;
00518             indxc[i]=icol;
00519             if (m_matrix[icol][icol] == 0)
00520                return false;
00521 
00522             pivinv = 1.0/m_matrix[icol][icol];
00523             m_matrix[icol][icol] = 1;
00524             for (l=0; l<4; l++) m_matrix[icol][l] *= pivinv;
00525             for (ll=0; ll<4; ll++)
00526                if (ll != icol)
00527                {
00528                   dum=m_matrix[ll][icol];
00529                   m_matrix[ll][icol] = 0;
00530                   for (l=0; l<4; l++) m_matrix[ll][l] -= m_matrix[icol][l]*dum;
00531                }
00532          }
00533          for (int lx=4; lx>0; --lx)
00534          {
00535             if (indxr[lx-1] != indxc[lx-1])
00536                for (k=0; k<4; k++) SGL_SWAP(m_matrix[k][indxr[lx-1]],
00537                                             m_matrix[k][indxc[lx-1]],temp);
00538          }
00539 
00540          return true;
00541       }
00542 
00543    template <class S>
00544    void inverseAffine(const sglMat4<S>& mat)
00545       {
00546          // |   R    p |'   |   R' -R'p |'
00547          // |          | -> |           |
00548          // | 0 0 0  1 |    | 0 0 0  1  |
00549          for (unsigned int i=0; i<3; i++)
00550          {
00551             m_matrix[i][3] = 0;
00552             m_matrix[3][i] = (T) -(mat[i][0]*mat[3][0] +
00553                                    mat[i][1]*mat[3][1] +
00554                                    mat[i][2]*mat[3][2]);
00555             for (unsigned int j=0; j<3; j++)
00556             {
00557                m_matrix[i][j] = (T)mat[j][i];
00558             }
00559          }
00560          m_matrix[3][3] = 1;
00561       }
00562 
00563 // stream operators
00564 #if defined(TEMPLATE_OSTREAM_HACK)
00565    template <class S>
00566    friend ostream& operator<<(ostream& ostrm, const sglMat4<S>& mat);
00567 
00568    template <class S>
00569    friend istream& operator>>(istream& istrm, sglMat4<S>& mat);
00570 #else
00571    friend ostream& operator<<(ostream& ostrm, const sglMat4<T>& mat)
00572         {
00573             ostrm << '(';
00574             for (unsigned int i=0; i<4; i++)
00575                 ostrm << '(' << mat[i][0] << ',' << mat[i][1] << ','
00576                       << mat[i][2] << ',' << mat[i][3] << ')';
00577             return (ostrm << ')');
00578         }
00579 
00580    friend istream& operator>>(istream& istrm, sglMat4<T>& mat)
00581         {
00582             for (unsigned int i=0;i<4;i++)
00583                 for (unsigned int j=0;j<4;j++)
00584                     istrm >> mat[i][j];
00585             return istrm;
00586         }
00587 #endif
00588 
00589 private:
00590    T m_matrix[4][4];
00591 };
00592 
00593 //----------------------------------------------------------------------------
00594 #if defined(TEMPLATE_OSTREAM_HACK)
00595 template <class T>
00596 inline ostream& operator<<(ostream& ostrm, const sglMat4<T>& mat)
00597 {
00598     return (ostrm << '(' << mat[0] << ", " << mat[1]
00599             << ", " << mat[2] << ", " << mat[3] << ')');
00600 }
00601 
00602 template <class T>
00603 inline istream& operator>>(istream& istrm, sglMat4<T>& mat)
00604 {
00605     for (unsigned int i=0;i<4;i++)
00606         for (unsigned int j=0;j<4;j++)
00607             istrm >> mat[i][j];
00608     return istrm;
00609 }
00610 #endif
00611 //----------------------------------------------------------------------------
00612 
00613 
00614 //----------------------------------------------------------------------------
00615 // vector/matrix4x4 operations
00616 //----------------------------------------------------------------------------
00617 
00618 //----------------------------------------------------------------------------
00619 // v' = v' * m (or v = m' * v)
00620 //----------------------------------------------------------------------------
00621 template <class S, class T>
00622 inline void mul(sglVec4<T>& vec, const sglMat4<S>& mat)
00623 {
00624    T x = (T)(vec[0]*mat[0][0] + vec[1]*mat[1][0] +
00625                                          vec[2]*mat[2][0] + vec[3]*mat[3][0]);
00626    T y = (T)(vec[0]*mat[0][1] + vec[1]*mat[1][1] +
00627                                          vec[2]*mat[2][1] + vec[3]*mat[3][1]);
00628    T z = (T)(vec[0]*mat[0][2] + vec[1]*mat[1][2] +
00629                                          vec[2]*mat[2][2] + vec[3]*mat[3][2]);
00630    vec[3] = (T)(vec[0]*mat[0][3] + vec[1]*mat[1][3] +
00631                                          vec[2]*mat[2][3] + vec[3]*mat[3][3]);
00632    vec[0] = x; vec[1] = y; vec[2] = z;
00633 }
00634 
00635 //----------------------------------------------------------------------------
00636 // dest' = v' * m (or dest = m' * v)
00637 //----------------------------------------------------------------------------
00638 template <class R, class S, class T>
00639 inline void mul(const sglVec4<R>& vec, const sglMat4<S>& mat,
00640                 sglVec4<T>& dest)
00641 {
00642    T x = (T)(vec[0]*mat[0][0] + vec[1]*mat[1][0] + vec[2]*mat[2][0] +
00643              vec[3]*mat[3][0]);
00644    T y = (T)(vec[0]*mat[0][1] + vec[1]*mat[1][1] + vec[2]*mat[2][1] +
00645              vec[3]*mat[3][1]);
00646    T z = (T)(vec[0]*mat[0][2] + vec[1]*mat[1][2] + vec[2]*mat[2][2] +
00647              vec[3]*mat[3][2]);
00648    dest[3] = (T)(vec[0]*mat[0][3] + vec[1]*mat[1][3] + vec[2]*mat[2][3] +
00649                  vec[3]*mat[3][3]);
00650    dest[0] = x; dest[1] = y; dest[2] = z;
00651 }
00652 
00653 //----------------------------------------------------------------------------
00654 // dest' = [v[0] v[1] v[2] w] * m (or dest = m' * v)
00655 //----------------------------------------------------------------------------
00656 template <class R, class S, class T>
00657 inline void mul(const sglVec3<R>& vec, R w, const sglMat4<S>& mat,
00658                 sglVec4<T>& dest)
00659 {
00660    dest[0] = (T)(vec[0]*mat[0][0] + vec[1]*mat[1][0] + vec[2]*mat[2][0] +
00661                  w*mat[3][0]);
00662    dest[1] = (T)(vec[0]*mat[0][1] + vec[1]*mat[1][1] + vec[2]*mat[2][1] +
00663                  w*mat[3][1]);
00664    dest[2] = (T)(vec[0]*mat[0][2] + vec[1]*mat[1][2] + vec[2]*mat[2][2] +
00665                  w*mat[3][2]);
00666    dest[3] = (T)(vec[0]*mat[0][3] + vec[1]*mat[1][3] + vec[2]*mat[2][3] +
00667                  w*mat[3][3]);
00668 }
00669 
00670 //----------------------------------------------------------------------------
00671 // v' = v' * m' (or v = m * v)
00672 //----------------------------------------------------------------------------
00673 template <class S, class T>
00674 inline void mulTranspose(sglVec4<T>& vec, const sglMat4<S>& mat)
00675 {
00676    T x =    (T)(vec[0]*mat[0][0] + vec[1]*mat[0][1] + vec[2]*mat[0][2] +
00677                 vec[3]*mat[0][3]);
00678    T y =    (T)(vec[0]*mat[1][0] + vec[1]*mat[1][1] + vec[2]*mat[1][2] +
00679                 vec[3]*mat[1][3]);
00680    T z =    (T)(vec[0]*mat[2][0] + vec[1]*mat[2][1] + vec[2]*mat[2][2] +
00681                 vec[3]*mat[2][3]);
00682    vec[3] = (T)(vec[0]*mat[3][0] + vec[1]*mat[3][1] + vec[2]*mat[3][2] +
00683                 vec[3]*mat[3][3]);
00684    vec[0] = x; vec[1] = y; vec[2] = z;
00685 }
00686 
00687 //----------------------------------------------------------------------------
00688 // dest' = v' * m' (or dest = m * v)
00689 //----------------------------------------------------------------------------
00690 template <class R, class S, class T>
00691 inline void mulTranspose(const sglVec4<R> &vec, const sglMat4<S>& mat,
00692                          sglVec4<T>& dest)
00693 {
00694    T x = (T)(vec[0]*mat[0][0] + vec[1]*mat[0][1] + vec[2]*mat[0][2] +
00695              vec[3]*mat[0][3]);
00696    T y = (T)(vec[0]*mat[1][0] + vec[1]*mat[1][1] + vec[2]*mat[1][2] +
00697              vec[3]*mat[1][3]);
00698    T z = (T)(vec[0]*mat[2][0] + vec[1]*mat[2][1] + vec[2]*mat[2][2] +
00699              vec[3]*mat[2][3]);
00700    dest[3] = (T)(vec[0]*mat[3][0] + vec[1]*mat[3][1] + vec[2]*mat[3][2] +
00701                      vec[3]*mat[3][3]);
00702    dest[0] = x; dest[1] = y; dest[2] = z;
00703 }
00704 
00705 //----------------------------------------------------------------------------
00706 //  dest' = [v[0] v[1] v[2] w] * m' (or dest = m * v)
00707 //----------------------------------------------------------------------------
00708 template <class R, class S, class T>
00709 inline void mulTranspose(const sglVec3<R>& vec, R w, const sglMat4<S>& mat,
00710                          sglVec4<T>& dest)
00711 {
00712    dest[0] = (T)(vec[0]*mat[0][0] + vec[1]*mat[0][1] + vec[2]*mat[0][2] +
00713                      w*mat[0][3]);
00714    dest[1] = (T)(vec[0]*mat[1][0] + vec[1]*mat[1][1] + vec[2]*mat[1][2] +
00715                      w*mat[1][3]);
00716    dest[2] = (T)(vec[0]*mat[2][0] + vec[1]*mat[2][1] + vec[2]*mat[2][2] +
00717                      w*mat[2][3]);
00718    dest[3] = (T)(vec[0]*mat[3][0] + vec[1]*mat[3][1] + vec[2]*mat[3][2] +
00719                      w*mat[3][3]);
00720 }
00721 
00722 #endif