src/rcrt/math/Matrix4D.h

Go to the documentation of this file.

#ifndef MATRIX4D_H_
#define MATRIX4D_H_

#include <cmath>
#include <iostream>
#include "Vec3D.h"
#include "Vec4D.h"
#include "HCoord4D.h"

namespace rcrt
{

class Matrix4D
{
        friend std::ostream& operator<<(std::ostream& out, const Matrix4D& m);
        
private:
        float values[4][4];//[row][column]
        void assign(const Matrix4D& mat);
        inline static void swap(float mat[4][4], int& i, int& ci);
        inline static void div(float mat[4][4], int& i, float& factor);
        inline static void res(float mat[4][4], int& i, int& b, float& factor);
        
public:
        Matrix4D();
        Matrix4D(float x);
        Matrix4D(const Matrix4D& mat);
        virtual ~Matrix4D();
        
        const Matrix4D& operator= (const Matrix4D& mat);
        const Matrix4D& operator= (float mat[4][4]);
        Matrix4D operator+ (const Matrix4D& mat) const;
        Matrix4D operator- (const Matrix4D& mat) const;
        Matrix4D operator* (const Matrix4D& mat) const;
        HCoord4D operator* (const HCoord4D& coords) const;
        Point3D operator* (const Point3D& p) const;
        Vec3D operator* (const Vec3D& p) const;
        Matrix4D operator* (const float f) const;
        Matrix4D operator/ (const float f) const;
        
        Vec4D operator[] (const int& i) const;
        
        Matrix4D transpose() const;
        Matrix4D inverse() const;
        
        static Matrix4D identity()
        {
                Matrix4D m;
                for (int i = 0; i < 4; ++i)
                        for (int j = 0; j < 4; ++j)
                                m.values[i][j] = i == j ? 1 : 0;
                
                return m;
        }
        
        static Matrix4D rotation(float angle, const Vec3D& axis)
        {
                Matrix4D m;
                float c = cos(angle);
                float s = sin(angle);
                float cm = 1-c;
                float xycm = axis.x()*axis.y()*cm;
                float xzcm = axis.x()*axis.z()*cm;
                float yzcm = axis.z()*axis.y()*cm;
                float zs = axis.z() * s;
                float ys = axis.y() * s;
                float xs = axis.x() * s;
                m.values[0][0] = axis.x()*axis.x()*cm + c;
                m.values[0][1] = xycm - zs;
                m.values[0][2] = xzcm + ys;
                m.values[0][3] = 0;
                
                m.values[1][0] = xycm + zs;
                m.values[1][1] = axis.y()*axis.y()*cm + c;
                m.values[1][2] = yzcm - xs;
                m.values[1][3] = 0;
                
                m.values[2][0] = xzcm - ys;
                m.values[2][1] = yzcm + xs;
                m.values[2][2] = axis.z()*axis.z()*cm + c;
                m.values[2][3] = 0;
                
                for(int i = 0; i < 3; i++)
                        m.values[3][i] = 0;
                
                m.values[3][3] = 1;
                
                return m;
        }
        
        static Matrix4D translation(const Vec3D& vec)
        {
                Matrix4D m;
                for (int i = 0; i < 3; ++i)
                        for (int j = 0; j < 3; ++j)
                                m.values[i][j] = i == j ? 1 : 0;
                m.values[0][3] = vec.x();
                m.values[1][3] = vec.y();
                m.values[2][3] = vec.z();
                m.values[3][3] = 1;
                
                return m;
        }
        
        static Matrix4D getOrthoNormalBasis(const Vec3D& n, const Vec3D& v) {
                Matrix4D onb;
                
                onb.values[0][2] = n.x();
                onb.values[1][2] = n.y();
                onb.values[2][2] = n.z();
                
                Vec3D vp = v*(v*n)/(n*n);
                Vec3D t = (v-vp).normalize();
                
                Vec3D b = (t.crossP(n)).normalize();
                //Vector3.cross(onb.w, onb.u, onb.v);
                
                onb.values[0][1] = t.x();
                onb.values[1][1] = t.y();
                onb.values[2][1] = t.z();
                
                onb.values[0][0] = b.x();
                onb.values[1][0] = b.y();
                onb.values[2][0] = b.z();
                onb.values[3][3] = 1;
                
                return onb;
        }
        
};

inline std::ostream& operator<<(std::ostream& out, const Matrix4D& m)
{
        out<<"Matrix4D:"<<std::endl;
        out<<"/ "<<m.values[0][0]<<" "<<m.values[0][1]<<" "<<m.values[0][2]<<" "<<m.values[0][3]<<" \\"<<std::endl;;
        out<<"| "<<m.values[1][0]<<" "<<m.values[1][1]<<" "<<m.values[1][2]<<" "<<m.values[1][3]<<" |"<<std::endl;;
        out<<"| "<<m.values[2][0]<<" "<<m.values[2][1]<<" "<<m.values[2][2]<<" "<<m.values[2][3]<<" |"<<std::endl;;
        out<<"\\ "<<m.values[3][0]<<" "<<m.values[3][1]<<" "<<m.values[3][2]<<" "<<m.values[3][3]<<" /"<<std::endl;;
        return out;
}

}

#endif /*MATRIX4D_H_*/

---