mat4.h

Go to the documentation of this file.

/*
 * $Id: mat4.h,v 1.6 2003/10/09 08:27:51 jmk Exp $
 */
///////////////////////////////////////////////////////////////////////////
//              _____________  ______________________    __   __  _____
//             /  ________  |  |   ___   ________   /   |  \ /  \   |
//            |  |       |  |_ |  |_ |  |       /  /    \__  |      |
//            |  |  ___  |  || |  || |  |      /  /        | |      |
//            |  | |   \ |  || |  || |  |     /  /      \__/ \__/ __|__
//            |  | |_@  ||  || |  || |  |    /  /          Institute
//            |  |___/  ||  ||_|  || |  |   /  /_____________________
//             \_______/  \______/ | |__|  /___________________________
//                        |  |__|  |
//                         \______/
//                    University of Utah       
//                           2002
///////////////////////////////////////////////////////////////////////////

// mat4.h

#ifndef mat4_h
#define mat4_h

#include "mm.h"
#include "mat3.h"

namespace gutz
{
template <class T>
class vec3;

template <class T>
class vec4;

template <class T>
class quat;

template <class T>
class mat4
{  
public:
  // defauts to identity  
  inline mat4 ();
    inline mat4 (const mat4<T>& M);
    inline mat4 (const T);
    inline mat4 (T*);
    inline mat4 (T, T, T, T,
                                 T, T, T, T,
                                 T, T, T, T,
                                 T, T, T, T);
    inline mat4 (const vec4<T>&, const vec4<T>&,
                                 const vec4<T>&, const vec4<T>&);
    // inline mat4 (const mat4<T>&); use default
    inline mat4 (const mat3<T>&, const vec3<T>&);
    inline ~mat4 () {}

    // inline mat4& operator= (const mat4&); use default
    inline mat4& operator+= (const mat4&);
    inline mat4& operator-= (const mat4&);
    inline mat4& operator*= (const T);
    inline mat4& operator/= (const T);
    inline mat4& operator*= (const mat4&);
  
    inline mat4 operator+ (const mat4&) const;
    inline mat4 operator- (const mat4&) const;
    inline mat4 operator- () const;
  
    inline mat4 operator* (const T) const;
    inline mat4 operator/ (const T) const;
    inline mat4 operator* (const mat4&) const;
    inline vec4<T> operator* (const vec4<T>&) const;

    inline T* v ();
        inline const T* v () const;
    inline void set(const T, const T, const T, const T,
                            const T, const T, const T, const T,
                            const T, const T, const T, const T,
                            const T, const T, const T, const T);
    inline T& operator[] (const int);
    inline T  operator[] (const int) const;

    inline vec4<T>& x() const;
    inline vec4<T>& y() const;
    inline vec4<T>& z() const;
    inline vec4<T>& w() const;

    inline void zero ();
    inline void eye ();
  
    mat4 inv () const;
    mat4 rinv () const;

    inline mat4 tp () const;

    inline mat3<T> rot () const;
    inline void rot (const mat3<T>&);

    inline vec3<T> trans () const;
    inline void trans (const vec3<T>&);

    inline vec3<T> tpoint (const vec3<T>&) const;
    inline vec3<T> tdir (const vec3<T>&) const;

public:
    T m[16];
};

//-----------------------------------------------------------------------------
// Typedefs
//-----------------------------------------------------------------------------

typedef mat4<uchar>         mat4ub;
typedef mat4<int>              mat4i;
typedef mat4<uint>          mat4ui;
typedef mat4<float>         mat4f;
typedef mat4<double>        mat4d;

//-----------------------------------------------------------------------------
// External functions
//-----------------------------------------------------------------------------
template <class T>
inline
mat4<T>
operator* (const T, const mat4<T>&);

template <class T>
mat4<T>
diag (const vec4<T>&); // diagonal matrix

template <class T>
std::ostream&
operator<< (std::ostream& os, const mat4<T>& M)
{
    os << M.m[ 0] << " " << M.m[ 4] << " " 
        << M.m[ 8] << " " << M.m[12] << "\n"
       << M.m[ 1] << " " << M.m[ 5] << " "
        << M.m[ 9] << " " << M.m[13] << "\n"
       << M.m[ 2] << " " << M.m[ 6] << " "
             << M.m[10] << " " << M.m[14] << "\n"
       << M.m[ 3] << " " << M.m[ 7] << " "
        << M.m[11] << " " << M.m[15];

    return os;
}

template <class T>
std::istream&
operator>> (std::istream& is, mat4<T>& M)
{
   is >> M.m[ 0] >> M.m[ 4] >> M.m[ 8] >> M.m[12];
   is >> M.m[ 1] >> M.m[ 5] >> M.m[ 9] >> M.m[13];
   is >> M.m[ 2] >> M.m[ 6] >> M.m[10] >> M.m[14];
   is >> M.m[ 3] >> M.m[ 7] >> M.m[11] >> M.m[15];
   return is;
}

template <class T>
bool
str2mat4 (const std::string& s, mat4<T>& v)
{
    std::istringstream z(s);
    mat4<T> w;
    bool status = z.eof();
    for (int i=0; i<15; i++)
    {
        z >> w.m[i];
        status = status || z.fail() || z.eof();
    }
    z >> w.m[15];
    status = status || z.fail();
    if (!status) v = w.tp();
    return status;
}

template <class T>
bool
str2mat4 (const std::string*, mat4<T>*);

//-----------------------------------------------------------------------------
// Implementation
//-----------------------------------------------------------------------------

template <class T>
mat4<T>::mat4 ()
{
    m[ 0] = 1; m[ 1] = 0; m[ 2] = 0; m[ 3] = 0;
    m[ 4] = 0; m[ 5] = 1; m[ 6] = 0; m[ 7] = 0;
    m[ 8] = 0; m[ 9] = 0; m[10] = 1; m[11] = 0;
    m[12] = 0; m[13] = 0; m[14] = 0; m[15] = 1;
}

template <class T>
mat4<T>::mat4 (const T c)
{
    m[ 0] = c; m[ 1] = c; m[ 2] = c; m[ 3] = c;
    m[ 4] = c; m[ 5] = c; m[ 6] = c; m[ 7] = c;
    m[ 8] = c; m[ 9] = c; m[10] = c; m[11] = c;
    m[12] = c; m[13] = c; m[14] = c; m[15] = c;
}

template <class T>
mat4<T>::mat4 (T* m0)
{
    m[ 0] = m0[ 0]; m[ 1] = m0[ 1]; m[ 2] = m0[ 2]; m[ 3] = m0[ 3];
    m[ 4] = m0[ 4]; m[ 5] = m0[ 5]; m[ 6] = m0[ 6]; m[ 7] = m0[ 7];
    m[ 8] = m0[ 8]; m[ 9] = m0[ 9]; m[10] = m0[10]; m[11] = m0[11];
    m[12] = m0[12]; m[13] = m0[13]; m[14] = m0[14]; m[15] = m0[15];
}

template <class T>
mat4<T>::mat4 (T m0, T m1, T m2, T m3,
                    T m4, T m5, T m6, T m7,
                    T m8, T m9, T m10, T m11,
                    T m12, T m13, T m14, T m15)
{
    m[ 0] = m0; m[ 1] = m1; m[ 2] = m2; m[ 3] = m3;
    m[ 4] = m4; m[ 5] = m5; m[ 6] = m6; m[ 7] = m7;
    m[ 8] = m8; m[ 9] = m9; m[10] = m10; m[11] = m11;
    m[12] = m12; m[13] = m13; m[14] = m14; m[15] = m15;
}

template <class T>
mat4<T>::mat4 (const vec4<T>& c0, const vec4<T>& c1,
               const vec4<T>& c2, const vec4<T>& c3)
{
    m[ 0] = c0.x; m[ 1] = c0.y; m[ 2] = c0.z; m[ 3] = c0.w;
    m[ 4] = c1.x; m[ 5] = c1.y; m[ 6] = c1.z; m[ 7] = c1.w;
    m[ 8] = c2.x; m[ 9] = c2.y; m[10] = c2.z; m[11] = c2.w;
    m[12] = c3.x; m[13] = c3.y; m[14] = c3.z; m[15] = c3.w;
}

template <class T>
mat4<T>::mat4(const mat3<T>& M, const vec3<T>& v)
{
    m[ 0] = M.m[0]; m[ 1] = M.m[1]; m[ 2] = M.m[2]; m[ 3] = 0;
    m[ 4] = M.m[3]; m[ 5] = M.m[4]; m[ 6] = M.m[5]; m[ 7] = 0;
    m[ 8] = M.m[6]; m[ 9] = M.m[7]; m[10] = M.m[8]; m[11] = 0;
    m[12] = v.x; m[13] = v.y; m[14] = v.z; m[15] = 1;
}

template <class T>
mat4<T>::mat4 (const mat4<T>& M)
{
    m[ 0] = M.m[ 0]; m[ 1] = M.m[ 1]; m[ 2] = M.m[ 2]; m[ 3] = M.m[ 3];
    m[ 4] = M.m[ 4]; m[ 5] = M.m[ 5]; m[ 6] = M.m[ 6]; m[ 7] = M.m[ 7];
    m[ 8] = M.m[ 8]; m[ 9] = M.m[ 9]; m[10] = M.m[10]; m[11] = M.m[11];
    m[12] = M.m[12]; m[13] = M.m[13]; m[14] = M.m[14]; m[15] = M.m[15];
}

template <class T>
mat4<T>&
mat4<T>::operator+= (const mat4<T>& M)
{    
    m[ 0] += M.m[ 0]; m[ 1] += M.m[ 1]; m[ 2] += M.m[ 2]; m[ 3] += M.m[ 3];
    m[ 4] += M.m[ 4]; m[ 5] += M.m[ 5]; m[ 6] += M.m[ 6]; m[ 7] += M.m[ 7];
    m[ 8] += M.m[ 8]; m[ 9] += M.m[ 9]; m[10] += M.m[10]; m[11] += M.m[11];
    m[12] += M.m[12]; m[13] += M.m[13]; m[14] += M.m[14]; m[15] += M.m[15];
    return *this;
}

template <class T>
mat4<T>&
mat4<T>::operator-= (const mat4<T>& M)
{
    m[ 0] -= M.m[ 0]; m[ 1] -= M.m[ 1]; m[ 2] -= M.m[ 2]; m[ 3] -= M.m[ 3];
    m[ 4] -= M.m[ 4]; m[ 5] -= M.m[ 5]; m[ 6] -= M.m[ 6]; m[ 7] -= M.m[ 7];
    m[ 8] -= M.m[ 8]; m[ 9] -= M.m[ 9]; m[10] -= M.m[10]; m[11] -= M.m[11];
    m[12] -= M.m[12]; m[13] -= M.m[13]; m[14] -= M.m[14]; m[15] -= M.m[15];
    return *this;
}

template <class T>
mat4<T>&
mat4<T>::operator*= (T c)
{
    m[ 0] *= c; m[ 1] *= c; m[ 2] *= c; m[ 3] *= c;
    m[ 4] *= c; m[ 5] *= c; m[ 6] *= c; m[ 7] *= c;
    m[ 8] *= c; m[ 9] *= c; m[10] *= c; m[11] *= c;
    m[12] *= c; m[13] *= c; m[14] *= c; m[15] *= c;
    return *this;
}

template <class T>
mat4<T>&
mat4<T>::operator/= (T c)
{
    T d = 1/c;
    m[ 0] *= d; m[ 1] *= d; m[ 2] *= d; m[ 3] *= d;
    m[ 4] *= d; m[ 5] *= d; m[ 6] *= d; m[ 7] *= d;
    m[ 8] *= d; m[ 9] *= d; m[10] *= d; m[11] *= d;
    m[12] *= d; m[13] *= d; m[14] *= d; m[15] *= d;
    return *this;
}

template <class T>
mat4<T>&
mat4<T>::operator*= (const mat4<T>& M)
{
    mat4<T> temp(*this);

    m[ 0] = 
        temp.m[ 0] * M.m[ 0] + temp.m[ 4] * M.m[ 1] + 
        temp.m[ 8] * M.m[ 2] + temp.m[12] * M.m[ 3];
    m[ 1] = 
        temp.m[ 1] * M.m[ 0] + temp.m[ 5] * M.m[ 1] + 
        temp.m[ 9] * M.m[ 2] + temp.m[13] * M.m[ 3];
    m[ 2] = 
        temp.m[ 2] * M.m[ 0] + temp.m[ 6] * M.m[ 1] + 
        temp.m[10] * M.m[ 2] + temp.m[14] * M.m[ 3];
    m[ 3] = 
        temp.m[ 3] * M.m[ 0] + temp.m[ 7] * M.m[ 1] + 
        temp.m[11] * M.m[ 2] + temp.m[15] * M.m[ 3];

    m[ 4] = 
        temp.m[ 0] * M.m[ 4] + temp.m[ 4] * M.m[ 5] + 
        temp.m[ 8] * M.m[ 6] + temp.m[12] * M.m[ 7];
    m[ 5] = 
        temp.m[ 1] * M.m[ 4] + temp.m[ 5] * M.m[ 5] + 
        temp.m[ 9] * M.m[ 6] + temp.m[13] * M.m[ 7];
    m[ 6] = 
        temp.m[ 2] * M.m[ 4] + temp.m[ 6] * M.m[ 5] + 
        temp.m[10] * M.m[ 6] + temp.m[14] * M.m[ 7];
    m[ 7] = 
        temp.m[ 3] * M.m[ 4] + temp.m[ 7] * M.m[ 5] + 
        temp.m[11] * M.m[ 6] + temp.m[15] * M.m[ 7];

    m[ 8] = 
        temp.m[ 0] * M.m[ 8] + temp.m[ 4] * M.m[ 9] + 
        temp.m[ 8] * M.m[10] + temp.m[12] * M.m[11];
    m[ 9] = 
        temp.m[ 1] * M.m[ 8] + temp.m[ 5] * M.m[ 9] + 
        temp.m[ 9] * M.m[10] + temp.m[13] * M.m[11];
    m[10] = 
        temp.m[ 2] * M.m[ 8] + temp.m[ 6] * M.m[ 9] + 
        temp.m[10] * M.m[10] + temp.m[14] * M.m[11];
    m[11] = 
        temp.m[ 3] * M.m[ 8] + temp.m[ 7] * M.m[ 9] + 
        temp.m[11] * M.m[10] + temp.m[15] * M.m[11];

    m[12] = 
        temp.m[ 0] * M.m[12] + temp.m[ 4] * M.m[13] + 
        temp.m[ 8] * M.m[14] + temp.m[12] * M.m[15];
    m[13] = 
        temp.m[ 1] * M.m[12] + temp.m[ 5] * M.m[13] + 
        temp.m[ 9] * M.m[14] + temp.m[13] * M.m[15];
    m[14] = 
        temp.m[ 2] * M.m[12] + temp.m[ 6] * M.m[13] + 
        temp.m[10] * M.m[14] + temp.m[14] * M.m[15];
    m[15] = 
        temp.m[ 3] * M.m[12] + temp.m[ 7] * M.m[13] + 
        temp.m[11] * M.m[14] + temp.m[15] * M.m[15];

    return *this;
}

template <class T>
mat4<T>
mat4<T>::operator+ (const mat4<T>& M) const
{
    return mat4<T>(m[ 0] + M.m[ 0], m[ 1] + M.m[ 1],
                   m[ 2] + M.m[ 2], m[ 3] + M.m[ 3],
                   m[ 4] + M.m[ 4], m[ 5] + M.m[ 5],
                   m[ 6] + M.m[ 6], m[ 7] + M.m[ 7],
                   m[ 8] + M.m[ 8], m[ 9] + M.m[ 9],
                   m[10] + M.m[10], m[11] + M.m[11],
                   m[12] + M.m[12], m[13] + M.m[13],
                   m[14] + M.m[14], m[15] + M.m[15]);
}

template <class T>
mat4<T>
mat4<T>::operator- (const mat4<T>& M) const
{
    return mat4<T>(m[ 0] - M.m[ 0], m[ 1] - M.m[ 1],
                   m[ 2] - M.m[ 2], m[ 3] - M.m[ 3],
                   m[ 4] - M.m[ 4], m[ 5] - M.m[ 5],
                   m[ 6] - M.m[ 6], m[ 7] - M.m[ 7],
                   m[ 8] - M.m[ 8], m[ 9] - M.m[ 9],
                   m[10] - M.m[10], m[11] - M.m[11],
                   m[12] - M.m[12], m[13] - M.m[13],
                   m[14] - M.m[14], m[15] - M.m[15]);
}

template <class T>
mat4<T>
mat4<T>::operator- () const
{
    return mat4<T>(-m[ 0], -m[ 1], -m[ 2], -m[ 3],
                   -m[ 4], -m[ 5], -m[ 6], -m[ 7],
                   -m[ 8], -m[ 9], -m[10], -m[11],
                   -m[12], -m[13], -m[14], -m[15]);
}

template <class T>
mat4<T>
mat4<T>::operator* (T c) const
{
    return mat4<T>(m[ 0]*c, m[ 1]*c, m[ 2]*c, m[ 3]*c,
                   m[ 4]*c, m[ 5]*c, m[ 6]*c, m[ 7]*c,
                   m[ 8]*c, m[ 9]*c, m[10]*c, m[11]*c,
                   m[12]*c, m[13]*c, m[14]*c, m[15]*c);
}

template <class T>
mat4<T>
mat4<T>::operator/ (T c) const
{
    T d = 1/c;
    return mat4<T>(m[ 0]*d, m[ 1]*d, m[ 2]*d, m[ 3]*d,
                   m[ 4]*d, m[ 5]*d, m[ 6]*d, m[ 7]*d,
                   m[ 8]*d, m[ 9]*d, m[10]*d, m[11]*d,
                   m[12]*d, m[13]*d, m[14]*d, m[15]*d);
}

template <class T>
mat4<T>
mat4<T>::operator* (const mat4<T>& M) const
{
    return mat4<T>(m[ 0]*M.m[ 0] + m[ 4]*M.m[ 1] + m[ 8]*M.m[ 2] + m[12]*M.m[ 3],
                   m[ 1]*M.m[ 0] + m[ 5]*M.m[ 1] + m[ 9]*M.m[ 2] + m[13]*M.m[ 3],
                   m[ 2]*M.m[ 0] + m[ 6]*M.m[ 1] + m[10]*M.m[ 2] + m[14]*M.m[ 3],
                   m[ 3]*M.m[ 0] + m[ 7]*M.m[ 1] + m[11]*M.m[ 2] + m[15]*M.m[ 3],
                    
                   m[ 0]*M.m[ 4] + m[ 4]*M.m[ 5] + m[ 8]*M.m[ 6] + m[12]*M.m[ 7],
                   m[ 1]*M.m[ 4] + m[ 5]*M.m[ 5] + m[ 9]*M.m[ 6] + m[13]*M.m[ 7],
                   m[ 2]*M.m[ 4] + m[ 6]*M.m[ 5] + m[10]*M.m[ 6] + m[14]*M.m[ 7],
                   m[ 3]*M.m[ 4] + m[ 7]*M.m[ 5] + m[11]*M.m[ 6] + m[15]*M.m[ 7],
                    
                   m[ 0]*M.m[ 8] + m[ 4]*M.m[ 9] + m[ 8]*M.m[10] + m[12]*M.m[11],
                   m[ 1]*M.m[ 8] + m[ 5]*M.m[ 9] + m[ 9]*M.m[10] + m[13]*M.m[11],
                   m[ 2]*M.m[ 8] + m[ 6]*M.m[ 9] + m[10]*M.m[10] + m[14]*M.m[11],
                   m[ 3]*M.m[ 8] + m[ 7]*M.m[ 9] + m[11]*M.m[10] + m[15]*M.m[11],
                    
                   m[ 0]*M.m[12] + m[ 4]*M.m[13] + m[ 8]*M.m[14] + m[12]*M.m[15],
                   m[ 1]*M.m[12] + m[ 5]*M.m[13] + m[ 9]*M.m[14] + m[13]*M.m[15],
                   m[ 2]*M.m[12] + m[ 6]*M.m[13] + m[10]*M.m[14] + m[14]*M.m[15],
                   m[ 3]*M.m[12] + m[ 7]*M.m[13] + m[11]*M.m[14] + m[15]*M.m[15]);
}

template <class T>
vec4<T>
mat4<T>::operator* (const vec4<T>& v) const
{
    return vec4<T>(     m[ 0]*v.x + m[ 4]*v.y + m[ 8]*v.z + m[12]*v.w,
                                        m[ 1]*v.x + m[ 5]*v.y + m[ 9]*v.z + m[13]*v.w,
                                        m[ 2]*v.x + m[ 6]*v.y + m[10]*v.z + m[14]*v.w,
                                        m[ 3]*v.x + m[ 7]*v.y + m[11]*v.z + m[15]*v.w);
}

template <class T>
T*
mat4<T>::v ()
{
    return m;
}

template <class T>
const T*
mat4<T>::v () const
{
    return m;
}

template <class T>
void
mat4<T>::set (const T m0, const T m1, const T m2, const T m3,
              const T m4, const T m5, const T m6, const T m7,
              const T m8, const T m9, const T m10, const T m11,
              const T m12, const T m13, const T m14, const T m15)
{
    m[ 0] = m0; m[ 1] = m1; m[ 2] = m2; m[ 3] = m3;
    m[ 4] = m4; m[ 5] = m5; m[ 6] = m6; m[ 7] = m7;
    m[ 8] = m8; m[ 9] = m9; m[10] = m10; m[11] = m11;
    m[12] = m12; m[13] = m13; m[14] = m14; m[15] = m15;
}

template <class T>
T&
mat4<T>::operator[] (const int i)
{
    return v()[i];
}

template <class T>
T
mat4<T>::operator[] (const int i) const
{
    return v()[i];
}

template <class T>
void
mat4<T>::zero ()
{
    m[ 0] = 0; m[ 1] = 0; m[ 2] = 0; m[ 3] = 0;
    m[ 4] = 0; m[ 5] = 0; m[ 6] = 0; m[ 7] = 0;
    m[ 8] = 0; m[ 9] = 0; m[10] = 0; m[11] = 0;
    m[12] = 0; m[13] = 0; m[14] = 0; m[15] = 0;
}

template <class T>
void
mat4<T>::eye ()
{
    m[ 0] = 1; m[ 1] = 0; m[ 2] = 0; m[ 3] = 0;
    m[ 4] = 0; m[ 5] = 1; m[ 6] = 0; m[ 7] = 0;
    m[ 8] = 0; m[ 9] = 0; m[10] = 1; m[11] = 0;
    m[12] = 0; m[13] = 0; m[14] = 0; m[15] = 1;
}

template <class T>
vec4<T>&
mat4<T>::x () const
{
    return (vec4<T>&)m[0];
}

template <class T>
vec4<T>&
mat4<T>::y () const
{
    return (vec4<T>&)m[4];
}

template <class T>
vec4<T>&
mat4<T>::z () const
{
    return (vec4<T>&)m[8];
}

template <class T>
vec4<T>&
mat4<T>::w () const
{
    return (vec4<T>&)m[12];
}

template <class T>
mat4<T>
mat4<T>::inv () const
{
    T det00 = 
        m[ 5]*(m[10]*m[15] - m[11]*m[14]) - 
        m[ 9]*(m[ 6]*m[15] - m[ 7]*m[14]) + 
        m[13]*(m[ 6]*m[11] - m[ 7]*m[10]);
    T det01 = 
        m[ 1]*(m[10]*m[15] - m[11]*m[14]) - 
        m[ 9]*(m[ 2]*m[15] - m[ 3]*m[14]) + 
        m[13]*(m[ 2]*m[11] - m[ 3]*m[10]);
    T det02 = 
        m[ 1]*(m[ 6]*m[15] - m[ 7]*m[14]) - 
        m[ 5]*(m[ 2]*m[15] - m[ 3]*m[14]) + 
        m[13]*(m[ 2]*m[ 7] - m[ 3]*m[ 6]);
    T det03 = 
        m[ 1]*(m[ 6]*m[11] - m[ 7]*m[10]) - 
        m[ 5]*(m[ 2]*m[11] - m[ 3]*m[10]) + 
        m[ 9]*(m[ 2]*m[ 7] - m[ 3]*m[ 6]);
  
    T det10 = 
        m[ 4]*(m[10]*m[15] - m[11]*m[14]) - 
        m[ 8]*(m[ 6]*m[15] - m[ 7]*m[14]) + 
        m[12]*(m[ 6]*m[11] - m[ 7]*m[10]);
    T det11 = 
        m[ 0]*(m[10]*m[15] - m[11]*m[14]) - 
        m[ 8]*(m[ 2]*m[15] - m[ 3]*m[14]) + 
        m[12]*(m[ 2]*m[11] - m[ 3]*m[10]);
    T det12 = 
        m[ 0]*(m[ 6]*m[15] - m[ 7]*m[14]) - 
        m[ 4]*(m[ 2]*m[15] - m[ 3]*m[14]) + 
        m[12]*(m[ 2]*m[ 7] - m[ 3]*m[ 6]);
    T det13 = 
        m[ 0]*(m[ 6]*m[11] - m[ 7]*m[10]) - 
        m[ 4]*(m[ 2]*m[11] - m[ 3]*m[10]) + 
        m[ 8]*(m[ 2]*m[ 7] - m[ 3]*m[ 6]);
  
    T det20 = 
        m[ 4]*(m[ 9]*m[15] - m[11]*m[13]) - 
        m[ 8]*(m[ 5]*m[15] - m[ 7]*m[13]) + 
        m[12]*(m[ 5]*m[11] - m[ 7]*m[ 9]);
    T det21 = 
        m[ 0]*(m[ 9]*m[15] - m[11]*m[13]) - 
        m[ 8]*(m[ 1]*m[15] - m[ 3]*m[13]) + 
        m[12]*(m[ 1]*m[11] - m[ 3]*m[ 9]);
    T det22 = 
        m[ 0]*(m[ 5]*m[15] - m[ 7]*m[13]) - 
        m[ 4]*(m[ 1]*m[15] - m[ 3]*m[13]) + 
        m[12]*(m[ 1]*m[ 7] - m[ 3]*m[ 5]);
    T det23 = 
        m[ 0]*(m[ 5]*m[11] - m[ 7]*m[ 9]) - 
        m[ 4]*(m[ 1]*m[11] - m[ 3]*m[ 9]) + 
        m[ 8]*(m[ 1]*m[ 7] - m[ 3]*m[ 5]);
  
    T det30 = 
        m[ 4]*(m[ 9]*m[14] - m[10]*m[13]) - 
        m[ 8]*(m[ 5]*m[14] - m[ 6]*m[13]) + 
        m[12]*(m[ 5]*m[10] - m[ 6]*m[ 9]);
    T det31 = 
        m[ 0]*(m[ 9]*m[14] - m[10]*m[13]) - 
        m[ 8]*(m[ 1]*m[14] - m[ 2]*m[13]) + 
        m[12]*(m[ 1]*m[10] - m[ 2]*m[ 9]);
    T det32 = 
        m[ 0]*(m[ 5]*m[14] - m[ 6]*m[13]) - 
        m[ 4]*(m[ 1]*m[14] - m[ 2]*m[13]) + 
        m[12]*(m[ 1]*m[ 6] - m[ 2]*m[ 5]);
    T det33 = 
        m[ 0]*(m[ 5]*m[10] - m[ 6]*m[ 9]) - 
        m[ 4]*(m[ 1]*m[ 9] - m[ 2]*m[ 9]) + 
        m[ 8]*(m[ 1]*m[ 6] - m[ 2]*m[ 5]);

    T d = 1 / (det00*m[ 0] - det01*m[ 4] + det02*m[ 8] - det03*m[12]);

    return mat4<T>( det00*d, -det01*d,  det02*d, -det03*d,
                   -det10*d,  det11*d, -det12*d,  det13*d,
                    det20*d, -det21*d,  det22*d, -det23*d,
                   -det30*d,  det31*d, -det32*d,  det33*d);
}

template <class T>
mat4<T>
mat4<T>::rinv () const
{
    mat3<T> R = rot().tp();
    vec3<T> p = -(R*trans());
    return mat4<T>(R, p);
}

template <class T>
mat4<T>
mat4<T>::tp () const
{
    return mat4<T>(m[ 0], m[ 4], m[ 8], m[12],
                   m[ 1], m[ 5], m[ 9], m[13],
                   m[ 2], m[ 6], m[10], m[14],
                   m[ 3], m[ 7], m[11], m[15]);
}

template <class T>
mat3<T>
mat4<T>::rot () const
{
    return mat3<T>(m[ 0], m[ 1], m[ 2],
                   m[ 4], m[ 5], m[ 6],
                   m[ 8], m[ 9], m[10]);
}

template <class T>
void
mat4<T>::rot (const mat3<T>& M)
{
    m[ 0] = M.m[0]; m[ 1] = M.m[1]; m[ 2] = M.m[2];
    m[ 4] = M.m[3]; m[ 5] = M.m[4]; m[ 6] = M.m[5];
    m[ 8] = M.m[6]; m[ 9] = M.m[7]; m[10] = M.m[8];
}

template <class T>
vec3<T>
mat4<T>::trans () const
{
    return vec3<T>(m[12], m[13], m[14]);
}

template <class T>
void
mat4<T>::trans (const vec3<T>& v)
{
    m[12] = v.x; m[13] = v.y; m[14] = v.z;
}

template <class T>
vec3<T>
mat4<T>::tpoint (const vec3<T>& p) const
{
    return vec3<T>(m[ 0]*p.x + m[ 4]*p.y + m[ 8]*p.z + m[12],
                           m[ 1]*p.x + m[ 5]*p.y + m[ 9]*p.z + m[13],
                           m[ 2]*p.x + m[ 6]*p.y + m[10]*p.z + m[14]);
}

template <class T>
vec3<T>
mat4<T>::tdir (const vec3<T>& p) const
{
    return vec3<T>(m[ 0]*p.x + m[ 4]*p.y + m[ 8]*p.z,
                           m[ 1]*p.x + m[ 5]*p.y + m[ 9]*p.z,
                           m[ 2]*p.x + m[ 6]*p.y + m[10]*p.z);
}

template <class T>
mat4<T>
operator* (const T c, const mat4<T>& M)
{
    return mat4<T>(M.m[ 0]*c, M.m[ 1]*c, M.m[ 2]*c, M.m[ 3]*c,
                   M.m[ 4]*c, M.m[ 5]*c, M.m[ 6]*c, M.m[ 7]*c,
                   M.m[ 8]*c, M.m[ 9]*c, M.m[10]*c, M.m[11]*c,
                   M.m[12]*c, M.m[13]*c, M.m[14]*c, M.m[15]*c);
}

template <class T>
mat4<T>
diag (const vec4<T>& v)
{
    return mat4<T>(v.x, 0, 0, 0,
                   0, v.y, 0, 0,
                   0, 0, v.z, 0,
                   0, 0, 0, v.w);
}

template <class T>
bool
str2mat4 (const std::string* s, mat4<T>& M)
{
    if (!s) return true;
    return str2mat4(*s, M);
}

//-----------------------------------------------------------------------------

typedef mat4<float>  mat4f;
typedef mat4<double> mat4d;

#define str2mat4f str2mat4<float>
#define str2mat4d str2mat4<double>

#if __win32
inline
mat4d
mm_castf2d (mat4f M)
{
    return mat4d(M.m[ 0], M.m[ 1], M.m[ 2], M.m[ 3],
                 M.m[ 4], M.m[ 5], M.m[ 6], M.m[ 7],
                 M.m[ 8], M.m[ 9], M.m[10], M.m[11],
                 M.m[12], M.m[13], M.m[14], M.m[15]);
}

inline
mat4f
mm_castd2f (mat4d M)
{
    return mat4f((float)M.m[ 0], (float)M.m[ 1], (float)M.m[ 2], (float)M.m[ 3],
                 (float)M.m[ 4], (float)M.m[ 5], (float)M.m[ 6], (float)M.m[ 7],
                 (float)M.m[ 8], (float)M.m[ 9], (float)M.m[10], (float)M.m[11],
                 (float)M.m[12], (float)M.m[13], (float)M.m[14], (float)M.m[15]);
}
#else //__real_os
template <class T1, class T2>
mat4<T2>
mm_cast (mat4<T1> v)
{
    return mat4<T2>(M.m[ 0], M.m[ 1], M.m[ 2], M.m[ 3],
                    M.m[ 4], M.m[ 5], M.m[ 6], M.m[ 7],
                    M.m[ 8], M.m[ 9], M.m[10], M.m[11],
                    M.m[12], M.m[13], M.m[14], M.m[15]);
}
#endif

} // End of namespace gutz

#endif // mat4_h

---

Send questions, comments, and bug reports to: