#include <cstdlib>
#include <iostream>
#include <math.h>
using namespace std;
struct Matrix
{
float e[16];
inline float& operator[](int index) {return e[index];}
inline const float& operator[](int index) const {return e[index];}
};
struct Vector
{
float e[4];
inline float& operator[](int index) {return e[index];}
inline const float& operator[](int index) const {return e[index];}
};
ostream& operator<<(ostream& out, Matrix &m);
ostream& operator<<(ostream& out, Vector &v);
Matrix buildidentitymatrix(void);
Matrix buildrotationmatrix(float theta, float x, float y, float z);
Matrix buildscalematrix(float x, float y, float z);
Matrix buildtranslationmatrix(float x, float y, float z);
Vector buildvector(float x, float y, float z);
Matrix operator*(const Matrix &a, const Matrix &b);
void operator*=(Matrix &a, const Matrix &b);
Vector operator*(const Matrix &a, const Vector &v);
ostream& operator<<(ostream& out, const Matrix &m)
{
out << m[0] << " " << m[4] << " " << m[8] << " " << m[12] << endl;
out << m[1] << " " << m[5] << " " << m[9] << " " << m[13] << endl;
out << m[2] << " " << m[6] << " " << m[10] << " " << m[14] << endl;
out << m[3] << " " << m[7] << " " << m[11] << " " << m[15];
return out;
}
ostream& operator<<(ostream& out, Vector &v)
{
out << "(" << v[0] << ", " << v[1] << ", " << v[2] << ")";
return out;
}
Matrix buildidentitymatrix(void)
{
Matrix m;
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;
return m;
}
Matrix buildrotationmatrix(float theta, float x, float y, float z)
{
Matrix m;
float c = cos(theta);
float s = sin(theta);
float t = 1 - cos(theta);
m[0] = t*x*x + c;
m[1] = t*x*y + s*z;
m[2] = t*x*z - s*y;
m[3] = 0;
m[4] = t*x*y - s*z;
m[5] = t*y*y + c;
m[6] = t*y*z + s*x;
m[7] = 0;
m[8] = t*x*z + s*y;
m[9] = t*y*z - s*x;
m[10] = t*z*z + c;
m[11] = 0;
m[12] = 0;
m[13] = 0;
m[14] = 0;
m[15] = 1;
return m;
}
Matrix buildscalematrix(float x, float y, float z)
{
Matrix m;
m[0] = x;
m[1] = 0;
m[2] = 0;
m[3] = 0;
m[4] = 0;
m[5] = y;
m[6] = 0;
m[7] = 0;
m[8] = 0;
m[9] = 0;
m[10] = z;
m[11] = 0;
m[12] = 0;
m[13] = 0;
m[14] = 0;
m[15] = 1;
return m;
}
Matrix buildtranslationmatrix(float x, float y, float z)
{
Matrix m;
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] = x;
m[13] = y;
m[14] = z;
m[15] = 1;
return m;
}
Matrix operator*(const Matrix &a, const Matrix &b)
{
Matrix m;
m[0] = b[0]*a[0] + b[1]*a[4] + b[2]*a[8] + b[3]*a[12];
m[1] = b[0]*a[1] + b[1]*a[5] + b[2]*a[9] + b[3]*a[13];
m[2] = b[0]*a[2] + b[1]*a[6] + b[2]*a[10] + b[3]*a[14];
m[3] = b[0]*a[3] + b[1]*a[7] + b[2]*a[11] + b[3]*a[15];
m[4] = b[4]*a[0] + b[5]*a[4] + b[6]*a[8] + b[7]*a[12];
m[5] = b[4]*a[1] + b[5]*a[5] + b[6]*a[9] + b[7]*a[13];
m[6] = b[4]*a[2] + b[5]*a[6] + b[6]*a[10] + b[7]*a[14];
m[7] = b[4]*a[3] + b[5]*a[7] + b[6]*a[11] + b[7]*a[15];
m[8] = b[8]*a[0] + b[9]*a[4] + b[10]*a[8] + b[11]*a[12];
m[9] = b[8]*a[1] + b[9]*a[5] + b[10]*a[9] + b[11]*a[13];
m[10] = b[8]*a[2] + b[9]*a[6] + b[10]*a[10] + b[11]*a[14];
m[11] = b[8]*a[3] + b[9]*a[7] + b[10]*a[11] + b[11]*a[15];
m[12] = b[12]*a[0] + b[13]*a[4] + b[14]*a[8] + b[15]*a[12];
m[13] = b[12]*a[1] + b[13]*a[5] + b[14]*a[9] + b[15]*a[13];
m[14] = b[12]*a[2] + b[13]*a[6] + b[14]*a[10] + b[15]*a[14];
m[15] = b[12]*a[3] + b[13]*a[7] + b[14]*a[11] + b[15]*a[15];
return m;
}
void operator*=(Matrix &a, const Matrix &b)
{
a = a*b;
}
Vector buildvector(float x, float y, float z)
{
Vector v;
v[0] = x;
v[1] = y;
v[2] = z;
v[3] = 1;
return v;
}
Vector operator*(const Matrix &a, const Vector &v)
{
Vector vnew;
vnew[0] = v[0]*a[0] + v[1]*a[4] + v[2]*a[8] + v[3]*a[12];
vnew[1] = v[0]*a[1] + v[1]*a[5] + v[2]*a[9] + v[3]*a[13];
vnew[2] = v[0]*a[2] + v[1]*a[6] + v[2]*a[10] + v[3]*a[14];
vnew[3] = v[0]*a[3] + v[1]*a[7] + v[2]*a[11] + v[3]*a[15];
return vnew;
}