src/rcrt/math/Matrix4D.cpp

Go to the documentation of this file.

#include "Matrix4D.h"

using namespace std;

namespace rcrt
{


Matrix4D::Matrix4D()
{
        for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                        values[i][j] = 0;
}


Matrix4D::Matrix4D(float x)
{
        for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                        values[i][j] = x;
}


Matrix4D::Matrix4D(const Matrix4D& mat)
{
        assign(mat);
}


void Matrix4D::assign(const Matrix4D& mat)
{
        if(&mat != this){
                for (int i = 0; i < 4; i++)
                        for (int j = 0; j < 4; j++)
                                values[i][j] = mat.values[i][j];
        }
}

const Matrix4D& Matrix4D::operator= (const Matrix4D& mat)
{
        assign(mat);
        
        return *this;
}

const Matrix4D& Matrix4D::operator= (float mat[4][4])
{
        for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                        values[i][j] = mat[i][j];
                
        return *this;
}


Matrix4D Matrix4D::operator+ (const Matrix4D& mat) const
{
        Matrix4D m;
        for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                        m.values[i][j] = values[i][j] + mat.values[i][j];
        
        return m;
}


Matrix4D Matrix4D::operator- (const Matrix4D& mat) const
{
        Matrix4D m;
        for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                        m.values[i][j] = values[i][j] - mat.values[i][j];
        
        return m;
}


Matrix4D Matrix4D::operator* (const Matrix4D& mat) const
{
        Matrix4D m;
        for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                        m.values[i][j] = values[i][0] * mat.values[0][j] +
                                                values[i][1] * mat.values[1][j] +
                                                values[i][2] * mat.values[2][j] +
                                                values[i][3] * mat.values[3][j];
        
        return m;
}


HCoord4D Matrix4D::operator* (const HCoord4D& coord) const
{
        float v[4];
        for(int i = 0; i < 4; i++)
                v[i] = values[i][0] * coord.x() +
                        values[i][1] * coord.y() +
                        values[i][2] * coord.z() +
                        values[i][3] * coord.w();
        
        return HCoord4D(v[0], v[1], v[2], v[3]);
}

Point3D Matrix4D::operator* (const Point3D& p) const
{
        float v[4];
        for(int i = 0; i < 4; i++)
                v[i] = values[i][0] * p.x() +
                        values[i][1] * p.y() +
                        values[i][2] * p.z() +
                        values[i][3];
        float invW = 1/v[3];
        return Point3D(v[0]*invW, v[1]*invW, v[2]*invW);
}

Vec3D Matrix4D::operator* (const Vec3D& vec) const
{
        float v[4];
        for(int i = 0; i < 4; i++)
                v[i] = values[i][0] * vec.x() +
                        values[i][1] * vec.y() +
                        values[i][2] * vec.z();
        return Vec3D(v[0], v[1], v[2]);
}

Matrix4D Matrix4D::operator* (float f) const
{
        Matrix4D m;
        for (int i = 0; i < 4; ++i)
                for (int j = 0; j < 4; ++j)
                        m.values[i][j] = values[i][j] * f;
        
        return m;
}


Matrix4D Matrix4D::operator/ (float f) const
{
        float g = 1/f;
        
        return *this * g;
}



Vec4D Matrix4D::operator[] (const int& i) const
{
        return Vec4D(values[0][i], values[1][i], values[2][i], values[3][i]);
}


Matrix4D::~Matrix4D()
{
}


Matrix4D Matrix4D::transpose() const
{
        Matrix4D m;
        for (int i = 0; i < 4; ++i)
                for (int j = 0; j < 4; ++j)
                        m.values[i][j] = values[j][i];
        
        return m;
}

void Matrix4D::swap(float mat[4][4], int& i, int& ci)
{
        for(int k = 0; k < 4; k++){
                float temp = mat[i][k];
                mat[i][k] = mat[ci][k];
                mat[ci][k] = temp;
        }
}

void Matrix4D::div(float mat[4][4], int& i, float& factor)
{
        if(factor == 0){
                cout<<"Error could not invert Matrix"<<endl;
                exit(1);
        }
        float inv = 1.0f/factor;
        for(int k = 0; k < 4; k++){
                mat[i][k] = mat[i][k]*inv;
        }
}

void Matrix4D::res(float mat[4][4], int& i, int& b, float& factor)
{
        int j;
        for(j=0;j<4;++j)
        {
                mat[i][j] -= mat[b][j]*factor;
        }
}

Matrix4D Matrix4D::inverse() const
{
        Matrix4D ident = Matrix4D::identity();
        Matrix4D result(*this);
        for(int i=0;i<4;++i)
        {
                float max=0;
                int ci=0;
                for(int k=i;k<4;++k)
                {
                        if(fabs(result.values[k][i])>max)
                        {
                                max=fabs(result.values[k][i]);
                                ci=k;
                        }
                }
                if(max==0)
                {
                        cout<<"Error could not invert Matrix"<<endl;
                        exit(1);
                }
                swap(result.values,i,ci);
                swap(ident.values,i,ci);
                float factor=result.values[i][i];
                div(result.values,i,factor);
                div(ident.values,i,factor);
                for(int k=0;k<4;++k)
                {
                        if(k!=i)
                        {
                                factor=result.values[k][i];
                                res(result.values,k,i,factor);
                                res(ident.values,k,i,factor);
                        }
                }
        }
        for(int i=0;i<4;++i)
                for(int j=0;j<4;++j)
                        result.values[i][j]=ident.values[i][j];
        return result;
}

}

---